Actual source code: rvector.c
petsc-3.7.5 2017-01-01
2: /*
3: Provides the interface functions for vector operations that have PetscScalar/PetscReal in the signature
4: These are the vector functions the user calls.
5: */
6: #include <petsc/private/vecimpl.h> /*I "petscvec.h" I*/
7: static PetscInt VecGetSubVectorSavedStateId = -1;
10: if ((x)->map->N != (y)->map->N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths %d != %d", (x)->map->N, (y)->map->N); \
11: if ((x)->map->n != (y)->map->n) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths %d != %d", (x)->map->n, (y)->map->n);
15: PETSC_EXTERN PetscErrorCode VecValidValues(Vec vec,PetscInt argnum,PetscBool begin)
16: {
17: #if defined(PETSC_USE_DEBUG)
18: PetscErrorCode ierr;
19: PetscInt n,i;
20: const PetscScalar *x;
23: #if defined(PETSC_HAVE_CUSP)
24: if ((vec->petscnative || vec->ops->getarray) && (vec->valid_GPU_array == PETSC_CUSP_CPU || vec->valid_GPU_array == PETSC_CUSP_BOTH)) {
25: #elif defined(PETSC_HAVE_VECCUDA)
26: if ((vec->petscnative || vec->ops->getarray) && (vec->valid_GPU_array == PETSC_CUDA_CPU || vec->valid_GPU_array == PETSC_CUDA_BOTH)) {
27: #else
28: if (vec->petscnative || vec->ops->getarray) {
29: #endif
30: VecGetLocalSize(vec,&n);
31: VecGetArrayRead(vec,&x);
32: for (i=0; i<n; i++) {
33: if (begin) {
34: if (PetscIsInfOrNanScalar(x[i])) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_FP,"Vec entry at local location %D is not-a-number or infinite at beginning of function: Parameter number %D",i,argnum);
35: } else {
36: if (PetscIsInfOrNanScalar(x[i])) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_FP,"Vec entry at local location %D is not-a-number or infinite at end of function: Parameter number %D",i,argnum);
37: }
38: }
39: VecRestoreArrayRead(vec,&x);
40: }
41: return(0);
42: #else
43: return 0;
44: #endif
45: }
49: /*@
50: VecMaxPointwiseDivide - Computes the maximum of the componentwise division max = max_i abs(x_i/y_i).
52: Logically Collective on Vec
54: Input Parameters:
55: . x, y - the vectors
57: Output Parameter:
58: . max - the result
60: Level: advanced
62: Notes: x and y may be the same vector
63: if a particular y_i is zero, it is treated as 1 in the above formula
65: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMax(), VecPointwiseMin(), VecPointwiseMaxAbs()
66: @*/
67: PetscErrorCode VecMaxPointwiseDivide(Vec x,Vec y,PetscReal *max)
68: {
80: (*x->ops->maxpointwisedivide)(x,y,max);
81: return(0);
82: }
86: /*@
87: VecDot - Computes the vector dot product.
89: Collective on Vec
91: Input Parameters:
92: . x, y - the vectors
94: Output Parameter:
95: . val - the dot product
97: Performance Issues:
98: $ per-processor memory bandwidth
99: $ interprocessor latency
100: $ work load inbalance that causes certain processes to arrive much earlier than others
102: Notes for Users of Complex Numbers:
103: For complex vectors, VecDot() computes
104: $ val = (x,y) = y^H x,
105: where y^H denotes the conjugate transpose of y. Note that this corresponds to the usual "mathematicians" complex
106: inner product where the SECOND argument gets the complex conjugate. Since the BLASdot() complex conjugates the first
107: first argument we call the BLASdot() with the arguments reversed.
109: Use VecTDot() for the indefinite form
110: $ val = (x,y) = y^T x,
111: where y^T denotes the transpose of y.
113: Level: intermediate
115: Concepts: inner product
116: Concepts: vector^inner product
118: .seealso: VecMDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecDotRealPart()
119: @*/
120: PetscErrorCode VecDot(Vec x,Vec y,PetscScalar *val)
121: {
133: PetscLogEventBarrierBegin(VEC_DotBarrier,x,y,0,0,PetscObjectComm((PetscObject)x));
134: (*x->ops->dot)(x,y,val);
135: PetscLogEventBarrierEnd(VEC_DotBarrier,x,y,0,0,PetscObjectComm((PetscObject)x));
136: return(0);
137: }
141: /*@
142: VecDotRealPart - Computes the real part of the vector dot product.
144: Collective on Vec
146: Input Parameters:
147: . x, y - the vectors
149: Output Parameter:
150: . val - the real part of the dot product;
152: Performance Issues:
153: $ per-processor memory bandwidth
154: $ interprocessor latency
155: $ work load inbalance that causes certain processes to arrive much earlier than others
157: Notes for Users of Complex Numbers:
158: See VecDot() for more details on the definition of the dot product for complex numbers
160: For real numbers this returns the same value as VecDot()
162: For complex numbers in C^n (that is a vector of n components with a complex number for each component) this is equal to the usual real dot product on the
163: the space R^{2n} (that is a vector of 2n components with the real or imaginary part of the complex numbers for components)
165: Developer Note: This is not currently optimized to compute only the real part of the dot product.
167: Level: intermediate
169: Concepts: inner product
170: Concepts: vector^inner product
172: .seealso: VecMDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecDot(), VecDotNorm2()
173: @*/
174: PetscErrorCode VecDotRealPart(Vec x,Vec y,PetscReal *val)
175: {
177: PetscScalar fdot;
180: VecDot(x,y,&fdot);
181: *val = PetscRealPart(fdot);
182: return(0);
183: }
187: /*@
188: VecNorm - Computes the vector norm.
190: Collective on Vec
192: Input Parameters:
193: + x - the vector
194: - type - one of NORM_1, NORM_2, NORM_INFINITY. Also available
195: NORM_1_AND_2, which computes both norms and stores them
196: in a two element array.
198: Output Parameter:
199: . val - the norm
201: Notes:
202: $ NORM_1 denotes sum_i |x_i|
203: $ NORM_2 denotes sqrt(sum_i (x_i)^2)
204: $ NORM_INFINITY denotes max_i |x_i|
206: For complex numbers NORM_1 will return the traditional 1 norm of the 2 norm of the complex numbers; that is the 1
207: norm of the absolutely values of the complex entries. In PETSc 3.6 and earlier releases it returned the 1 norm of
208: the 1 norm of the complex entries (what is returned by the BLAS routine asum()). Both are valid norms but most
209: people expect the former.
211: Level: intermediate
213: Performance Issues:
214: $ per-processor memory bandwidth
215: $ interprocessor latency
216: $ work load inbalance that causes certain processes to arrive much earlier than others
218: Concepts: norm
219: Concepts: vector^norm
221: .seealso: VecDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecNormAvailable(),
222: VecNormBegin(), VecNormEnd()
224: @*/
225: PetscErrorCode VecNorm(Vec x,NormType type,PetscReal *val)
226: {
227: PetscBool flg;
234: if (((PetscObject)x)->precision != sizeof(PetscReal)) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Wrong precision of input argument");
236: /*
237: * Cached data?
238: */
239: if (type!=NORM_1_AND_2) {
240: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[type],*val,flg);
241: if (flg) return(0);
242: }
243: PetscLogEventBarrierBegin(VEC_NormBarrier,x,0,0,0,PetscObjectComm((PetscObject)x));
244: (*x->ops->norm)(x,type,val);
245: PetscLogEventBarrierEnd(VEC_NormBarrier,x,0,0,0,PetscObjectComm((PetscObject)x));
247: if (type!=NORM_1_AND_2) {
248: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[type],*val);
249: }
250: return(0);
251: }
255: /*@
256: VecNormAvailable - Returns the vector norm if it is already known.
258: Not Collective
260: Input Parameters:
261: + x - the vector
262: - type - one of NORM_1, NORM_2, NORM_INFINITY. Also available
263: NORM_1_AND_2, which computes both norms and stores them
264: in a two element array.
266: Output Parameter:
267: + available - PETSC_TRUE if the val returned is valid
268: - val - the norm
270: Notes:
271: $ NORM_1 denotes sum_i |x_i|
272: $ NORM_2 denotes sqrt(sum_i (x_i)^2)
273: $ NORM_INFINITY denotes max_i |x_i|
275: Level: intermediate
277: Performance Issues:
278: $ per-processor memory bandwidth
279: $ interprocessor latency
280: $ work load inbalance that causes certain processes to arrive much earlier than others
282: Compile Option:
283: PETSC_HAVE_SLOW_BLAS_NORM2 will cause a C (loop unrolled) version of the norm to be used, rather
284: than the BLAS. This should probably only be used when one is using the FORTRAN BLAS routines
285: (as opposed to vendor provided) because the FORTRAN BLAS NRM2() routine is very slow.
287: Concepts: norm
288: Concepts: vector^norm
290: .seealso: VecDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecNorm()
291: VecNormBegin(), VecNormEnd()
293: @*/
294: PetscErrorCode VecNormAvailable(Vec x,NormType type,PetscBool *available,PetscReal *val)
295: {
303: *available = PETSC_FALSE;
304: if (type!=NORM_1_AND_2) {
305: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[type],*val,*available);
306: }
307: return(0);
308: }
312: /*@
313: VecNormalize - Normalizes a vector by 2-norm.
315: Collective on Vec
317: Input Parameters:
318: + x - the vector
320: Output Parameter:
321: . x - the normalized vector
322: - val - the vector norm before normalization
324: Level: intermediate
326: Concepts: vector^normalizing
327: Concepts: normalizing^vector
329: @*/
330: PetscErrorCode VecNormalize(Vec x,PetscReal *val)
331: {
333: PetscReal norm;
338: PetscLogEventBegin(VEC_Normalize,x,0,0,0);
339: VecNorm(x,NORM_2,&norm);
340: if (norm == 0.0) {
341: PetscInfo(x,"Vector of zero norm can not be normalized; Returning only the zero norm\n");
342: } else if (norm != 1.0) {
343: PetscScalar tmp = 1.0/norm;
344: VecScale(x,tmp);
345: }
346: if (val) *val = norm;
347: PetscLogEventEnd(VEC_Normalize,x,0,0,0);
348: return(0);
349: }
353: /*@C
354: VecMax - Determines the maximum vector component and its location.
356: Collective on Vec
358: Input Parameter:
359: . x - the vector
361: Output Parameters:
362: + val - the maximum component
363: - p - the location of val (pass NULL if you don't want this)
365: Notes:
366: Returns the value PETSC_MIN_REAL and p = -1 if the vector is of length 0.
368: Returns the smallest index with the maximum value
369: Level: intermediate
371: Concepts: maximum^of vector
372: Concepts: vector^maximum value
374: .seealso: VecNorm(), VecMin()
375: @*/
376: PetscErrorCode VecMax(Vec x,PetscInt *p,PetscReal *val)
377: {
384: PetscLogEventBegin(VEC_Max,x,0,0,0);
385: (*x->ops->max)(x,p,val);
386: PetscLogEventEnd(VEC_Max,x,0,0,0);
387: return(0);
388: }
392: /*@
393: VecMin - Determines the minimum vector component and its location.
395: Collective on Vec
397: Input Parameters:
398: . x - the vector
400: Output Parameter:
401: + val - the minimum component
402: - p - the location of val (pass NULL if you don't want this location)
404: Level: intermediate
406: Notes:
407: Returns the value PETSC_MAX_REAL and p = -1 if the vector is of length 0.
409: This returns the smallest index with the minumum value
411: Concepts: minimum^of vector
412: Concepts: vector^minimum entry
414: .seealso: VecMax()
415: @*/
416: PetscErrorCode VecMin(Vec x,PetscInt *p,PetscReal *val)
417: {
424: PetscLogEventBegin(VEC_Min,x,0,0,0);
425: (*x->ops->min)(x,p,val);
426: PetscLogEventEnd(VEC_Min,x,0,0,0);
427: return(0);
428: }
432: /*@
433: VecTDot - Computes an indefinite vector dot product. That is, this
434: routine does NOT use the complex conjugate.
436: Collective on Vec
438: Input Parameters:
439: . x, y - the vectors
441: Output Parameter:
442: . val - the dot product
444: Notes for Users of Complex Numbers:
445: For complex vectors, VecTDot() computes the indefinite form
446: $ val = (x,y) = y^T x,
447: where y^T denotes the transpose of y.
449: Use VecDot() for the inner product
450: $ val = (x,y) = y^H x,
451: where y^H denotes the conjugate transpose of y.
453: Level: intermediate
455: Concepts: inner product^non-Hermitian
456: Concepts: vector^inner product
457: Concepts: non-Hermitian inner product
459: .seealso: VecDot(), VecMTDot()
460: @*/
461: PetscErrorCode VecTDot(Vec x,Vec y,PetscScalar *val)
462: {
474: PetscLogEventBegin(VEC_TDot,x,y,0,0);
475: (*x->ops->tdot)(x,y,val);
476: PetscLogEventEnd(VEC_TDot,x,y,0,0);
477: return(0);
478: }
482: /*@
483: VecScale - Scales a vector.
485: Not collective on Vec
487: Input Parameters:
488: + x - the vector
489: - alpha - the scalar
491: Output Parameter:
492: . x - the scaled vector
494: Note:
495: For a vector with n components, VecScale() computes
496: $ x[i] = alpha * x[i], for i=1,...,n.
498: Level: intermediate
500: Concepts: vector^scaling
501: Concepts: scaling^vector
503: @*/
504: PetscErrorCode VecScale(Vec x, PetscScalar alpha)
505: {
506: PetscReal norms[4] = {0.0,0.0,0.0, 0.0};
507: PetscBool flgs[4];
509: PetscInt i;
514: if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
515: PetscLogEventBegin(VEC_Scale,x,0,0,0);
516: if (alpha != (PetscScalar)1.0) {
517: /* get current stashed norms */
518: for (i=0; i<4; i++) {
519: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[i],norms[i],flgs[i]);
520: }
521: (*x->ops->scale)(x,alpha);
522: PetscObjectStateIncrease((PetscObject)x);
523: /* put the scaled stashed norms back into the Vec */
524: for (i=0; i<4; i++) {
525: if (flgs[i]) {
526: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[i],PetscAbsScalar(alpha)*norms[i]);
527: }
528: }
529: }
530: PetscLogEventEnd(VEC_Scale,x,0,0,0);
531: return(0);
532: }
536: /*@
537: VecSet - Sets all components of a vector to a single scalar value.
539: Logically Collective on Vec
541: Input Parameters:
542: + x - the vector
543: - alpha - the scalar
545: Output Parameter:
546: . x - the vector
548: Note:
549: For a vector of dimension n, VecSet() computes
550: $ x[i] = alpha, for i=1,...,n,
551: so that all vector entries then equal the identical
552: scalar value, alpha. Use the more general routine
553: VecSetValues() to set different vector entries.
555: You CANNOT call this after you have called VecSetValues() but before you call
556: VecAssemblyBegin/End().
558: Level: beginner
560: .seealso VecSetValues(), VecSetValuesBlocked(), VecSetRandom()
562: Concepts: vector^setting to constant
564: @*/
565: PetscErrorCode VecSet(Vec x,PetscScalar alpha)
566: {
567: PetscReal val;
573: if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"You cannot call this after you have called VecSetValues() but\n before you have called VecAssemblyBegin/End()");
575: VecLocked(x,1);
577: PetscLogEventBegin(VEC_Set,x,0,0,0);
578: (*x->ops->set)(x,alpha);
579: PetscLogEventEnd(VEC_Set,x,0,0,0);
580: PetscObjectStateIncrease((PetscObject)x);
582: /* norms can be simply set (if |alpha|*N not too large) */
583: val = PetscAbsScalar(alpha);
584: if (x->map->N == 0) {
585: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_1],0.0l);
586: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],0.0);
587: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_2],0.0);
588: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_FROBENIUS],0.0);
589: } else if (val > PETSC_MAX_REAL/x->map->N) {
590: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],val);
591: } else {
592: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_1],x->map->N * val);
593: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],val);
594: val = PetscSqrtReal((PetscReal)x->map->N) * val;
595: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_2],val);
596: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_FROBENIUS],val);
597: }
598: return(0);
599: }
604: /*@
605: VecAXPY - Computes y = alpha x + y.
607: Logically Collective on Vec
609: Input Parameters:
610: + alpha - the scalar
611: - x, y - the vectors
613: Output Parameter:
614: . y - output vector
616: Level: intermediate
618: Notes: x and y MUST be different vectors
620: Concepts: vector^BLAS
621: Concepts: BLAS
623: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY()
624: @*/
625: PetscErrorCode VecAXPY(Vec y,PetscScalar alpha,Vec x)
626: {
636: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y cannot be the same vector");
638: VecLocked(y,1);
640: VecLockPush(x);
641: PetscLogEventBegin(VEC_AXPY,x,y,0,0);
642: (*y->ops->axpy)(y,alpha,x);
643: PetscLogEventEnd(VEC_AXPY,x,y,0,0);
644: VecLockPop(x);
645: PetscObjectStateIncrease((PetscObject)y);
646: return(0);
647: }
651: /*@
652: VecAXPBY - Computes y = alpha x + beta y.
654: Logically Collective on Vec
656: Input Parameters:
657: + alpha,beta - the scalars
658: - x, y - the vectors
660: Output Parameter:
661: . y - output vector
663: Level: intermediate
665: Notes: x and y MUST be different vectors
667: Concepts: BLAS
668: Concepts: vector^BLAS
670: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPY()
671: @*/
672: PetscErrorCode VecAXPBY(Vec y,PetscScalar alpha,PetscScalar beta,Vec x)
673: {
683: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y cannot be the same vector");
687: PetscLogEventBegin(VEC_AXPY,x,y,0,0);
688: (*y->ops->axpby)(y,alpha,beta,x);
689: PetscLogEventEnd(VEC_AXPY,x,y,0,0);
690: PetscObjectStateIncrease((PetscObject)y);
691: return(0);
692: }
696: /*@
697: VecAXPBYPCZ - Computes z = alpha x + beta y + gamma z
699: Logically Collective on Vec
701: Input Parameters:
702: + alpha,beta, gamma - the scalars
703: - x, y, z - the vectors
705: Output Parameter:
706: . z - output vector
708: Level: intermediate
710: Notes: x, y and z must be different vectors
712: Developer Note: alpha = 1 or gamma = 1 or gamma = 0.0 are handled as special cases
714: Concepts: BLAS
715: Concepts: vector^BLAS
717: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPY()
718: @*/
719: PetscErrorCode VecAXPBYPCZ(Vec z,PetscScalar alpha,PetscScalar beta,PetscScalar gamma,Vec x,Vec y)
720: {
734: if (x == y || x == z) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x, y, and z must be different vectors");
735: if (y == z) SETERRQ(PetscObjectComm((PetscObject)y),PETSC_ERR_ARG_IDN,"x, y, and z must be different vectors");
740: PetscLogEventBegin(VEC_AXPBYPCZ,x,y,z,0);
741: (*y->ops->axpbypcz)(z,alpha,beta,gamma,x,y);
742: PetscLogEventEnd(VEC_AXPBYPCZ,x,y,z,0);
743: PetscObjectStateIncrease((PetscObject)z);
744: return(0);
745: }
749: /*@
750: VecAYPX - Computes y = x + alpha y.
752: Logically Collective on Vec
754: Input Parameters:
755: + alpha - the scalar
756: - x, y - the vectors
758: Output Parameter:
759: . y - output vector
761: Level: intermediate
763: Notes: x and y MUST be different vectors
765: Concepts: vector^BLAS
766: Concepts: BLAS
768: .seealso: VecAXPY(), VecWAXPY()
769: @*/
770: PetscErrorCode VecAYPX(Vec y,PetscScalar alpha,Vec x)
771: {
779: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y must be different vectors");
782: PetscLogEventBegin(VEC_AYPX,x,y,0,0);
783: (*y->ops->aypx)(y,alpha,x);
784: PetscLogEventEnd(VEC_AYPX,x,y,0,0);
785: PetscObjectStateIncrease((PetscObject)y);
786: return(0);
787: }
792: /*@
793: VecWAXPY - Computes w = alpha x + y.
795: Logically Collective on Vec
797: Input Parameters:
798: + alpha - the scalar
799: - x, y - the vectors
801: Output Parameter:
802: . w - the result
804: Level: intermediate
806: Notes: w cannot be either x or y, but x and y can be the same
808: Concepts: vector^BLAS
809: Concepts: BLAS
811: .seealso: VecAXPY(), VecAYPX(), VecAXPBY()
812: @*/
813: PetscErrorCode VecWAXPY(Vec w,PetscScalar alpha,Vec x,Vec y)
814: {
828: if (w == y) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Result vector w cannot be same as input vector y, suggest VecAXPY()");
829: if (w == x) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Result vector w cannot be same as input vector x, suggest VecAYPX()");
832: PetscLogEventBegin(VEC_WAXPY,x,y,w,0);
833: (*w->ops->waxpy)(w,alpha,x,y);
834: PetscLogEventEnd(VEC_WAXPY,x,y,w,0);
835: PetscObjectStateIncrease((PetscObject)w);
836: return(0);
837: }
842: /*@
843: VecSetValues - Inserts or adds values into certain locations of a vector.
845: Not Collective
847: Input Parameters:
848: + x - vector to insert in
849: . ni - number of elements to add
850: . ix - indices where to add
851: . y - array of values
852: - iora - either INSERT_VALUES or ADD_VALUES, where
853: ADD_VALUES adds values to any existing entries, and
854: INSERT_VALUES replaces existing entries with new values
856: Notes:
857: VecSetValues() sets x[ix[i]] = y[i], for i=0,...,ni-1.
859: Calls to VecSetValues() with the INSERT_VALUES and ADD_VALUES
860: options cannot be mixed without intervening calls to the assembly
861: routines.
863: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
864: MUST be called after all calls to VecSetValues() have been completed.
866: VecSetValues() uses 0-based indices in Fortran as well as in C.
868: If you call VecSetOption(x, VEC_IGNORE_NEGATIVE_INDICES,PETSC_TRUE),
869: negative indices may be passed in ix. These rows are
870: simply ignored. This allows easily inserting element load matrices
871: with homogeneous Dirchlet boundary conditions that you don't want represented
872: in the vector.
874: Level: beginner
876: Concepts: vector^setting values
878: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValuesLocal(),
879: VecSetValue(), VecSetValuesBlocked(), InsertMode, INSERT_VALUES, ADD_VALUES, VecGetValues()
880: @*/
881: PetscErrorCode VecSetValues(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
882: {
887: if (!ni) return(0);
891: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
892: (*x->ops->setvalues)(x,ni,ix,y,iora);
893: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
894: PetscObjectStateIncrease((PetscObject)x);
895: return(0);
896: }
900: /*@
901: VecGetValues - Gets values from certain locations of a vector. Currently
902: can only get values on the same processor
904: Not Collective
906: Input Parameters:
907: + x - vector to get values from
908: . ni - number of elements to get
909: - ix - indices where to get them from (in global 1d numbering)
911: Output Parameter:
912: . y - array of values
914: Notes:
915: The user provides the allocated array y; it is NOT allocated in this routine
917: VecGetValues() gets y[i] = x[ix[i]], for i=0,...,ni-1.
919: VecAssemblyBegin() and VecAssemblyEnd() MUST be called before calling this
921: VecGetValues() uses 0-based indices in Fortran as well as in C.
923: If you call VecSetOption(x, VEC_IGNORE_NEGATIVE_INDICES,PETSC_TRUE),
924: negative indices may be passed in ix. These rows are
925: simply ignored.
927: Level: beginner
929: Concepts: vector^getting values
931: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues()
932: @*/
933: PetscErrorCode VecGetValues(Vec x,PetscInt ni,const PetscInt ix[],PetscScalar y[])
934: {
939: if (!ni) return(0);
943: (*x->ops->getvalues)(x,ni,ix,y);
944: return(0);
945: }
949: /*@
950: VecSetValuesBlocked - Inserts or adds blocks of values into certain locations of a vector.
952: Not Collective
954: Input Parameters:
955: + x - vector to insert in
956: . ni - number of blocks to add
957: . ix - indices where to add in block count, rather than element count
958: . y - array of values
959: - iora - either INSERT_VALUES or ADD_VALUES, where
960: ADD_VALUES adds values to any existing entries, and
961: INSERT_VALUES replaces existing entries with new values
963: Notes:
964: VecSetValuesBlocked() sets x[bs*ix[i]+j] = y[bs*i+j],
965: for j=0,...,bs-1, for i=0,...,ni-1. where bs was set with VecSetBlockSize().
967: Calls to VecSetValuesBlocked() with the INSERT_VALUES and ADD_VALUES
968: options cannot be mixed without intervening calls to the assembly
969: routines.
971: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
972: MUST be called after all calls to VecSetValuesBlocked() have been completed.
974: VecSetValuesBlocked() uses 0-based indices in Fortran as well as in C.
976: Negative indices may be passed in ix, these rows are
977: simply ignored. This allows easily inserting element load matrices
978: with homogeneous Dirchlet boundary conditions that you don't want represented
979: in the vector.
981: Level: intermediate
983: Concepts: vector^setting values blocked
985: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValuesBlockedLocal(),
986: VecSetValues()
987: @*/
988: PetscErrorCode VecSetValuesBlocked(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
989: {
997: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
998: (*x->ops->setvaluesblocked)(x,ni,ix,y,iora);
999: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
1000: PetscObjectStateIncrease((PetscObject)x);
1001: return(0);
1002: }
1007: /*@
1008: VecSetValuesLocal - Inserts or adds values into certain locations of a vector,
1009: using a local ordering of the nodes.
1011: Not Collective
1013: Input Parameters:
1014: + x - vector to insert in
1015: . ni - number of elements to add
1016: . ix - indices where to add
1017: . y - array of values
1018: - iora - either INSERT_VALUES or ADD_VALUES, where
1019: ADD_VALUES adds values to any existing entries, and
1020: INSERT_VALUES replaces existing entries with new values
1022: Level: intermediate
1024: Notes:
1025: VecSetValuesLocal() sets x[ix[i]] = y[i], for i=0,...,ni-1.
1027: Calls to VecSetValues() with the INSERT_VALUES and ADD_VALUES
1028: options cannot be mixed without intervening calls to the assembly
1029: routines.
1031: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
1032: MUST be called after all calls to VecSetValuesLocal() have been completed.
1034: VecSetValuesLocal() uses 0-based indices in Fortran as well as in C.
1036: Concepts: vector^setting values with local numbering
1038: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetLocalToGlobalMapping(),
1039: VecSetValuesBlockedLocal()
1040: @*/
1041: PetscErrorCode VecSetValuesLocal(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
1042: {
1044: PetscInt lixp[128],*lix = lixp;
1048: if (!ni) return(0);
1053: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
1054: if (!x->ops->setvalueslocal) {
1055: if (!x->map->mapping) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Local to global never set with VecSetLocalToGlobalMapping()");
1056: if (ni > 128) {
1057: PetscMalloc1(ni,&lix);
1058: }
1059: ISLocalToGlobalMappingApply(x->map->mapping,ni,(PetscInt*)ix,lix);
1060: (*x->ops->setvalues)(x,ni,lix,y,iora);
1061: if (ni > 128) {
1062: PetscFree(lix);
1063: }
1064: } else {
1065: (*x->ops->setvalueslocal)(x,ni,ix,y,iora);
1066: }
1067: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
1068: PetscObjectStateIncrease((PetscObject)x);
1069: return(0);
1070: }
1074: /*@
1075: VecSetValuesBlockedLocal - Inserts or adds values into certain locations of a vector,
1076: using a local ordering of the nodes.
1078: Not Collective
1080: Input Parameters:
1081: + x - vector to insert in
1082: . ni - number of blocks to add
1083: . ix - indices where to add in block count, not element count
1084: . y - array of values
1085: - iora - either INSERT_VALUES or ADD_VALUES, where
1086: ADD_VALUES adds values to any existing entries, and
1087: INSERT_VALUES replaces existing entries with new values
1089: Level: intermediate
1091: Notes:
1092: VecSetValuesBlockedLocal() sets x[bs*ix[i]+j] = y[bs*i+j],
1093: for j=0,..bs-1, for i=0,...,ni-1, where bs has been set with VecSetBlockSize().
1095: Calls to VecSetValuesBlockedLocal() with the INSERT_VALUES and ADD_VALUES
1096: options cannot be mixed without intervening calls to the assembly
1097: routines.
1099: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
1100: MUST be called after all calls to VecSetValuesBlockedLocal() have been completed.
1102: VecSetValuesBlockedLocal() uses 0-based indices in Fortran as well as in C.
1105: Concepts: vector^setting values blocked with local numbering
1107: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetValuesBlocked(),
1108: VecSetLocalToGlobalMapping()
1109: @*/
1110: PetscErrorCode VecSetValuesBlockedLocal(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
1111: {
1113: PetscInt lixp[128],*lix = lixp;
1120: if (ni > 128) {
1121: PetscMalloc1(ni,&lix);
1122: }
1124: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
1125: ISLocalToGlobalMappingApplyBlock(x->map->mapping,ni,(PetscInt*)ix,lix);
1126: (*x->ops->setvaluesblocked)(x,ni,lix,y,iora);
1127: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
1128: if (ni > 128) {
1129: PetscFree(lix);
1130: }
1131: PetscObjectStateIncrease((PetscObject)x);
1132: return(0);
1133: }
1137: /*@
1138: VecMTDot - Computes indefinite vector multiple dot products.
1139: That is, it does NOT use the complex conjugate.
1141: Collective on Vec
1143: Input Parameters:
1144: + x - one vector
1145: . nv - number of vectors
1146: - y - array of vectors. Note that vectors are pointers
1148: Output Parameter:
1149: . val - array of the dot products
1151: Notes for Users of Complex Numbers:
1152: For complex vectors, VecMTDot() computes the indefinite form
1153: $ val = (x,y) = y^T x,
1154: where y^T denotes the transpose of y.
1156: Use VecMDot() for the inner product
1157: $ val = (x,y) = y^H x,
1158: where y^H denotes the conjugate transpose of y.
1160: Level: intermediate
1162: Concepts: inner product^multiple
1163: Concepts: vector^multiple inner products
1165: .seealso: VecMDot(), VecTDot()
1166: @*/
1167: PetscErrorCode VecMTDot(Vec x,PetscInt nv,const Vec y[],PetscScalar val[])
1168: {
1181: PetscLogEventBegin(VEC_MTDot,x,*y,0,0);
1182: (*x->ops->mtdot)(x,nv,y,val);
1183: PetscLogEventEnd(VEC_MTDot,x,*y,0,0);
1184: return(0);
1185: }
1189: /*@
1190: VecMDot - Computes vector multiple dot products.
1192: Collective on Vec
1194: Input Parameters:
1195: + x - one vector
1196: . nv - number of vectors
1197: - y - array of vectors.
1199: Output Parameter:
1200: . val - array of the dot products (does not allocate the array)
1202: Notes for Users of Complex Numbers:
1203: For complex vectors, VecMDot() computes
1204: $ val = (x,y) = y^H x,
1205: where y^H denotes the conjugate transpose of y.
1207: Use VecMTDot() for the indefinite form
1208: $ val = (x,y) = y^T x,
1209: where y^T denotes the transpose of y.
1211: Level: intermediate
1213: Concepts: inner product^multiple
1214: Concepts: vector^multiple inner products
1216: .seealso: VecMTDot(), VecDot()
1217: @*/
1218: PetscErrorCode VecMDot(Vec x,PetscInt nv,const Vec y[],PetscScalar val[])
1219: {
1224: if (!nv) return(0);
1225: if (nv < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Number of vectors (given %D) cannot be negative",nv);
1234: PetscLogEventBarrierBegin(VEC_MDotBarrier,x,*y,0,0,PetscObjectComm((PetscObject)x));
1235: (*x->ops->mdot)(x,nv,y,val);
1236: PetscLogEventBarrierEnd(VEC_MDotBarrier,x,*y,0,0,PetscObjectComm((PetscObject)x));
1237: return(0);
1238: }
1242: /*@
1243: VecMAXPY - Computes y = y + sum alpha[j] x[j]
1245: Logically Collective on Vec
1247: Input Parameters:
1248: + nv - number of scalars and x-vectors
1249: . alpha - array of scalars
1250: . y - one vector
1251: - x - array of vectors
1253: Level: intermediate
1255: Notes: y cannot be any of the x vectors
1257: Concepts: BLAS
1259: .seealso: VecAXPY(), VecWAXPY(), VecAYPX()
1260: @*/
1261: PetscErrorCode VecMAXPY(Vec y,PetscInt nv,const PetscScalar alpha[],Vec x[])
1262: {
1264: PetscInt i;
1268: if (!nv) return(0);
1269: if (nv < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Number of vectors (given %D) cannot be negative",nv);
1279: PetscLogEventBegin(VEC_MAXPY,*x,y,0,0);
1280: (*y->ops->maxpy)(y,nv,alpha,x);
1281: PetscLogEventEnd(VEC_MAXPY,*x,y,0,0);
1282: PetscObjectStateIncrease((PetscObject)y);
1283: return(0);
1284: }
1288: /*@
1289: VecGetSubVector - Gets a vector representing part of another vector
1291: Collective on IS (and Vec if nonlocal entries are needed)
1293: Input Arguments:
1294: + X - vector from which to extract a subvector
1295: - is - index set representing portion of X to extract
1297: Output Arguments:
1298: . Y - subvector corresponding to is
1300: Level: advanced
1302: Notes:
1303: The subvector Y should be returned with VecRestoreSubVector().
1305: This function may return a subvector without making a copy, therefore it is not safe to use the original vector while
1306: modifying the subvector. Other non-overlapping subvectors can still be obtained from X using this function.
1308: .seealso: MatGetSubMatrix()
1309: @*/
1310: PetscErrorCode VecGetSubVector(Vec X,IS is,Vec *Y)
1311: {
1312: PetscErrorCode ierr;
1313: Vec Z;
1319: if (X->ops->getsubvector) {
1320: (*X->ops->getsubvector)(X,is,&Z);
1321: } else { /* Default implementation currently does no caching */
1322: PetscInt gstart,gend,start;
1323: PetscBool contiguous,gcontiguous;
1324: VecGetOwnershipRange(X,&gstart,&gend);
1325: ISContiguousLocal(is,gstart,gend,&start,&contiguous);
1326: MPIU_Allreduce(&contiguous,&gcontiguous,1,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)is));
1327: if (gcontiguous) { /* We can do a no-copy implementation */
1328: PetscInt n,N,bs;
1329: PetscMPIInt size;
1330: PetscInt state;
1332: ISGetLocalSize(is,&n);
1333: VecGetBlockSize(X,&bs);
1334: if (n%bs || bs == 1) bs = -1; /* Do not decide block size if we do not have to */
1335: MPI_Comm_size(PetscObjectComm((PetscObject)X),&size);
1336: VecLockGet(X,&state);
1337: if (state) {
1338: const PetscScalar *x;
1339: VecGetArrayRead(X,&x);
1340: if (size == 1) {
1341: VecCreateSeqWithArray(PetscObjectComm((PetscObject)X),bs,n,(PetscScalar*)x+start,&Z);
1342: } else {
1343: ISGetSize(is,&N);
1344: VecCreateMPIWithArray(PetscObjectComm((PetscObject)X),bs,n,N,(PetscScalar*)x+start,&Z);
1345: }
1346: VecRestoreArrayRead(X,&x);
1347: VecLockPush(Z);
1348: } else {
1349: PetscScalar *x;
1350: VecGetArray(X,&x);
1351: if (size == 1) {
1352: VecCreateSeqWithArray(PetscObjectComm((PetscObject)X),bs,n,x+start,&Z);
1353: } else {
1354: ISGetSize(is,&N);
1355: VecCreateMPIWithArray(PetscObjectComm((PetscObject)X),bs,n,N,x+start,&Z);
1356: }
1357: VecRestoreArray(X,&x);
1358: }
1359: } else { /* Have to create a scatter and do a copy */
1360: VecScatter scatter;
1361: PetscInt n,N;
1362: ISGetLocalSize(is,&n);
1363: ISGetSize(is,&N);
1364: VecCreate(PetscObjectComm((PetscObject)is),&Z);
1365: VecSetSizes(Z,n,N);
1366: VecSetType(Z,((PetscObject)X)->type_name);
1367: VecScatterCreate(X,is,Z,NULL,&scatter);
1368: VecScatterBegin(scatter,X,Z,INSERT_VALUES,SCATTER_FORWARD);
1369: VecScatterEnd(scatter,X,Z,INSERT_VALUES,SCATTER_FORWARD);
1370: PetscObjectCompose((PetscObject)Z,"VecGetSubVector_Scatter",(PetscObject)scatter);
1371: VecScatterDestroy(&scatter);
1372: }
1373: }
1374: /* Record the state when the subvector was gotten so we know whether its values need to be put back */
1375: if (VecGetSubVectorSavedStateId < 0) {PetscObjectComposedDataRegister(&VecGetSubVectorSavedStateId);}
1376: PetscObjectComposedDataSetInt((PetscObject)Z,VecGetSubVectorSavedStateId,1);
1377: *Y = Z;
1378: return(0);
1379: }
1383: /*@
1384: VecRestoreSubVector - Restores a subvector extracted using VecGetSubVector()
1386: Collective on IS (and Vec if nonlocal entries need to be written)
1388: Input Arguments:
1389: + X - vector from which subvector was obtained
1390: . is - index set representing the subset of X
1391: - Y - subvector being restored
1393: Level: advanced
1395: .seealso: VecGetSubVector()
1396: @*/
1397: PetscErrorCode VecRestoreSubVector(Vec X,IS is,Vec *Y)
1398: {
1406: if (X->ops->restoresubvector) {
1407: (*X->ops->restoresubvector)(X,is,Y);
1408: } else {
1409: PETSC_UNUSED PetscObjectState dummystate = 0;
1410: PetscBool valid;
1411: PetscObjectComposedDataGetInt((PetscObject)*Y,VecGetSubVectorSavedStateId,dummystate,valid);
1412: if (!valid) {
1413: VecScatter scatter;
1415: PetscObjectQuery((PetscObject)*Y,"VecGetSubVector_Scatter",(PetscObject*)&scatter);
1416: if (scatter) {
1417: VecScatterBegin(scatter,*Y,X,INSERT_VALUES,SCATTER_REVERSE);
1418: VecScatterEnd(scatter,*Y,X,INSERT_VALUES,SCATTER_REVERSE);
1419: }
1420: }
1421: VecDestroy(Y);
1422: }
1423: return(0);
1424: }
1426: /*@
1427: VecGetLocalVectorRead - Maps the local portion of a vector into a
1428: vector. You must call VecRestoreLocalVectorRead() when the local
1429: vector is no longer needed.
1431: Not collective.
1433: Input parameter:
1434: . v - The vector for which the local vector is desired.
1436: Output parameter:
1437: . w - Upon exit this contains the local vector.
1439: Level: beginner
1440:
1441: Notes:
1442: This function is similar to VecGetArrayRead() which maps the local
1443: portion into a raw pointer. VecGetLocalVectorRead() is usually
1444: almost as efficient as VecGetArrayRead() but in certain circumstances
1445: VecGetLocalVectorRead() can be much more efficient than
1446: VecGetArrayRead(). This is because the construction of a contiguous
1447: array representing the vector data required by VecGetArrayRead() can
1448: be an expensive operation for certain vector types. For example, for
1449: GPU vectors VecGetArrayRead() requires that the data between device
1450: and host is synchronized.
1452: Unlike VecGetLocalVector(), this routine is not collective and
1453: preserves cached information.
1455: .seealso: VecRestoreLocalVectorRead(), VecGetLocalVector(), VecGetArrayRead(), VecGetArray()
1456: @*/
1459: PetscErrorCode VecGetLocalVectorRead(Vec v,Vec w)
1460: {
1462: PetscScalar *a;
1463: PetscInt m1,m2;
1468: VecGetLocalSize(v,&m1);
1469: VecGetLocalSize(w,&m2);
1470: if (m1 != m2) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Vectors of different local sizes.");
1471: if (v->ops->getlocalvectorread) {
1472: (*v->ops->getlocalvectorread)(v,w);
1473: } else {
1474: VecGetArrayRead(v,(const PetscScalar**)&a);
1475: VecPlaceArray(w,a);
1476: }
1477: return(0);
1478: }
1480: /*@
1481: VecRestoreLocalVectorRead - Unmaps the local portion of a vector
1482: previously mapped into a vector using VecGetLocalVectorRead().
1484: Not collective.
1486: Input parameter:
1487: . v - The local portion of this vector was previously mapped into w using VecGetLocalVectorRead().
1488: . w - The vector into which the local portion of v was mapped.
1490: Level: beginner
1492: .seealso: VecGetLocalVectorRead(), VecGetLocalVector(), VecGetArrayRead(), VecGetArray()
1493: @*/
1496: PetscErrorCode VecRestoreLocalVectorRead(Vec v,Vec w)
1497: {
1499: PetscScalar *a;
1504: if (v->ops->restorelocalvectorread) {
1505: (*v->ops->restorelocalvectorread)(v,w);
1506: } else {
1507: VecGetArrayRead(w,(const PetscScalar**)&a);
1508: VecRestoreArrayRead(v,(const PetscScalar**)&a);
1509: VecResetArray(w);
1510: }
1511: return(0);
1512: }
1514: /*@
1515: VecGetLocalVector - Maps the local portion of a vector into a
1516: vector.
1518: Collective on v, not collective on w.
1520: Input parameter:
1521: . v - The vector for which the local vector is desired.
1523: Output parameter:
1524: . w - Upon exit this contains the local vector.
1526: Level: beginner
1527:
1528: Notes:
1529: This function is similar to VecGetArray() which maps the local
1530: portion into a raw pointer. VecGetLocalVector() is usually about as
1531: efficient as VecGetArray() but in certain circumstances
1532: VecGetLocalVector() can be much more efficient than VecGetArray().
1533: This is because the construction of a contiguous array representing
1534: the vector data required by VecGetArray() can be an expensive
1535: operation for certain vector types. For example, for GPU vectors
1536: VecGetArray() requires that the data between device and host is
1537: synchronized.
1539: .seealso: VecRestoreLocalVector(), VecGetLocalVectorRead(), VecGetArrayRead(), VecGetArray()
1540: @*/
1543: PetscErrorCode VecGetLocalVector(Vec v,Vec w)
1544: {
1546: PetscScalar *a;
1547: PetscInt m1,m2;
1552: VecGetLocalSize(v,&m1);
1553: VecGetLocalSize(w,&m2);
1554: if (m1 != m2) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Vectors of different local sizes.");
1555: if (v->ops->getlocalvector) {
1556: (*v->ops->getlocalvector)(v,w);
1557: } else {
1558: VecGetArray(v,&a);
1559: VecPlaceArray(w,a);
1560: }
1561: return(0);
1562: }
1564: /*@
1565: VecRestoreLocalVector - Unmaps the local portion of a vector
1566: previously mapped into a vector using VecGetLocalVector().
1568: Logically collective.
1570: Input parameter:
1571: . v - The local portion of this vector was previously mapped into w using VecGetLocalVector().
1572: . w - The vector into which the local portion of v was mapped.
1574: Level: beginner
1576: .seealso: VecGetLocalVector(), VecGetLocalVectorRead(), VecRestoreLocalVectorRead(), LocalVectorRead(), VecGetArrayRead(), VecGetArray()
1577: @*/
1580: PetscErrorCode VecRestoreLocalVector(Vec v,Vec w)
1581: {
1583: PetscScalar *a;
1588: if (v->ops->restorelocalvector) {
1589: (*v->ops->restorelocalvector)(v,w);
1590: } else {
1591: VecGetArray(w,&a);
1592: VecRestoreArray(v,&a);
1593: VecResetArray(w);
1594: }
1595: return(0);
1596: }
1600: /*@C
1601: VecGetArray - Returns a pointer to a contiguous array that contains this
1602: processor's portion of the vector data. For the standard PETSc
1603: vectors, VecGetArray() returns a pointer to the local data array and
1604: does not use any copies. If the underlying vector data is not stored
1605: in a contiquous array this routine will copy the data to a contiquous
1606: array and return a pointer to that. You MUST call VecRestoreArray()
1607: when you no longer need access to the array.
1609: Logically Collective on Vec
1611: Input Parameter:
1612: . x - the vector
1614: Output Parameter:
1615: . a - location to put pointer to the array
1617: Fortran Note:
1618: This routine is used differently from Fortran 77
1619: $ Vec x
1620: $ PetscScalar x_array(1)
1621: $ PetscOffset i_x
1622: $ PetscErrorCode ierr
1623: $ call VecGetArray(x,x_array,i_x,ierr)
1624: $
1625: $ Access first local entry in vector with
1626: $ value = x_array(i_x + 1)
1627: $
1628: $ ...... other code
1629: $ call VecRestoreArray(x,x_array,i_x,ierr)
1630: For Fortran 90 see VecGetArrayF90()
1632: See the Fortran chapter of the users manual and
1633: petsc/src/snes/examples/tutorials/ex5f.F for details.
1635: Level: beginner
1637: Concepts: vector^accessing local values
1639: .seealso: VecRestoreArray(), VecGetArrayRead(), VecGetArrays(), VecGetArrayF90(), VecGetArrayReadF90(), VecPlaceArray(), VecGetArray2d()
1640: @*/
1641: PetscErrorCode VecGetArray(Vec x,PetscScalar **a)
1642: {
1647: VecLocked(x,1);
1648: if (x->petscnative) {
1649: #if defined(PETSC_HAVE_CUSP)
1650: if (x->valid_GPU_array == PETSC_CUSP_GPU) {
1651: VecCUSPCopyFromGPU(x);
1652: } else if (x->valid_GPU_array == PETSC_CUSP_UNALLOCATED) {
1653: VecCUSPAllocateCheckHost(x);
1654: }
1655: #elif defined(PETSC_HAVE_VIENNACL)
1656: if (x->valid_GPU_array == PETSC_VIENNACL_GPU) {
1657: VecViennaCLCopyFromGPU(x);
1658: }
1659: #elif defined(PETSC_HAVE_VECCUDA)
1660: if (x->valid_GPU_array == PETSC_CUDA_GPU) {
1661: VecCUDACopyFromGPU(x);
1662: } else if (x->valid_GPU_array == PETSC_CUDA_UNALLOCATED) {
1663: VecCUDAAllocateCheckHost(x);
1664: }
1665: #endif
1666: *a = *((PetscScalar**)x->data);
1667: } else {
1668: (*x->ops->getarray)(x,a);
1669: }
1670: return(0);
1671: }
1675: /*@C
1676: VecGetArrayRead - Get read-only pointer to contiguous array containing this processor's portion of the vector data.
1678: Not Collective
1680: Input Parameters:
1681: . x - the vector
1683: Output Parameter:
1684: . a - the array
1686: Level: beginner
1688: Notes:
1689: The array must be returned using a matching call to VecRestoreArrayRead().
1691: Unlike VecGetArray(), this routine is not collective and preserves cached information like vector norms.
1693: Standard PETSc vectors use contiguous storage so that this routine does not perform a copy. Other vector
1694: implementations may require a copy, but must such implementations should cache the contiguous representation so that
1695: only one copy is performed when this routine is called multiple times in sequence.
1697: .seealso: VecGetArray(), VecRestoreArray()
1698: @*/
1699: PetscErrorCode VecGetArrayRead(Vec x,const PetscScalar **a)
1700: {
1705: if (x->petscnative) {
1706: #if defined(PETSC_HAVE_CUSP)
1707: if (x->valid_GPU_array == PETSC_CUSP_GPU) {
1708: VecCUSPCopyFromGPU(x);
1709: }
1710: #elif defined(PETSC_HAVE_VIENNACL)
1711: if (x->valid_GPU_array == PETSC_VIENNACL_GPU) {
1712: VecViennaCLCopyFromGPU(x);
1713: }
1714: #elif defined(PETSC_HAVE_VECCUDA)
1715: if (x->valid_GPU_array == PETSC_CUDA_GPU) {
1716: VecCUDACopyFromGPU(x);
1717: }
1718: #endif
1719: *a = *((PetscScalar **)x->data);
1720: } else if (x->ops->getarrayread) {
1721: (*x->ops->getarrayread)(x,a);
1722: } else {
1723: (*x->ops->getarray)(x,(PetscScalar**)a);
1724: }
1725: return(0);
1726: }
1730: /*@C
1731: VecGetArrays - Returns a pointer to the arrays in a set of vectors
1732: that were created by a call to VecDuplicateVecs(). You MUST call
1733: VecRestoreArrays() when you no longer need access to the array.
1735: Logically Collective on Vec
1737: Input Parameter:
1738: + x - the vectors
1739: - n - the number of vectors
1741: Output Parameter:
1742: . a - location to put pointer to the array
1744: Fortran Note:
1745: This routine is not supported in Fortran.
1747: Level: intermediate
1749: .seealso: VecGetArray(), VecRestoreArrays()
1750: @*/
1751: PetscErrorCode VecGetArrays(const Vec x[],PetscInt n,PetscScalar **a[])
1752: {
1754: PetscInt i;
1755: PetscScalar **q;
1761: if (n <= 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must get at least one array n = %D",n);
1762: PetscMalloc1(n,&q);
1763: for (i=0; i<n; ++i) {
1764: VecGetArray(x[i],&q[i]);
1765: }
1766: *a = q;
1767: return(0);
1768: }
1772: /*@C
1773: VecRestoreArrays - Restores a group of vectors after VecGetArrays()
1774: has been called.
1776: Logically Collective on Vec
1778: Input Parameters:
1779: + x - the vector
1780: . n - the number of vectors
1781: - a - location of pointer to arrays obtained from VecGetArrays()
1783: Notes:
1784: For regular PETSc vectors this routine does not involve any copies. For
1785: any special vectors that do not store local vector data in a contiguous
1786: array, this routine will copy the data back into the underlying
1787: vector data structure from the arrays obtained with VecGetArrays().
1789: Fortran Note:
1790: This routine is not supported in Fortran.
1792: Level: intermediate
1794: .seealso: VecGetArrays(), VecRestoreArray()
1795: @*/
1796: PetscErrorCode VecRestoreArrays(const Vec x[],PetscInt n,PetscScalar **a[])
1797: {
1799: PetscInt i;
1800: PetscScalar **q = *a;
1807: for (i=0; i<n; ++i) {
1808: VecRestoreArray(x[i],&q[i]);
1809: }
1810: PetscFree(q);
1811: return(0);
1812: }
1816: /*@C
1817: VecRestoreArray - Restores a vector after VecGetArray() has been called.
1819: Logically Collective on Vec
1821: Input Parameters:
1822: + x - the vector
1823: - a - location of pointer to array obtained from VecGetArray()
1825: Level: beginner
1827: Notes:
1828: For regular PETSc vectors this routine does not involve any copies. For
1829: any special vectors that do not store local vector data in a contiguous
1830: array, this routine will copy the data back into the underlying
1831: vector data structure from the array obtained with VecGetArray().
1833: This routine actually zeros out the a pointer. This is to prevent accidental
1834: us of the array after it has been restored. If you pass null for a it will
1835: not zero the array pointer a.
1837: Fortran Note:
1838: This routine is used differently from Fortran 77
1839: $ Vec x
1840: $ PetscScalar x_array(1)
1841: $ PetscOffset i_x
1842: $ PetscErrorCode ierr
1843: $ call VecGetArray(x,x_array,i_x,ierr)
1844: $
1845: $ Access first local entry in vector with
1846: $ value = x_array(i_x + 1)
1847: $
1848: $ ...... other code
1849: $ call VecRestoreArray(x,x_array,i_x,ierr)
1851: See the Fortran chapter of the users manual and
1852: petsc/src/snes/examples/tutorials/ex5f.F for details.
1853: For Fortran 90 see VecRestoreArrayF90()
1855: .seealso: VecGetArray(), VecRestoreArrayRead(), VecRestoreArrays(), VecRestoreArrayF90(), VecRestoreArrayReadF90(), VecPlaceArray(), VecRestoreArray2d()
1856: @*/
1857: PetscErrorCode VecRestoreArray(Vec x,PetscScalar **a)
1858: {
1863: if (x->petscnative) {
1864: #if defined(PETSC_HAVE_CUSP)
1865: x->valid_GPU_array = PETSC_CUSP_CPU;
1866: #elif defined(PETSC_HAVE_VIENNACL)
1867: x->valid_GPU_array = PETSC_VIENNACL_CPU;
1868: #elif defined(PETSC_HAVE_VECCUDA)
1869: x->valid_GPU_array = PETSC_CUDA_CPU;
1870: #endif
1871: } else {
1872: (*x->ops->restorearray)(x,a);
1873: }
1874: if (a) *a = NULL;
1875: PetscObjectStateIncrease((PetscObject)x);
1876: return(0);
1877: }
1881: /*@C
1882: VecRestoreArrayRead - Restore array obtained with VecGetArrayRead()
1884: Not Collective
1886: Input Parameters:
1887: + vec - the vector
1888: - array - the array
1890: Level: beginner
1892: .seealso: VecGetArray(), VecRestoreArray()
1893: @*/
1894: PetscErrorCode VecRestoreArrayRead(Vec x,const PetscScalar **a)
1895: {
1900: if (x->petscnative) {
1901: #if defined(PETSC_HAVE_VIENNACL)
1902: x->valid_GPU_array = PETSC_VIENNACL_CPU;
1903: #endif
1904: } else if (x->ops->restorearrayread) {
1905: (*x->ops->restorearrayread)(x,a);
1906: } else {
1907: (*x->ops->restorearray)(x,(PetscScalar**)a);
1908: }
1909: if (a) *a = NULL;
1910: return(0);
1911: }
1915: /*@
1916: VecPlaceArray - Allows one to replace the array in a vector with an
1917: array provided by the user. This is useful to avoid copying an array
1918: into a vector.
1920: Not Collective
1922: Input Parameters:
1923: + vec - the vector
1924: - array - the array
1926: Notes:
1927: You can return to the original array with a call to VecResetArray()
1929: Level: developer
1931: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecResetArray()
1933: @*/
1934: PetscErrorCode VecPlaceArray(Vec vec,const PetscScalar array[])
1935: {
1942: if (vec->ops->placearray) {
1943: (*vec->ops->placearray)(vec,array);
1944: } else SETERRQ(PetscObjectComm((PetscObject)vec),PETSC_ERR_SUP,"Cannot place array in this type of vector");
1945: PetscObjectStateIncrease((PetscObject)vec);
1946: return(0);
1947: }
1951: /*@C
1952: VecReplaceArray - Allows one to replace the array in a vector with an
1953: array provided by the user. This is useful to avoid copying an array
1954: into a vector.
1956: Not Collective
1958: Input Parameters:
1959: + vec - the vector
1960: - array - the array
1962: Notes:
1963: This permanently replaces the array and frees the memory associated
1964: with the old array.
1966: The memory passed in MUST be obtained with PetscMalloc() and CANNOT be
1967: freed by the user. It will be freed when the vector is destroy.
1969: Not supported from Fortran
1971: Level: developer
1973: .seealso: VecGetArray(), VecRestoreArray(), VecPlaceArray(), VecResetArray()
1975: @*/
1976: PetscErrorCode VecReplaceArray(Vec vec,const PetscScalar array[])
1977: {
1983: if (vec->ops->replacearray) {
1984: (*vec->ops->replacearray)(vec,array);
1985: } else SETERRQ(PetscObjectComm((PetscObject)vec),PETSC_ERR_SUP,"Cannot replace array in this type of vector");
1986: PetscObjectStateIncrease((PetscObject)vec);
1987: return(0);
1988: }
1990: /*MC
1991: VecDuplicateVecsF90 - Creates several vectors of the same type as an existing vector
1992: and makes them accessible via a Fortran90 pointer.
1994: Synopsis:
1995: VecDuplicateVecsF90(Vec x,PetscInt n,{Vec, pointer :: y(:)},integer ierr)
1997: Collective on Vec
1999: Input Parameters:
2000: + x - a vector to mimic
2001: - n - the number of vectors to obtain
2003: Output Parameters:
2004: + y - Fortran90 pointer to the array of vectors
2005: - ierr - error code
2007: Example of Usage:
2008: .vb
2009: Vec x
2010: Vec, pointer :: y(:)
2011: ....
2012: call VecDuplicateVecsF90(x,2,y,ierr)
2013: call VecSet(y(2),alpha,ierr)
2014: call VecSet(y(2),alpha,ierr)
2015: ....
2016: call VecDestroyVecsF90(2,y,ierr)
2017: .ve
2019: Notes:
2020: Not yet supported for all F90 compilers
2022: Use VecDestroyVecsF90() to free the space.
2024: Level: beginner
2026: .seealso: VecDestroyVecsF90(), VecDuplicateVecs()
2028: M*/
2030: /*MC
2031: VecRestoreArrayF90 - Restores a vector to a usable state after a call to
2032: VecGetArrayF90().
2034: Synopsis:
2035: VecRestoreArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2037: Logically Collective on Vec
2039: Input Parameters:
2040: + x - vector
2041: - xx_v - the Fortran90 pointer to the array
2043: Output Parameter:
2044: . ierr - error code
2046: Example of Usage:
2047: .vb
2048: PetscScalar, pointer :: xx_v(:)
2049: ....
2050: call VecGetArrayF90(x,xx_v,ierr)
2051: xx_v(3) = a
2052: call VecRestoreArrayF90(x,xx_v,ierr)
2053: .ve
2055: Level: beginner
2057: .seealso: VecGetArrayF90(), VecGetArray(), VecRestoreArray(), UsingFortran, VecRestoreArrayReadF90()
2059: M*/
2061: /*MC
2062: VecDestroyVecsF90 - Frees a block of vectors obtained with VecDuplicateVecsF90().
2064: Synopsis:
2065: VecDestroyVecsF90(PetscInt n,{Vec, pointer :: x(:)},PetscErrorCode ierr)
2067: Collective on Vec
2069: Input Parameters:
2070: + n - the number of vectors previously obtained
2071: - x - pointer to array of vector pointers
2073: Output Parameter:
2074: . ierr - error code
2076: Notes:
2077: Not yet supported for all F90 compilers
2079: Level: beginner
2081: .seealso: VecDestroyVecs(), VecDuplicateVecsF90()
2083: M*/
2085: /*MC
2086: VecGetArrayF90 - Accesses a vector array from Fortran90. For default PETSc
2087: vectors, VecGetArrayF90() returns a pointer to the local data array. Otherwise,
2088: this routine is implementation dependent. You MUST call VecRestoreArrayF90()
2089: when you no longer need access to the array.
2091: Synopsis:
2092: VecGetArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2094: Logically Collective on Vec
2096: Input Parameter:
2097: . x - vector
2099: Output Parameters:
2100: + xx_v - the Fortran90 pointer to the array
2101: - ierr - error code
2103: Example of Usage:
2104: .vb
2105: PetscScalar, pointer :: xx_v(:)
2106: ....
2107: call VecGetArrayF90(x,xx_v,ierr)
2108: xx_v(3) = a
2109: call VecRestoreArrayF90(x,xx_v,ierr)
2110: .ve
2112: If you ONLY intend to read entries from the array and not change any entries you should use VecGetArrayReadF90().
2114: Level: beginner
2116: .seealso: VecRestoreArrayF90(), VecGetArray(), VecRestoreArray(), VecGetArrayReadF90(), UsingFortran
2118: M*/
2120: /*MC
2121: VecGetArrayReadF90 - Accesses a read only array from Fortran90. For default PETSc
2122: vectors, VecGetArrayF90() returns a pointer to the local data array. Otherwise,
2123: this routine is implementation dependent. You MUST call VecRestoreArrayReadF90()
2124: when you no longer need access to the array.
2126: Synopsis:
2127: VecGetArrayReadF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2129: Logically Collective on Vec
2131: Input Parameter:
2132: . x - vector
2134: Output Parameters:
2135: + xx_v - the Fortran90 pointer to the array
2136: - ierr - error code
2138: Example of Usage:
2139: .vb
2140: PetscScalar, pointer :: xx_v(:)
2141: ....
2142: call VecGetArrayReadF90(x,xx_v,ierr)
2143: a = xx_v(3)
2144: call VecRestoreArrayReadF90(x,xx_v,ierr)
2145: .ve
2147: If you intend to write entries into the array you must use VecGetArrayF90().
2149: Level: beginner
2151: .seealso: VecRestoreArrayReadF90(), VecGetArray(), VecRestoreArray(), VecGetArrayRead(), VecRestoreArrayRead(), VecGetArrayF90(), UsingFortran
2153: M*/
2155: /*MC
2156: VecRestoreArrayReadF90 - Restores a readonly vector to a usable state after a call to
2157: VecGetArrayReadF90().
2159: Synopsis:
2160: VecRestoreArrayReadF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2162: Logically Collective on Vec
2164: Input Parameters:
2165: + x - vector
2166: - xx_v - the Fortran90 pointer to the array
2168: Output Parameter:
2169: . ierr - error code
2171: Example of Usage:
2172: .vb
2173: PetscScalar, pointer :: xx_v(:)
2174: ....
2175: call VecGetArrayReadF90(x,xx_v,ierr)
2176: a = xx_v(3)
2177: call VecRestoreArrayReadF90(x,xx_v,ierr)
2178: .ve
2180: Level: beginner
2182: .seealso: VecGetArrayReadF90(), VecGetArray(), VecRestoreArray(), VecGetArrayRead(), VecRestoreArrayRead(),UsingFortran, VecRestoreArrayF90()
2184: M*/
2188: /*@C
2189: VecGetArray2d - Returns a pointer to a 2d contiguous array that contains this
2190: processor's portion of the vector data. You MUST call VecRestoreArray2d()
2191: when you no longer need access to the array.
2193: Logically Collective
2195: Input Parameter:
2196: + x - the vector
2197: . m - first dimension of two dimensional array
2198: . n - second dimension of two dimensional array
2199: . mstart - first index you will use in first coordinate direction (often 0)
2200: - nstart - first index in the second coordinate direction (often 0)
2202: Output Parameter:
2203: . a - location to put pointer to the array
2205: Level: developer
2207: Notes:
2208: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
2209: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2210: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2211: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
2213: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2215: Concepts: vector^accessing local values as 2d array
2217: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2218: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2219: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2220: @*/
2221: PetscErrorCode VecGetArray2d(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2222: {
2224: PetscInt i,N;
2225: PetscScalar *aa;
2231: VecGetLocalSize(x,&N);
2232: if (m*n != N) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Local array size %D does not match 2d array dimensions %D by %D",N,m,n);
2233: VecGetArray(x,&aa);
2235: PetscMalloc1(m,a);
2236: for (i=0; i<m; i++) (*a)[i] = aa + i*n - nstart;
2237: *a -= mstart;
2238: return(0);
2239: }
2243: /*@C
2244: VecRestoreArray2d - Restores a vector after VecGetArray2d() has been called.
2246: Logically Collective
2248: Input Parameters:
2249: + x - the vector
2250: . m - first dimension of two dimensional array
2251: . n - second dimension of the two dimensional array
2252: . mstart - first index you will use in first coordinate direction (often 0)
2253: . nstart - first index in the second coordinate direction (often 0)
2254: - a - location of pointer to array obtained from VecGetArray2d()
2256: Level: developer
2258: Notes:
2259: For regular PETSc vectors this routine does not involve any copies. For
2260: any special vectors that do not store local vector data in a contiguous
2261: array, this routine will copy the data back into the underlying
2262: vector data structure from the array obtained with VecGetArray().
2264: This routine actually zeros out the a pointer.
2266: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2267: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2268: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2269: @*/
2270: PetscErrorCode VecRestoreArray2d(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2271: {
2273: void *dummy;
2279: dummy = (void*)(*a + mstart);
2280: PetscFree(dummy);
2281: VecRestoreArray(x,NULL);
2282: return(0);
2283: }
2287: /*@C
2288: VecGetArray1d - Returns a pointer to a 1d contiguous array that contains this
2289: processor's portion of the vector data. You MUST call VecRestoreArray1d()
2290: when you no longer need access to the array.
2292: Logically Collective
2294: Input Parameter:
2295: + x - the vector
2296: . m - first dimension of two dimensional array
2297: - mstart - first index you will use in first coordinate direction (often 0)
2299: Output Parameter:
2300: . a - location to put pointer to the array
2302: Level: developer
2304: Notes:
2305: For a vector obtained from DMCreateLocalVector() mstart are likely
2306: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2307: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
2309: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2311: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2312: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2313: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2314: @*/
2315: PetscErrorCode VecGetArray1d(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2316: {
2318: PetscInt N;
2324: VecGetLocalSize(x,&N);
2325: if (m != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Local array size %D does not match 1d array dimensions %D",N,m);
2326: VecGetArray(x,a);
2327: *a -= mstart;
2328: return(0);
2329: }
2333: /*@C
2334: VecRestoreArray1d - Restores a vector after VecGetArray1d() has been called.
2336: Logically Collective
2338: Input Parameters:
2339: + x - the vector
2340: . m - first dimension of two dimensional array
2341: . mstart - first index you will use in first coordinate direction (often 0)
2342: - a - location of pointer to array obtained from VecGetArray21()
2344: Level: developer
2346: Notes:
2347: For regular PETSc vectors this routine does not involve any copies. For
2348: any special vectors that do not store local vector data in a contiguous
2349: array, this routine will copy the data back into the underlying
2350: vector data structure from the array obtained with VecGetArray1d().
2352: This routine actually zeros out the a pointer.
2354: Concepts: vector^accessing local values as 1d array
2356: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2357: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2358: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
2359: @*/
2360: PetscErrorCode VecRestoreArray1d(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2361: {
2367: VecRestoreArray(x,NULL);
2368: return(0);
2369: }
2374: /*@C
2375: VecGetArray3d - Returns a pointer to a 3d contiguous array that contains this
2376: processor's portion of the vector data. You MUST call VecRestoreArray3d()
2377: when you no longer need access to the array.
2379: Logically Collective
2381: Input Parameter:
2382: + x - the vector
2383: . m - first dimension of three dimensional array
2384: . n - second dimension of three dimensional array
2385: . p - third dimension of three dimensional array
2386: . mstart - first index you will use in first coordinate direction (often 0)
2387: . nstart - first index in the second coordinate direction (often 0)
2388: - pstart - first index in the third coordinate direction (often 0)
2390: Output Parameter:
2391: . a - location to put pointer to the array
2393: Level: developer
2395: Notes:
2396: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
2397: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2398: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2399: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
2401: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2403: Concepts: vector^accessing local values as 3d array
2405: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2406: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2407: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2408: @*/
2409: PetscErrorCode VecGetArray3d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
2410: {
2412: PetscInt i,N,j;
2413: PetscScalar *aa,**b;
2419: VecGetLocalSize(x,&N);
2420: if (m*n*p != N) SETERRQ4(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Local array size %D does not match 3d array dimensions %D by %D by %D",N,m,n,p);
2421: VecGetArray(x,&aa);
2423: PetscMalloc1(m*sizeof(PetscScalar**)+m*n,a);
2424: b = (PetscScalar**)((*a) + m);
2425: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
2426: for (i=0; i<m; i++)
2427: for (j=0; j<n; j++)
2428: b[i*n+j] = aa + i*n*p + j*p - pstart;
2430: *a -= mstart;
2431: return(0);
2432: }
2436: /*@C
2437: VecRestoreArray3d - Restores a vector after VecGetArray3d() has been called.
2439: Logically Collective
2441: Input Parameters:
2442: + x - the vector
2443: . m - first dimension of three dimensional array
2444: . n - second dimension of the three dimensional array
2445: . p - third dimension of the three dimensional array
2446: . mstart - first index you will use in first coordinate direction (often 0)
2447: . nstart - first index in the second coordinate direction (often 0)
2448: . pstart - first index in the third coordinate direction (often 0)
2449: - a - location of pointer to array obtained from VecGetArray3d()
2451: Level: developer
2453: Notes:
2454: For regular PETSc vectors this routine does not involve any copies. For
2455: any special vectors that do not store local vector data in a contiguous
2456: array, this routine will copy the data back into the underlying
2457: vector data structure from the array obtained with VecGetArray().
2459: This routine actually zeros out the a pointer.
2461: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2462: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2463: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
2464: @*/
2465: PetscErrorCode VecRestoreArray3d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
2466: {
2468: void *dummy;
2474: dummy = (void*)(*a + mstart);
2475: PetscFree(dummy);
2476: VecRestoreArray(x,NULL);
2477: return(0);
2478: }
2482: /*@C
2483: VecGetArray4d - Returns a pointer to a 4d contiguous array that contains this
2484: processor's portion of the vector data. You MUST call VecRestoreArray4d()
2485: when you no longer need access to the array.
2487: Logically Collective
2489: Input Parameter:
2490: + x - the vector
2491: . m - first dimension of four dimensional array
2492: . n - second dimension of four dimensional array
2493: . p - third dimension of four dimensional array
2494: . q - fourth dimension of four dimensional array
2495: . mstart - first index you will use in first coordinate direction (often 0)
2496: . nstart - first index in the second coordinate direction (often 0)
2497: . pstart - first index in the third coordinate direction (often 0)
2498: - qstart - first index in the fourth coordinate direction (often 0)
2500: Output Parameter:
2501: . a - location to put pointer to the array
2503: Level: beginner
2505: Notes:
2506: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
2507: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2508: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2509: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
2511: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2513: Concepts: vector^accessing local values as 3d array
2515: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2516: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2517: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2518: @*/
2519: PetscErrorCode VecGetArray4d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
2520: {
2522: PetscInt i,N,j,k;
2523: PetscScalar *aa,***b,**c;
2529: VecGetLocalSize(x,&N);
2530: if (m*n*p*q != N) SETERRQ5(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Local array size %D does not match 4d array dimensions %D by %D by %D by %D",N,m,n,p,q);
2531: VecGetArray(x,&aa);
2533: PetscMalloc1(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p,a);
2534: b = (PetscScalar***)((*a) + m);
2535: c = (PetscScalar**)(b + m*n);
2536: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
2537: for (i=0; i<m; i++)
2538: for (j=0; j<n; j++)
2539: b[i*n+j] = c + i*n*p + j*p - pstart;
2540: for (i=0; i<m; i++)
2541: for (j=0; j<n; j++)
2542: for (k=0; k<p; k++)
2543: c[i*n*p+j*p+k] = aa + i*n*p*q + j*p*q + k*q - qstart;
2544: *a -= mstart;
2545: return(0);
2546: }
2550: /*@C
2551: VecRestoreArray4d - Restores a vector after VecGetArray3d() has been called.
2553: Logically Collective
2555: Input Parameters:
2556: + x - the vector
2557: . m - first dimension of four dimensional array
2558: . n - second dimension of the four dimensional array
2559: . p - third dimension of the four dimensional array
2560: . q - fourth dimension of the four dimensional array
2561: . mstart - first index you will use in first coordinate direction (often 0)
2562: . nstart - first index in the second coordinate direction (often 0)
2563: . pstart - first index in the third coordinate direction (often 0)
2564: . qstart - first index in the fourth coordinate direction (often 0)
2565: - a - location of pointer to array obtained from VecGetArray4d()
2567: Level: beginner
2569: Notes:
2570: For regular PETSc vectors this routine does not involve any copies. For
2571: any special vectors that do not store local vector data in a contiguous
2572: array, this routine will copy the data back into the underlying
2573: vector data structure from the array obtained with VecGetArray().
2575: This routine actually zeros out the a pointer.
2577: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2578: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2579: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
2580: @*/
2581: PetscErrorCode VecRestoreArray4d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
2582: {
2584: void *dummy;
2590: dummy = (void*)(*a + mstart);
2591: PetscFree(dummy);
2592: VecRestoreArray(x,NULL);
2593: return(0);
2594: }
2598: /*@C
2599: VecGetArray2dRead - Returns a pointer to a 2d contiguous array that contains this
2600: processor's portion of the vector data. You MUST call VecRestoreArray2dRead()
2601: when you no longer need access to the array.
2603: Logically Collective
2605: Input Parameter:
2606: + x - the vector
2607: . m - first dimension of two dimensional array
2608: . n - second dimension of two dimensional array
2609: . mstart - first index you will use in first coordinate direction (often 0)
2610: - nstart - first index in the second coordinate direction (often 0)
2612: Output Parameter:
2613: . a - location to put pointer to the array
2615: Level: developer
2617: Notes:
2618: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
2619: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2620: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2621: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
2623: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2625: Concepts: vector^accessing local values as 2d array
2627: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2628: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2629: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2630: @*/
2631: PetscErrorCode VecGetArray2dRead(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2632: {
2633: PetscErrorCode ierr;
2634: PetscInt i,N;
2635: const PetscScalar *aa;
2641: VecGetLocalSize(x,&N);
2642: if (m*n != N) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Local array size %D does not match 2d array dimensions %D by %D",N,m,n);
2643: VecGetArrayRead(x,&aa);
2645: PetscMalloc1(m,a);
2646: for (i=0; i<m; i++) (*a)[i] = (PetscScalar*) aa + i*n - nstart;
2647: *a -= mstart;
2648: return(0);
2649: }
2653: /*@C
2654: VecRestoreArray2dRead - Restores a vector after VecGetArray2dRead() has been called.
2656: Logically Collective
2658: Input Parameters:
2659: + x - the vector
2660: . m - first dimension of two dimensional array
2661: . n - second dimension of the two dimensional array
2662: . mstart - first index you will use in first coordinate direction (often 0)
2663: . nstart - first index in the second coordinate direction (often 0)
2664: - a - location of pointer to array obtained from VecGetArray2d()
2666: Level: developer
2668: Notes:
2669: For regular PETSc vectors this routine does not involve any copies. For
2670: any special vectors that do not store local vector data in a contiguous
2671: array, this routine will copy the data back into the underlying
2672: vector data structure from the array obtained with VecGetArray().
2674: This routine actually zeros out the a pointer.
2676: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2677: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2678: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2679: @*/
2680: PetscErrorCode VecRestoreArray2dRead(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2681: {
2683: void *dummy;
2689: dummy = (void*)(*a + mstart);
2690: PetscFree(dummy);
2691: VecRestoreArrayRead(x,NULL);
2692: return(0);
2693: }
2697: /*@C
2698: VecGetArray1dRead - Returns a pointer to a 1d contiguous array that contains this
2699: processor's portion of the vector data. You MUST call VecRestoreArray1dRead()
2700: when you no longer need access to the array.
2702: Logically Collective
2704: Input Parameter:
2705: + x - the vector
2706: . m - first dimension of two dimensional array
2707: - mstart - first index you will use in first coordinate direction (often 0)
2709: Output Parameter:
2710: . a - location to put pointer to the array
2712: Level: developer
2714: Notes:
2715: For a vector obtained from DMCreateLocalVector() mstart are likely
2716: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2717: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
2719: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2721: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2722: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2723: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2724: @*/
2725: PetscErrorCode VecGetArray1dRead(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2726: {
2728: PetscInt N;
2734: VecGetLocalSize(x,&N);
2735: if (m != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Local array size %D does not match 1d array dimensions %D",N,m);
2736: VecGetArrayRead(x,(const PetscScalar**)a);
2737: *a -= mstart;
2738: return(0);
2739: }
2743: /*@C
2744: VecRestoreArray1dRead - Restores a vector after VecGetArray1dRead() has been called.
2746: Logically Collective
2748: Input Parameters:
2749: + x - the vector
2750: . m - first dimension of two dimensional array
2751: . mstart - first index you will use in first coordinate direction (often 0)
2752: - a - location of pointer to array obtained from VecGetArray21()
2754: Level: developer
2756: Notes:
2757: For regular PETSc vectors this routine does not involve any copies. For
2758: any special vectors that do not store local vector data in a contiguous
2759: array, this routine will copy the data back into the underlying
2760: vector data structure from the array obtained with VecGetArray1dRead().
2762: This routine actually zeros out the a pointer.
2764: Concepts: vector^accessing local values as 1d array
2766: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2767: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2768: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
2769: @*/
2770: PetscErrorCode VecRestoreArray1dRead(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2771: {
2777: VecRestoreArrayRead(x,NULL);
2778: return(0);
2779: }
2784: /*@C
2785: VecGetArray3dRead - Returns a pointer to a 3d contiguous array that contains this
2786: processor's portion of the vector data. You MUST call VecRestoreArray3dRead()
2787: when you no longer need access to the array.
2789: Logically Collective
2791: Input Parameter:
2792: + x - the vector
2793: . m - first dimension of three dimensional array
2794: . n - second dimension of three dimensional array
2795: . p - third dimension of three dimensional array
2796: . mstart - first index you will use in first coordinate direction (often 0)
2797: . nstart - first index in the second coordinate direction (often 0)
2798: - pstart - first index in the third coordinate direction (often 0)
2800: Output Parameter:
2801: . a - location to put pointer to the array
2803: Level: developer
2805: Notes:
2806: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
2807: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2808: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2809: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3dRead().
2811: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2813: Concepts: vector^accessing local values as 3d array
2815: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2816: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2817: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2818: @*/
2819: PetscErrorCode VecGetArray3dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
2820: {
2821: PetscErrorCode ierr;
2822: PetscInt i,N,j;
2823: const PetscScalar *aa;
2824: PetscScalar **b;
2830: VecGetLocalSize(x,&N);
2831: if (m*n*p != N) SETERRQ4(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Local array size %D does not match 3d array dimensions %D by %D by %D",N,m,n,p);
2832: VecGetArrayRead(x,&aa);
2834: PetscMalloc1(m*sizeof(PetscScalar**)+m*n,a);
2835: b = (PetscScalar**)((*a) + m);
2836: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
2837: for (i=0; i<m; i++)
2838: for (j=0; j<n; j++)
2839: b[i*n+j] = (PetscScalar *)aa + i*n*p + j*p - pstart;
2841: *a -= mstart;
2842: return(0);
2843: }
2847: /*@C
2848: VecRestoreArray3dRead - Restores a vector after VecGetArray3dRead() has been called.
2850: Logically Collective
2852: Input Parameters:
2853: + x - the vector
2854: . m - first dimension of three dimensional array
2855: . n - second dimension of the three dimensional array
2856: . p - third dimension of the three dimensional array
2857: . mstart - first index you will use in first coordinate direction (often 0)
2858: . nstart - first index in the second coordinate direction (often 0)
2859: . pstart - first index in the third coordinate direction (often 0)
2860: - a - location of pointer to array obtained from VecGetArray3dRead()
2862: Level: developer
2864: Notes:
2865: For regular PETSc vectors this routine does not involve any copies. For
2866: any special vectors that do not store local vector data in a contiguous
2867: array, this routine will copy the data back into the underlying
2868: vector data structure from the array obtained with VecGetArray().
2870: This routine actually zeros out the a pointer.
2872: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2873: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2874: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
2875: @*/
2876: PetscErrorCode VecRestoreArray3dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
2877: {
2879: void *dummy;
2885: dummy = (void*)(*a + mstart);
2886: PetscFree(dummy);
2887: VecRestoreArrayRead(x,NULL);
2888: return(0);
2889: }
2893: /*@C
2894: VecGetArray4dRead - Returns a pointer to a 4d contiguous array that contains this
2895: processor's portion of the vector data. You MUST call VecRestoreArray4dRead()
2896: when you no longer need access to the array.
2898: Logically Collective
2900: Input Parameter:
2901: + x - the vector
2902: . m - first dimension of four dimensional array
2903: . n - second dimension of four dimensional array
2904: . p - third dimension of four dimensional array
2905: . q - fourth dimension of four dimensional array
2906: . mstart - first index you will use in first coordinate direction (often 0)
2907: . nstart - first index in the second coordinate direction (often 0)
2908: . pstart - first index in the third coordinate direction (often 0)
2909: - qstart - first index in the fourth coordinate direction (often 0)
2911: Output Parameter:
2912: . a - location to put pointer to the array
2914: Level: beginner
2916: Notes:
2917: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
2918: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2919: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2920: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
2922: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2924: Concepts: vector^accessing local values as 3d array
2926: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2927: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2928: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2929: @*/
2930: PetscErrorCode VecGetArray4dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
2931: {
2932: PetscErrorCode ierr;
2933: PetscInt i,N,j,k;
2934: const PetscScalar *aa;
2935: PetscScalar ***b,**c;
2941: VecGetLocalSize(x,&N);
2942: if (m*n*p*q != N) SETERRQ5(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Local array size %D does not match 4d array dimensions %D by %D by %D by %D",N,m,n,p,q);
2943: VecGetArrayRead(x,&aa);
2945: PetscMalloc1(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p,a);
2946: b = (PetscScalar***)((*a) + m);
2947: c = (PetscScalar**)(b + m*n);
2948: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
2949: for (i=0; i<m; i++)
2950: for (j=0; j<n; j++)
2951: b[i*n+j] = c + i*n*p + j*p - pstart;
2952: for (i=0; i<m; i++)
2953: for (j=0; j<n; j++)
2954: for (k=0; k<p; k++)
2955: c[i*n*p+j*p+k] = (PetscScalar*) aa + i*n*p*q + j*p*q + k*q - qstart;
2956: *a -= mstart;
2957: return(0);
2958: }
2962: /*@C
2963: VecRestoreArray4dRead - Restores a vector after VecGetArray3d() has been called.
2965: Logically Collective
2967: Input Parameters:
2968: + x - the vector
2969: . m - first dimension of four dimensional array
2970: . n - second dimension of the four dimensional array
2971: . p - third dimension of the four dimensional array
2972: . q - fourth dimension of the four dimensional array
2973: . mstart - first index you will use in first coordinate direction (often 0)
2974: . nstart - first index in the second coordinate direction (often 0)
2975: . pstart - first index in the third coordinate direction (often 0)
2976: . qstart - first index in the fourth coordinate direction (often 0)
2977: - a - location of pointer to array obtained from VecGetArray4dRead()
2979: Level: beginner
2981: Notes:
2982: For regular PETSc vectors this routine does not involve any copies. For
2983: any special vectors that do not store local vector data in a contiguous
2984: array, this routine will copy the data back into the underlying
2985: vector data structure from the array obtained with VecGetArray().
2987: This routine actually zeros out the a pointer.
2989: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2990: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2991: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
2992: @*/
2993: PetscErrorCode VecRestoreArray4dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
2994: {
2996: void *dummy;
3002: dummy = (void*)(*a + mstart);
3003: PetscFree(dummy);
3004: VecRestoreArrayRead(x,NULL);
3005: return(0);
3006: }
3008: #if defined(PETSC_USE_DEBUG)
3012: /*@
3013: VecLockGet - Gets the current lock status of a vector
3015: Logically Collective on Vec
3017: Input Parameter:
3018: . x - the vector
3020: Output Parameter:
3021: . state - greater than zero indicates the vector is still locked
3023: Level: beginner
3025: Concepts: vector^accessing local values
3027: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockPush(), VecLockGet()
3028: @*/
3029: PetscErrorCode VecLockGet(Vec x,PetscInt *state)
3030: {
3033: *state = x->lock;
3034: return(0);
3035: }
3039: /*@
3040: VecLockPush - Lock a vector from writing
3042: Logically Collective on Vec
3044: Input Parameter:
3045: . x - the vector
3047: Notes: If this is set then calls to VecGetArray() or VecSetValues() or any other routines that change the vectors values will fail.
3049: Call VecLockPop() to remove the latest lock
3051: Level: beginner
3053: Concepts: vector^accessing local values
3055: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockPop(), VecLockGet()
3056: @*/
3057: PetscErrorCode VecLockPush(Vec x)
3058: {
3061: x->lock++;
3062: return(0);
3063: }
3067: /*@
3068: VecLockPop - Unlock a vector from writing
3070: Logically Collective on Vec
3072: Input Parameter:
3073: . x - the vector
3075: Level: beginner
3077: Concepts: vector^accessing local values
3079: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockPush(), VecLockGet()
3080: @*/
3081: PetscErrorCode VecLockPop(Vec x)
3082: {
3085: x->lock--;
3086: if (x->lock < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Vector has been unlocked too many times");
3087: return(0);
3088: }
3090: #endif