Actual source code: rvector.c
petsc-3.6.4 2016-04-12
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: #else
26: if (vec->petscnative || vec->ops->getarray) {
27: #endif
28: VecGetLocalSize(vec,&n);
29: VecGetArrayRead(vec,&x);
30: for (i=0; i<n; i++) {
31: if (begin) {
32: 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);
33: } else {
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 end of function: Parameter number %D",i,argnum);
35: }
36: }
37: VecRestoreArrayRead(vec,&x);
38: }
39: return(0);
40: #else
41: return 0;
42: #endif
43: }
47: /*@
48: VecMaxPointwiseDivide - Computes the maximum of the componentwise division max = max_i abs(x_i/y_i).
50: Logically Collective on Vec
52: Input Parameters:
53: . x, y - the vectors
55: Output Parameter:
56: . max - the result
58: Level: advanced
60: Notes: x and y may be the same vector
61: if a particular y_i is zero, it is treated as 1 in the above formula
63: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMax(), VecPointwiseMin(), VecPointwiseMaxAbs()
64: @*/
65: PetscErrorCode VecMaxPointwiseDivide(Vec x,Vec y,PetscReal *max)
66: {
78: (*x->ops->maxpointwisedivide)(x,y,max);
79: return(0);
80: }
84: /*@
85: VecDot - Computes the vector dot product.
87: Collective on Vec
89: Input Parameters:
90: . x, y - the vectors
92: Output Parameter:
93: . val - the dot product
95: Performance Issues:
96: $ per-processor memory bandwidth
97: $ interprocessor latency
98: $ work load inbalance that causes certain processes to arrive much earlier than others
100: Notes for Users of Complex Numbers:
101: For complex vectors, VecDot() computes
102: $ val = (x,y) = y^H x,
103: where y^H denotes the conjugate transpose of y. Note that this corresponds to the usual "mathematicians" complex
104: inner product where the SECOND argument gets the complex conjugate. Since the BLASdot() complex conjugates the first
105: first argument we call the BLASdot() with the arguments reversed.
107: Use VecTDot() for the indefinite form
108: $ val = (x,y) = y^T x,
109: where y^T denotes the transpose of y.
111: Level: intermediate
113: Concepts: inner product
114: Concepts: vector^inner product
116: .seealso: VecMDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecDotRealPart()
117: @*/
118: PetscErrorCode VecDot(Vec x,Vec y,PetscScalar *val)
119: {
131: PetscLogEventBarrierBegin(VEC_DotBarrier,x,y,0,0,PetscObjectComm((PetscObject)x));
132: (*x->ops->dot)(x,y,val);
133: PetscLogEventBarrierEnd(VEC_DotBarrier,x,y,0,0,PetscObjectComm((PetscObject)x));
134: return(0);
135: }
139: /*@
140: VecDotRealPart - Computes the real part of the vector dot product.
142: Collective on Vec
144: Input Parameters:
145: . x, y - the vectors
147: Output Parameter:
148: . val - the real part of the dot product;
150: Performance Issues:
151: $ per-processor memory bandwidth
152: $ interprocessor latency
153: $ work load inbalance that causes certain processes to arrive much earlier than others
155: Notes for Users of Complex Numbers:
156: See VecDot() for more details on the definition of the dot product for complex numbers
158: For real numbers this returns the same value as VecDot()
160: 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
161: the space R^{2n} (that is a vector of 2n components with the real or imaginary part of the complex numbers for components)
163: Developer Note: This is not currently optimized to compute only the real part of the dot product.
165: Level: intermediate
167: Concepts: inner product
168: Concepts: vector^inner product
170: .seealso: VecMDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecDot(), VecDotNorm2()
171: @*/
172: PetscErrorCode VecDotRealPart(Vec x,Vec y,PetscReal *val)
173: {
175: PetscScalar fdot;
178: VecDot(x,y,&fdot);
179: *val = PetscRealPart(fdot);
180: return(0);
181: }
185: /*@
186: VecNorm - Computes the vector norm.
188: Collective on Vec
190: Input Parameters:
191: + x - the vector
192: - type - one of NORM_1, NORM_2, NORM_INFINITY. Also available
193: NORM_1_AND_2, which computes both norms and stores them
194: in a two element array.
196: Output Parameter:
197: . val - the norm
199: Notes:
200: $ NORM_1 denotes sum_i |x_i|
201: $ NORM_2 denotes sqrt(sum_i (x_i)^2)
202: $ NORM_INFINITY denotes max_i |x_i|
204: Level: intermediate
206: Performance Issues:
207: $ per-processor memory bandwidth
208: $ interprocessor latency
209: $ work load inbalance that causes certain processes to arrive much earlier than others
211: Compile Option:
212: PETSC_HAVE_SLOW_BLAS_NORM2 will cause a C (loop unrolled) version of the norm to be used, rather
213: than the BLAS. This should probably only be used when one is using the FORTRAN BLAS routines
214: (as opposed to vendor provided) because the FORTRAN BLAS NRM2() routine is very slow.
216: Concepts: norm
217: Concepts: vector^norm
219: .seealso: VecDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecNormAvailable(),
220: VecNormBegin(), VecNormEnd()
222: @*/
223: PetscErrorCode VecNorm(Vec x,NormType type,PetscReal *val)
224: {
225: PetscBool flg;
232: if (((PetscObject)x)->precision != sizeof(PetscReal)) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Wrong precision of input argument");
234: /*
235: * Cached data?
236: */
237: if (type!=NORM_1_AND_2) {
238: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[type],*val,flg);
239: if (flg) return(0);
240: }
241: PetscLogEventBarrierBegin(VEC_NormBarrier,x,0,0,0,PetscObjectComm((PetscObject)x));
242: (*x->ops->norm)(x,type,val);
243: PetscLogEventBarrierEnd(VEC_NormBarrier,x,0,0,0,PetscObjectComm((PetscObject)x));
245: if (type!=NORM_1_AND_2) {
246: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[type],*val);
247: }
248: return(0);
249: }
253: /*@
254: VecNormAvailable - Returns the vector norm if it is already known.
256: Not Collective
258: Input Parameters:
259: + x - the vector
260: - type - one of NORM_1, NORM_2, NORM_INFINITY. Also available
261: NORM_1_AND_2, which computes both norms and stores them
262: in a two element array.
264: Output Parameter:
265: + available - PETSC_TRUE if the val returned is valid
266: - val - the norm
268: Notes:
269: $ NORM_1 denotes sum_i |x_i|
270: $ NORM_2 denotes sqrt(sum_i (x_i)^2)
271: $ NORM_INFINITY denotes max_i |x_i|
273: Level: intermediate
275: Performance Issues:
276: $ per-processor memory bandwidth
277: $ interprocessor latency
278: $ work load inbalance that causes certain processes to arrive much earlier than others
280: Compile Option:
281: PETSC_HAVE_SLOW_BLAS_NORM2 will cause a C (loop unrolled) version of the norm to be used, rather
282: than the BLAS. This should probably only be used when one is using the FORTRAN BLAS routines
283: (as opposed to vendor provided) because the FORTRAN BLAS NRM2() routine is very slow.
285: Concepts: norm
286: Concepts: vector^norm
288: .seealso: VecDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecNorm()
289: VecNormBegin(), VecNormEnd()
291: @*/
292: PetscErrorCode VecNormAvailable(Vec x,NormType type,PetscBool *available,PetscReal *val)
293: {
301: *available = PETSC_FALSE;
302: if (type!=NORM_1_AND_2) {
303: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[type],*val,*available);
304: }
305: return(0);
306: }
310: /*@
311: VecNormalize - Normalizes a vector by 2-norm.
313: Collective on Vec
315: Input Parameters:
316: + x - the vector
318: Output Parameter:
319: . x - the normalized vector
320: - val - the vector norm before normalization
322: Level: intermediate
324: Concepts: vector^normalizing
325: Concepts: normalizing^vector
327: @*/
328: PetscErrorCode VecNormalize(Vec x,PetscReal *val)
329: {
331: PetscReal norm;
336: PetscLogEventBegin(VEC_Normalize,x,0,0,0);
337: VecNorm(x,NORM_2,&norm);
338: if (norm == 0.0) {
339: PetscInfo(x,"Vector of zero norm can not be normalized; Returning only the zero norm\n");
340: } else if (norm != 1.0) {
341: PetscScalar tmp = 1.0/norm;
342: VecScale(x,tmp);
343: }
344: if (val) *val = norm;
345: PetscLogEventEnd(VEC_Normalize,x,0,0,0);
346: return(0);
347: }
351: /*@C
352: VecMax - Determines the maximum vector component and its location.
354: Collective on Vec
356: Input Parameter:
357: . x - the vector
359: Output Parameters:
360: + val - the maximum component
361: - p - the location of val (pass NULL if you don't want this)
363: Notes:
364: Returns the value PETSC_MIN_REAL and p = -1 if the vector is of length 0.
366: Returns the smallest index with the maximum value
367: Level: intermediate
369: Concepts: maximum^of vector
370: Concepts: vector^maximum value
372: .seealso: VecNorm(), VecMin()
373: @*/
374: PetscErrorCode VecMax(Vec x,PetscInt *p,PetscReal *val)
375: {
382: PetscLogEventBegin(VEC_Max,x,0,0,0);
383: (*x->ops->max)(x,p,val);
384: PetscLogEventEnd(VEC_Max,x,0,0,0);
385: return(0);
386: }
390: /*@
391: VecMin - Determines the minimum vector component and its location.
393: Collective on Vec
395: Input Parameters:
396: . x - the vector
398: Output Parameter:
399: + val - the minimum component
400: - p - the location of val (pass NULL if you don't want this location)
402: Level: intermediate
404: Notes:
405: Returns the value PETSC_MAX_REAL and p = -1 if the vector is of length 0.
407: This returns the smallest index with the minumum value
409: Concepts: minimum^of vector
410: Concepts: vector^minimum entry
412: .seealso: VecMax()
413: @*/
414: PetscErrorCode VecMin(Vec x,PetscInt *p,PetscReal *val)
415: {
422: PetscLogEventBegin(VEC_Min,x,0,0,0);
423: (*x->ops->min)(x,p,val);
424: PetscLogEventEnd(VEC_Min,x,0,0,0);
425: return(0);
426: }
430: /*@
431: VecTDot - Computes an indefinite vector dot product. That is, this
432: routine does NOT use the complex conjugate.
434: Collective on Vec
436: Input Parameters:
437: . x, y - the vectors
439: Output Parameter:
440: . val - the dot product
442: Notes for Users of Complex Numbers:
443: For complex vectors, VecTDot() computes the indefinite form
444: $ val = (x,y) = y^T x,
445: where y^T denotes the transpose of y.
447: Use VecDot() for the inner product
448: $ val = (x,y) = y^H x,
449: where y^H denotes the conjugate transpose of y.
451: Level: intermediate
453: Concepts: inner product^non-Hermitian
454: Concepts: vector^inner product
455: Concepts: non-Hermitian inner product
457: .seealso: VecDot(), VecMTDot()
458: @*/
459: PetscErrorCode VecTDot(Vec x,Vec y,PetscScalar *val)
460: {
472: PetscLogEventBegin(VEC_TDot,x,y,0,0);
473: (*x->ops->tdot)(x,y,val);
474: PetscLogEventEnd(VEC_TDot,x,y,0,0);
475: return(0);
476: }
480: /*@
481: VecScale - Scales a vector.
483: Not collective on Vec
485: Input Parameters:
486: + x - the vector
487: - alpha - the scalar
489: Output Parameter:
490: . x - the scaled vector
492: Note:
493: For a vector with n components, VecScale() computes
494: $ x[i] = alpha * x[i], for i=1,...,n.
496: Level: intermediate
498: Concepts: vector^scaling
499: Concepts: scaling^vector
501: @*/
502: PetscErrorCode VecScale(Vec x, PetscScalar alpha)
503: {
504: PetscReal norms[4] = {0.0,0.0,0.0, 0.0};
505: PetscBool flgs[4];
507: PetscInt i;
512: if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
513: PetscLogEventBegin(VEC_Scale,x,0,0,0);
514: if (alpha != (PetscScalar)1.0) {
515: /* get current stashed norms */
516: for (i=0; i<4; i++) {
517: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[i],norms[i],flgs[i]);
518: }
519: (*x->ops->scale)(x,alpha);
520: PetscObjectStateIncrease((PetscObject)x);
521: /* put the scaled stashed norms back into the Vec */
522: for (i=0; i<4; i++) {
523: if (flgs[i]) {
524: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[i],PetscAbsScalar(alpha)*norms[i]);
525: }
526: }
527: }
528: PetscLogEventEnd(VEC_Scale,x,0,0,0);
529: return(0);
530: }
534: /*@
535: VecSet - Sets all components of a vector to a single scalar value.
537: Logically Collective on Vec
539: Input Parameters:
540: + x - the vector
541: - alpha - the scalar
543: Output Parameter:
544: . x - the vector
546: Note:
547: For a vector of dimension n, VecSet() computes
548: $ x[i] = alpha, for i=1,...,n,
549: so that all vector entries then equal the identical
550: scalar value, alpha. Use the more general routine
551: VecSetValues() to set different vector entries.
553: You CANNOT call this after you have called VecSetValues() but before you call
554: VecAssemblyBegin/End().
556: Level: beginner
558: .seealso VecSetValues(), VecSetValuesBlocked(), VecSetRandom()
560: Concepts: vector^setting to constant
562: @*/
563: PetscErrorCode VecSet(Vec x,PetscScalar alpha)
564: {
565: PetscReal val;
571: 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()");
573: VecLocked(x,1);
575: PetscLogEventBegin(VEC_Set,x,0,0,0);
576: (*x->ops->set)(x,alpha);
577: PetscLogEventEnd(VEC_Set,x,0,0,0);
578: PetscObjectStateIncrease((PetscObject)x);
580: /* norms can be simply set (if |alpha|*N not too large) */
581: val = PetscAbsScalar(alpha);
582: if (x->map->N == 0) {
583: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_1],0.0l);
584: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],0.0);
585: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_2],0.0);
586: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_FROBENIUS],0.0);
587: } else if (val > PETSC_MAX_REAL/x->map->N) {
588: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],val);
589: } else {
590: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_1],x->map->N * val);
591: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],val);
592: val = PetscSqrtReal((PetscReal)x->map->N) * val;
593: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_2],val);
594: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_FROBENIUS],val);
595: }
596: return(0);
597: }
602: /*@
603: VecAXPY - Computes y = alpha x + y.
605: Logically Collective on Vec
607: Input Parameters:
608: + alpha - the scalar
609: - x, y - the vectors
611: Output Parameter:
612: . y - output vector
614: Level: intermediate
616: Notes: x and y MUST be different vectors
618: Concepts: vector^BLAS
619: Concepts: BLAS
621: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY()
622: @*/
623: PetscErrorCode VecAXPY(Vec y,PetscScalar alpha,Vec x)
624: {
634: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y cannot be the same vector");
636: VecLocked(y,1);
638: VecLockPush(x);
639: PetscLogEventBegin(VEC_AXPY,x,y,0,0);
640: (*y->ops->axpy)(y,alpha,x);
641: PetscLogEventEnd(VEC_AXPY,x,y,0,0);
642: VecLockPop(x);
643: PetscObjectStateIncrease((PetscObject)y);
644: return(0);
645: }
649: /*@
650: VecAXPBY - Computes y = alpha x + beta y.
652: Logically Collective on Vec
654: Input Parameters:
655: + alpha,beta - the scalars
656: - x, y - the vectors
658: Output Parameter:
659: . y - output vector
661: Level: intermediate
663: Notes: x and y MUST be different vectors
665: Concepts: BLAS
666: Concepts: vector^BLAS
668: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPY()
669: @*/
670: PetscErrorCode VecAXPBY(Vec y,PetscScalar alpha,PetscScalar beta,Vec x)
671: {
681: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y cannot be the same vector");
685: PetscLogEventBegin(VEC_AXPY,x,y,0,0);
686: (*y->ops->axpby)(y,alpha,beta,x);
687: PetscLogEventEnd(VEC_AXPY,x,y,0,0);
688: PetscObjectStateIncrease((PetscObject)y);
689: return(0);
690: }
694: /*@
695: VecAXPBYPCZ - Computes z = alpha x + beta y + gamma z
697: Logically Collective on Vec
699: Input Parameters:
700: + alpha,beta, gamma - the scalars
701: - x, y, z - the vectors
703: Output Parameter:
704: . z - output vector
706: Level: intermediate
708: Notes: x, y and z must be different vectors
710: Developer Note: alpha = 1 or gamma = 1 or gamma = 0.0 are handled as special cases
712: Concepts: BLAS
713: Concepts: vector^BLAS
715: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPY()
716: @*/
717: PetscErrorCode VecAXPBYPCZ(Vec z,PetscScalar alpha,PetscScalar beta,PetscScalar gamma,Vec x,Vec y)
718: {
732: if (x == y || x == z) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x, y, and z must be different vectors");
733: if (y == z) SETERRQ(PetscObjectComm((PetscObject)y),PETSC_ERR_ARG_IDN,"x, y, and z must be different vectors");
738: PetscLogEventBegin(VEC_AXPBYPCZ,x,y,z,0);
739: (*y->ops->axpbypcz)(z,alpha,beta,gamma,x,y);
740: PetscLogEventEnd(VEC_AXPBYPCZ,x,y,z,0);
741: PetscObjectStateIncrease((PetscObject)z);
742: return(0);
743: }
747: /*@
748: VecAYPX - Computes y = x + alpha y.
750: Logically Collective on Vec
752: Input Parameters:
753: + alpha - the scalar
754: - x, y - the vectors
756: Output Parameter:
757: . y - output vector
759: Level: intermediate
761: Notes: x and y MUST be different vectors
763: Concepts: vector^BLAS
764: Concepts: BLAS
766: .seealso: VecAXPY(), VecWAXPY()
767: @*/
768: PetscErrorCode VecAYPX(Vec y,PetscScalar alpha,Vec x)
769: {
777: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y must be different vectors");
780: PetscLogEventBegin(VEC_AYPX,x,y,0,0);
781: (*y->ops->aypx)(y,alpha,x);
782: PetscLogEventEnd(VEC_AYPX,x,y,0,0);
783: PetscObjectStateIncrease((PetscObject)y);
784: return(0);
785: }
790: /*@
791: VecWAXPY - Computes w = alpha x + y.
793: Logically Collective on Vec
795: Input Parameters:
796: + alpha - the scalar
797: - x, y - the vectors
799: Output Parameter:
800: . w - the result
802: Level: intermediate
804: Notes: w cannot be either x or y, but x and y can be the same
806: Concepts: vector^BLAS
807: Concepts: BLAS
809: .seealso: VecAXPY(), VecAYPX(), VecAXPBY()
810: @*/
811: PetscErrorCode VecWAXPY(Vec w,PetscScalar alpha,Vec x,Vec y)
812: {
826: if (w == y) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Result vector w cannot be same as input vector y, suggest VecAXPY()");
827: if (w == x) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Result vector w cannot be same as input vector x, suggest VecAYPX()");
830: PetscLogEventBegin(VEC_WAXPY,x,y,w,0);
831: (*w->ops->waxpy)(w,alpha,x,y);
832: PetscLogEventEnd(VEC_WAXPY,x,y,w,0);
833: PetscObjectStateIncrease((PetscObject)w);
834: return(0);
835: }
840: /*@
841: VecSetValues - Inserts or adds values into certain locations of a vector.
843: Not Collective
845: Input Parameters:
846: + x - vector to insert in
847: . ni - number of elements to add
848: . ix - indices where to add
849: . y - array of values
850: - iora - either INSERT_VALUES or ADD_VALUES, where
851: ADD_VALUES adds values to any existing entries, and
852: INSERT_VALUES replaces existing entries with new values
854: Notes:
855: VecSetValues() sets x[ix[i]] = y[i], for i=0,...,ni-1.
857: Calls to VecSetValues() with the INSERT_VALUES and ADD_VALUES
858: options cannot be mixed without intervening calls to the assembly
859: routines.
861: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
862: MUST be called after all calls to VecSetValues() have been completed.
864: VecSetValues() uses 0-based indices in Fortran as well as in C.
866: If you call VecSetOption(x, VEC_IGNORE_NEGATIVE_INDICES,PETSC_TRUE),
867: negative indices may be passed in ix. These rows are
868: simply ignored. This allows easily inserting element load matrices
869: with homogeneous Dirchlet boundary conditions that you don't want represented
870: in the vector.
872: Level: beginner
874: Concepts: vector^setting values
876: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValuesLocal(),
877: VecSetValue(), VecSetValuesBlocked(), InsertMode, INSERT_VALUES, ADD_VALUES, VecGetValues()
878: @*/
879: PetscErrorCode VecSetValues(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
880: {
885: if (!ni) return(0);
889: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
890: (*x->ops->setvalues)(x,ni,ix,y,iora);
891: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
892: PetscObjectStateIncrease((PetscObject)x);
893: return(0);
894: }
898: /*@
899: VecGetValues - Gets values from certain locations of a vector. Currently
900: can only get values on the same processor
902: Not Collective
904: Input Parameters:
905: + x - vector to get values from
906: . ni - number of elements to get
907: - ix - indices where to get them from (in global 1d numbering)
909: Output Parameter:
910: . y - array of values
912: Notes:
913: The user provides the allocated array y; it is NOT allocated in this routine
915: VecGetValues() gets y[i] = x[ix[i]], for i=0,...,ni-1.
917: VecAssemblyBegin() and VecAssemblyEnd() MUST be called before calling this
919: VecGetValues() uses 0-based indices in Fortran as well as in C.
921: If you call VecSetOption(x, VEC_IGNORE_NEGATIVE_INDICES,PETSC_TRUE),
922: negative indices may be passed in ix. These rows are
923: simply ignored.
925: Level: beginner
927: Concepts: vector^getting values
929: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecGetValuesLocal(),
930: VecGetValuesBlocked(), InsertMode, INSERT_VALUES, ADD_VALUES, VecSetValues()
931: @*/
932: PetscErrorCode VecGetValues(Vec x,PetscInt ni,const PetscInt ix[],PetscScalar y[])
933: {
938: if (!ni) return(0);
942: (*x->ops->getvalues)(x,ni,ix,y);
943: return(0);
944: }
948: /*@
949: VecSetValuesBlocked - Inserts or adds blocks of values into certain locations of a vector.
951: Not Collective
953: Input Parameters:
954: + x - vector to insert in
955: . ni - number of blocks to add
956: . ix - indices where to add in block count, rather than element count
957: . y - array of values
958: - iora - either INSERT_VALUES or ADD_VALUES, where
959: ADD_VALUES adds values to any existing entries, and
960: INSERT_VALUES replaces existing entries with new values
962: Notes:
963: VecSetValuesBlocked() sets x[bs*ix[i]+j] = y[bs*i+j],
964: for j=0,...,bs, for i=0,...,ni-1. where bs was set with VecSetBlockSize().
966: Calls to VecSetValuesBlocked() with the INSERT_VALUES and ADD_VALUES
967: options cannot be mixed without intervening calls to the assembly
968: routines.
970: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
971: MUST be called after all calls to VecSetValuesBlocked() have been completed.
973: VecSetValuesBlocked() uses 0-based indices in Fortran as well as in C.
975: Negative indices may be passed in ix, these rows are
976: simply ignored. This allows easily inserting element load matrices
977: with homogeneous Dirchlet boundary conditions that you don't want represented
978: in the vector.
980: Level: intermediate
982: Concepts: vector^setting values blocked
984: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValuesBlockedLocal(),
985: VecSetValues()
986: @*/
987: PetscErrorCode VecSetValuesBlocked(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
988: {
996: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
997: (*x->ops->setvaluesblocked)(x,ni,ix,y,iora);
998: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
999: PetscObjectStateIncrease((PetscObject)x);
1000: return(0);
1001: }
1006: /*@
1007: VecSetValuesLocal - Inserts or adds values into certain locations of a vector,
1008: using a local ordering of the nodes.
1010: Not Collective
1012: Input Parameters:
1013: + x - vector to insert in
1014: . ni - number of elements to add
1015: . ix - indices where to add
1016: . y - array of values
1017: - iora - either INSERT_VALUES or ADD_VALUES, where
1018: ADD_VALUES adds values to any existing entries, and
1019: INSERT_VALUES replaces existing entries with new values
1021: Level: intermediate
1023: Notes:
1024: VecSetValuesLocal() sets x[ix[i]] = y[i], for i=0,...,ni-1.
1026: Calls to VecSetValues() with the INSERT_VALUES and ADD_VALUES
1027: options cannot be mixed without intervening calls to the assembly
1028: routines.
1030: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
1031: MUST be called after all calls to VecSetValuesLocal() have been completed.
1033: VecSetValuesLocal() uses 0-based indices in Fortran as well as in C.
1035: Concepts: vector^setting values with local numbering
1037: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetLocalToGlobalMapping(),
1038: VecSetValuesBlockedLocal()
1039: @*/
1040: PetscErrorCode VecSetValuesLocal(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
1041: {
1043: PetscInt lixp[128],*lix = lixp;
1047: if (!ni) return(0);
1052: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
1053: if (!x->ops->setvalueslocal) {
1054: if (!x->map->mapping) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Local to global never set with VecSetLocalToGlobalMapping()");
1055: if (ni > 128) {
1056: PetscMalloc1(ni,&lix);
1057: }
1058: ISLocalToGlobalMappingApply(x->map->mapping,ni,(PetscInt*)ix,lix);
1059: (*x->ops->setvalues)(x,ni,lix,y,iora);
1060: if (ni > 128) {
1061: PetscFree(lix);
1062: }
1063: } else {
1064: (*x->ops->setvalueslocal)(x,ni,ix,y,iora);
1065: }
1066: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
1067: PetscObjectStateIncrease((PetscObject)x);
1068: return(0);
1069: }
1073: /*@
1074: VecSetValuesBlockedLocal - Inserts or adds values into certain locations of a vector,
1075: using a local ordering of the nodes.
1077: Not Collective
1079: Input Parameters:
1080: + x - vector to insert in
1081: . ni - number of blocks to add
1082: . ix - indices where to add in block count, not element count
1083: . y - array of values
1084: - iora - either INSERT_VALUES or ADD_VALUES, where
1085: ADD_VALUES adds values to any existing entries, and
1086: INSERT_VALUES replaces existing entries with new values
1088: Level: intermediate
1090: Notes:
1091: VecSetValuesBlockedLocal() sets x[bs*ix[i]+j] = y[bs*i+j],
1092: for j=0,..bs-1, for i=0,...,ni-1, where bs has been set with VecSetBlockSize().
1094: Calls to VecSetValuesBlockedLocal() with the INSERT_VALUES and ADD_VALUES
1095: options cannot be mixed without intervening calls to the assembly
1096: routines.
1098: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
1099: MUST be called after all calls to VecSetValuesBlockedLocal() have been completed.
1101: VecSetValuesBlockedLocal() uses 0-based indices in Fortran as well as in C.
1104: Concepts: vector^setting values blocked with local numbering
1106: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetValuesBlocked(),
1107: VecSetLocalToGlobalMapping()
1108: @*/
1109: PetscErrorCode VecSetValuesBlockedLocal(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
1110: {
1112: PetscInt lixp[128],*lix = lixp;
1119: if (ni > 128) {
1120: PetscMalloc1(ni,&lix);
1121: }
1123: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
1124: ISLocalToGlobalMappingApplyBlock(x->map->mapping,ni,(PetscInt*)ix,lix);
1125: (*x->ops->setvaluesblocked)(x,ni,lix,y,iora);
1126: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
1127: if (ni > 128) {
1128: PetscFree(lix);
1129: }
1130: PetscObjectStateIncrease((PetscObject)x);
1131: return(0);
1132: }
1136: /*@
1137: VecMTDot - Computes indefinite vector multiple dot products.
1138: That is, it does NOT use the complex conjugate.
1140: Collective on Vec
1142: Input Parameters:
1143: + x - one vector
1144: . nv - number of vectors
1145: - y - array of vectors. Note that vectors are pointers
1147: Output Parameter:
1148: . val - array of the dot products
1150: Notes for Users of Complex Numbers:
1151: For complex vectors, VecMTDot() computes the indefinite form
1152: $ val = (x,y) = y^T x,
1153: where y^T denotes the transpose of y.
1155: Use VecMDot() for the inner product
1156: $ val = (x,y) = y^H x,
1157: where y^H denotes the conjugate transpose of y.
1159: Level: intermediate
1161: Concepts: inner product^multiple
1162: Concepts: vector^multiple inner products
1164: .seealso: VecMDot(), VecTDot()
1165: @*/
1166: PetscErrorCode VecMTDot(Vec x,PetscInt nv,const Vec y[],PetscScalar val[])
1167: {
1180: PetscLogEventBegin(VEC_MTDot,x,*y,0,0);
1181: (*x->ops->mtdot)(x,nv,y,val);
1182: PetscLogEventEnd(VEC_MTDot,x,*y,0,0);
1183: return(0);
1184: }
1188: /*@
1189: VecMDot - Computes vector multiple dot products.
1191: Collective on Vec
1193: Input Parameters:
1194: + x - one vector
1195: . nv - number of vectors
1196: - y - array of vectors.
1198: Output Parameter:
1199: . val - array of the dot products (does not allocate the array)
1201: Notes for Users of Complex Numbers:
1202: For complex vectors, VecMDot() computes
1203: $ val = (x,y) = y^H x,
1204: where y^H denotes the conjugate transpose of y.
1206: Use VecMTDot() for the indefinite form
1207: $ val = (x,y) = y^T x,
1208: where y^T denotes the transpose of y.
1210: Level: intermediate
1212: Concepts: inner product^multiple
1213: Concepts: vector^multiple inner products
1215: .seealso: VecMTDot(), VecDot()
1216: @*/
1217: PetscErrorCode VecMDot(Vec x,PetscInt nv,const Vec y[],PetscScalar val[])
1218: {
1223: if (!nv) return(0);
1224: if (nv < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Number of vectors (given %D) cannot be negative",nv);
1233: PetscLogEventBarrierBegin(VEC_MDotBarrier,x,*y,0,0,PetscObjectComm((PetscObject)x));
1234: (*x->ops->mdot)(x,nv,y,val);
1235: PetscLogEventBarrierEnd(VEC_MDotBarrier,x,*y,0,0,PetscObjectComm((PetscObject)x));
1236: return(0);
1237: }
1241: /*@
1242: VecMAXPY - Computes y = y + sum alpha[j] x[j]
1244: Logically Collective on Vec
1246: Input Parameters:
1247: + nv - number of scalars and x-vectors
1248: . alpha - array of scalars
1249: . y - one vector
1250: - x - array of vectors
1252: Level: intermediate
1254: Notes: y cannot be any of the x vectors
1256: Concepts: BLAS
1258: .seealso: VecAXPY(), VecWAXPY(), VecAYPX()
1259: @*/
1260: PetscErrorCode VecMAXPY(Vec y,PetscInt nv,const PetscScalar alpha[],Vec x[])
1261: {
1263: PetscInt i;
1267: if (!nv) return(0);
1268: if (nv < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Number of vectors (given %D) cannot be negative",nv);
1278: PetscLogEventBegin(VEC_MAXPY,*x,y,0,0);
1279: (*y->ops->maxpy)(y,nv,alpha,x);
1280: PetscLogEventEnd(VEC_MAXPY,*x,y,0,0);
1281: PetscObjectStateIncrease((PetscObject)y);
1282: return(0);
1283: }
1287: /*@
1288: VecGetSubVector - Gets a vector representing part of another vector
1290: Collective on IS (and Vec if nonlocal entries are needed)
1292: Input Arguments:
1293: + X - vector from which to extract a subvector
1294: - is - index set representing portion of X to extract
1296: Output Arguments:
1297: . Y - subvector corresponding to is
1299: Level: advanced
1301: Notes:
1302: The subvector Y should be returned with VecRestoreSubVector().
1304: This function may return a subvector without making a copy, therefore it is not safe to use the original vector while
1305: modifying the subvector. Other non-overlapping subvectors can still be obtained from X using this function.
1307: .seealso: MatGetSubMatrix()
1308: @*/
1309: PetscErrorCode VecGetSubVector(Vec X,IS is,Vec *Y)
1310: {
1311: PetscErrorCode ierr;
1312: Vec Z;
1318: if (X->ops->getsubvector) {
1319: (*X->ops->getsubvector)(X,is,&Z);
1320: } else { /* Default implementation currently does no caching */
1321: PetscInt gstart,gend,start;
1322: PetscBool contiguous,gcontiguous;
1323: VecGetOwnershipRange(X,&gstart,&gend);
1324: ISContiguousLocal(is,gstart,gend,&start,&contiguous);
1325: MPI_Allreduce(&contiguous,&gcontiguous,1,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)is));
1326: if (gcontiguous) { /* We can do a no-copy implementation */
1327: PetscInt n,N,bs;
1328: PetscMPIInt size;
1329: PetscInt state;
1331: ISGetLocalSize(is,&n);
1332: VecGetBlockSize(X,&bs);
1333: if (n%bs || bs == 1) bs = -1; /* Do not decide block size if we do not have to */
1334: MPI_Comm_size(PetscObjectComm((PetscObject)X),&size);
1335: VecLockGet(X,&state);
1336: if (state) {
1337: const PetscScalar *x;
1338: VecGetArrayRead(X,&x);
1339: if (size == 1) {
1340: VecCreateSeqWithArray(PetscObjectComm((PetscObject)X),bs,n,(PetscScalar*)x+start,&Z);
1341: } else {
1342: ISGetSize(is,&N);
1343: VecCreateMPIWithArray(PetscObjectComm((PetscObject)X),bs,n,N,(PetscScalar*)x+start,&Z);
1344: }
1345: VecRestoreArrayRead(X,&x);
1346: VecLockPush(Z);
1347: } else {
1348: PetscScalar *x;
1349: VecGetArray(X,&x);
1350: if (size == 1) {
1351: VecCreateSeqWithArray(PetscObjectComm((PetscObject)X),bs,n,x+start,&Z);
1352: } else {
1353: ISGetSize(is,&N);
1354: VecCreateMPIWithArray(PetscObjectComm((PetscObject)X),bs,n,N,x+start,&Z);
1355: }
1356: VecRestoreArray(X,&x);
1357: }
1358: } else { /* Have to create a scatter and do a copy */
1359: VecScatter scatter;
1360: PetscInt n,N;
1361: ISGetLocalSize(is,&n);
1362: ISGetSize(is,&N);
1363: VecCreate(PetscObjectComm((PetscObject)is),&Z);
1364: VecSetSizes(Z,n,N);
1365: VecSetType(Z,((PetscObject)X)->type_name);
1366: VecScatterCreate(X,is,Z,NULL,&scatter);
1367: VecScatterBegin(scatter,X,Z,INSERT_VALUES,SCATTER_FORWARD);
1368: VecScatterEnd(scatter,X,Z,INSERT_VALUES,SCATTER_FORWARD);
1369: PetscObjectCompose((PetscObject)Z,"VecGetSubVector_Scatter",(PetscObject)scatter);
1370: VecScatterDestroy(&scatter);
1371: }
1372: }
1373: /* Record the state when the subvector was gotten so we know whether its values need to be put back */
1374: if (VecGetSubVectorSavedStateId < 0) {PetscObjectComposedDataRegister(&VecGetSubVectorSavedStateId);}
1375: PetscObjectComposedDataSetInt((PetscObject)Z,VecGetSubVectorSavedStateId,1);
1376: *Y = Z;
1377: return(0);
1378: }
1382: /*@
1383: VecRestoreSubVector - Restores a subvector extracted using VecGetSubVector()
1385: Collective on IS (and Vec if nonlocal entries need to be written)
1387: Input Arguments:
1388: + X - vector from which subvector was obtained
1389: . is - index set representing the subset of X
1390: - Y - subvector being restored
1392: Level: advanced
1394: .seealso: VecGetSubVector()
1395: @*/
1396: PetscErrorCode VecRestoreSubVector(Vec X,IS is,Vec *Y)
1397: {
1405: if (X->ops->restoresubvector) {
1406: (*X->ops->restoresubvector)(X,is,Y);
1407: } else {
1408: PETSC_UNUSED PetscObjectState dummystate = 0;
1409: PetscBool valid;
1410: PetscObjectComposedDataGetInt((PetscObject)*Y,VecGetSubVectorSavedStateId,dummystate,valid);
1411: if (!valid) {
1412: VecScatter scatter;
1414: PetscObjectQuery((PetscObject)*Y,"VecGetSubVector_Scatter",(PetscObject*)&scatter);
1415: if (scatter) {
1416: VecScatterBegin(scatter,*Y,X,INSERT_VALUES,SCATTER_REVERSE);
1417: VecScatterEnd(scatter,*Y,X,INSERT_VALUES,SCATTER_REVERSE);
1418: }
1419: }
1420: VecDestroy(Y);
1421: }
1422: return(0);
1423: }
1425: /*@C
1426: VecGetLocalVectorRead - Maps the local portion of a vector into a
1427: vector. You must call VecRestoreLocalVectorRead() when the local
1428: vector is no longer needed.
1430: Not collective.
1432: Input parameter:
1433: . v - The vector for which the local vector is desired.
1435: Output parameter:
1436: . w - Upon exit this contains the local vector.
1438: Level: beginner
1439:
1440: Notes:
1441: This function is similar to VecGetArrayRead() which maps the local
1442: portion into a raw pointer. VecGetLocalVectorRead() is usually
1443: almost as efficient as VecGetArrayRead() but in certain circumstances
1444: VecGetLocalVectorRead() can be much more efficient than
1445: VecGetArrayRead(). This is because the construction of a contiguous
1446: array representing the vector data required by VecGetArrayRead() can
1447: be an expensive operation for certain vector types. For example, for
1448: GPU vectors VecGetArrayRead() requires that the data between device
1449: and host is synchronized.
1451: Unlike VecGetLocalVector(), this routine is not collective and
1452: preserves cached information.
1454: .seealso: VecRestoreLocalVectorRead(), VecGetLocalVector(), VecGetArrayRead(), VecGetArray()
1455: @*/
1458: PetscErrorCode VecGetLocalVectorRead(Vec v,Vec w)
1459: {
1461: PetscScalar *a;
1462: PetscInt m1,m2;
1467: VecGetLocalSize(v,&m1);
1468: VecGetLocalSize(w,&m2);
1469: if (m1 != m2) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Vectors of different local sizes.");
1470: if (v->ops->getlocalvectorread) {
1471: (*v->ops->getlocalvectorread)(v,w);
1472: } else {
1473: VecGetArrayRead(v,(const PetscScalar**)&a);
1474: VecPlaceArray(w,a);
1475: }
1476: return(0);
1477: }
1479: /*@C
1480: VecRestoreLocalVectorRead - Unmaps the local portion of a vector
1481: previously mapped into a vector using VecGetLocalVectorRead().
1483: Not collective.
1485: Input parameter:
1486: . v - The local portion of this vector was previously mapped into w using VecGetLocalVectorRead().
1487: . w - The vector into which the local portion of v was mapped.
1489: Level: beginner
1491: .seealso: VecGetLocalVectorRead(), VecGetLocalVector(), VecGetArrayRead(), VecGetArray()
1492: @*/
1495: PetscErrorCode VecRestoreLocalVectorRead(Vec v,Vec w)
1496: {
1498: PetscScalar *a;
1503: if (v->ops->restorelocalvectorread) {
1504: (*v->ops->restorelocalvectorread)(v,w);
1505: } else {
1506: VecGetArrayRead(w,(const PetscScalar**)&a);
1507: VecRestoreArrayRead(v,(const PetscScalar**)&a);
1508: VecResetArray(w);
1509: }
1510: return(0);
1511: }
1513: /*@C
1514: VecGetLocalVector - Maps the local portion of a vector into a
1515: vector.
1517: Collective on v, not collective on w.
1519: Input parameter:
1520: . v - The vector for which the local vector is desired.
1522: Output parameter:
1523: . w - Upon exit this contains the local vector.
1525: Level: beginner
1526:
1527: Notes:
1528: This function is similar to VecGetArray() which maps the local
1529: portion into a raw pointer. VecGetLocalVector() is usually about as
1530: efficient as VecGetArray() but in certain circumstances
1531: VecGetLocalVector() can be much more efficient than VecGetArray().
1532: This is because the construction of a contiguous array representing
1533: the vector data required by VecGetArray() can be an expensive
1534: operation for certain vector types. For example, for GPU vectors
1535: VecGetArray() requires that the data between device and host is
1536: synchronized.
1538: .seealso: VecRestoreLocalVector(), VecGetLocalVectorRead(), VecGetArrayRead(), VecGetArray()
1539: @*/
1542: PetscErrorCode VecGetLocalVector(Vec v,Vec w)
1543: {
1545: PetscScalar *a;
1546: PetscInt m1,m2;
1551: VecGetLocalSize(v,&m1);
1552: VecGetLocalSize(w,&m2);
1553: if (m1 != m2) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Vectors of different local sizes.");
1554: if (v->ops->getlocalvector) {
1555: (*v->ops->getlocalvector)(v,w);
1556: } else {
1557: VecGetArray(v,&a);
1558: VecPlaceArray(w,a);
1559: }
1560: return(0);
1561: }
1563: /*@C
1564: VecRestoreLocalVector - Unmaps the local portion of a vector
1565: previously mapped into a vector using VecGetLocalVector().
1567: Logically collective.
1569: Input parameter:
1570: . v - The local portion of this vector was previously mapped into w using VecGetLocalVector().
1571: . w - The vector into which the local portion of v was mapped.
1573: Level: beginner
1575: .seealso: VecGetLocalVector(), VecGetLocalVectorRead(), VecRestoreLocalVectorRead(), LocalVectorRead(), VecGetArrayRead(), VecGetArray()
1576: @*/
1579: PetscErrorCode VecRestoreLocalVector(Vec v,Vec w)
1580: {
1582: PetscScalar *a;
1587: if (v->ops->restorelocalvector) {
1588: (*v->ops->restorelocalvector)(v,w);
1589: } else {
1590: VecGetArray(w,&a);
1591: VecRestoreArray(v,&a);
1592: VecResetArray(w);
1593: }
1594: return(0);
1595: }
1599: /*@C
1600: VecGetArray - Returns a pointer to a contiguous array that contains this
1601: processor's portion of the vector data. For the standard PETSc
1602: vectors, VecGetArray() returns a pointer to the local data array and
1603: does not use any copies. If the underlying vector data is not stored
1604: in a contiquous array this routine will copy the data to a contiquous
1605: array and return a pointer to that. You MUST call VecRestoreArray()
1606: when you no longer need access to the array.
1608: Logically Collective on Vec
1610: Input Parameter:
1611: . x - the vector
1613: Output Parameter:
1614: . a - location to put pointer to the array
1616: Fortran Note:
1617: This routine is used differently from Fortran 77
1618: $ Vec x
1619: $ PetscScalar x_array(1)
1620: $ PetscOffset i_x
1621: $ PetscErrorCode ierr
1622: $ call VecGetArray(x,x_array,i_x,ierr)
1623: $
1624: $ Access first local entry in vector with
1625: $ value = x_array(i_x + 1)
1626: $
1627: $ ...... other code
1628: $ call VecRestoreArray(x,x_array,i_x,ierr)
1629: For Fortran 90 see VecGetArrayF90()
1631: See the Fortran chapter of the users manual and
1632: petsc/src/snes/examples/tutorials/ex5f.F for details.
1634: Level: beginner
1636: Concepts: vector^accessing local values
1638: .seealso: VecRestoreArray(), VecGetArrayRead(), VecGetArrays(), VecGetArrayF90(), VecGetArrayReadF90(), VecPlaceArray(), VecGetArray2d()
1639: @*/
1640: PetscErrorCode VecGetArray(Vec x,PetscScalar **a)
1641: {
1646: VecLocked(x,1);
1647: if (x->petscnative) {
1648: #if defined(PETSC_HAVE_CUSP)
1649: if (x->valid_GPU_array == PETSC_CUSP_GPU) {
1650: VecCUSPCopyFromGPU(x);
1651: } else if (x->valid_GPU_array == PETSC_CUSP_UNALLOCATED) {
1652: VecCUSPAllocateCheckHost(x);
1653: }
1654: #endif
1655: #if defined(PETSC_HAVE_VIENNACL)
1656: if (x->valid_GPU_array == PETSC_VIENNACL_GPU) {
1657: VecViennaCLCopyFromGPU(x);
1658: }
1659: #endif
1660: *a = *((PetscScalar**)x->data);
1661: } else {
1662: (*x->ops->getarray)(x,a);
1663: }
1664: return(0);
1665: }
1669: /*@C
1670: VecGetArrayRead - Get read-only pointer to contiguous array containing this processor's portion of the vector data.
1672: Not Collective
1674: Input Parameters:
1675: . x - the vector
1677: Output Parameter:
1678: . a - the array
1680: Level: beginner
1682: Notes:
1683: The array must be returned using a matching call to VecRestoreArrayRead().
1685: Unlike VecGetArray(), this routine is not collective and preserves cached information like vector norms.
1687: Standard PETSc vectors use contiguous storage so that this routine does not perform a copy. Other vector
1688: implementations may require a copy, but must such implementations should cache the contiguous representation so that
1689: only one copy is performed when this routine is called multiple times in sequence.
1691: .seealso: VecGetArray(), VecRestoreArray()
1692: @*/
1693: PetscErrorCode VecGetArrayRead(Vec x,const PetscScalar **a)
1694: {
1699: if (x->petscnative) {
1700: #if defined(PETSC_HAVE_CUSP)
1701: if (x->valid_GPU_array == PETSC_CUSP_GPU) {
1702: VecCUSPCopyFromGPU(x);
1703: }
1704: #endif
1705: #if defined(PETSC_HAVE_VIENNACL)
1706: if (x->valid_GPU_array == PETSC_VIENNACL_GPU) {
1707: VecViennaCLCopyFromGPU(x);
1708: }
1709: #endif
1710: *a = *((PetscScalar **)x->data);
1711: } else if (x->ops->getarrayread) {
1712: (*x->ops->getarrayread)(x,a);
1713: } else {
1714: (*x->ops->getarray)(x,(PetscScalar**)a);
1715: }
1716: return(0);
1717: }
1721: /*@C
1722: VecGetArrays - Returns a pointer to the arrays in a set of vectors
1723: that were created by a call to VecDuplicateVecs(). You MUST call
1724: VecRestoreArrays() when you no longer need access to the array.
1726: Logically Collective on Vec
1728: Input Parameter:
1729: + x - the vectors
1730: - n - the number of vectors
1732: Output Parameter:
1733: . a - location to put pointer to the array
1735: Fortran Note:
1736: This routine is not supported in Fortran.
1738: Level: intermediate
1740: .seealso: VecGetArray(), VecRestoreArrays()
1741: @*/
1742: PetscErrorCode VecGetArrays(const Vec x[],PetscInt n,PetscScalar **a[])
1743: {
1745: PetscInt i;
1746: PetscScalar **q;
1752: if (n <= 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must get at least one array n = %D",n);
1753: PetscMalloc1(n,&q);
1754: for (i=0; i<n; ++i) {
1755: VecGetArray(x[i],&q[i]);
1756: }
1757: *a = q;
1758: return(0);
1759: }
1763: /*@C
1764: VecRestoreArrays - Restores a group of vectors after VecGetArrays()
1765: has been called.
1767: Logically Collective on Vec
1769: Input Parameters:
1770: + x - the vector
1771: . n - the number of vectors
1772: - a - location of pointer to arrays obtained from VecGetArrays()
1774: Notes:
1775: For regular PETSc vectors this routine does not involve any copies. For
1776: any special vectors that do not store local vector data in a contiguous
1777: array, this routine will copy the data back into the underlying
1778: vector data structure from the arrays obtained with VecGetArrays().
1780: Fortran Note:
1781: This routine is not supported in Fortran.
1783: Level: intermediate
1785: .seealso: VecGetArrays(), VecRestoreArray()
1786: @*/
1787: PetscErrorCode VecRestoreArrays(const Vec x[],PetscInt n,PetscScalar **a[])
1788: {
1790: PetscInt i;
1791: PetscScalar **q = *a;
1798: for (i=0; i<n; ++i) {
1799: VecRestoreArray(x[i],&q[i]);
1800: }
1801: PetscFree(q);
1802: return(0);
1803: }
1807: /*@C
1808: VecRestoreArray - Restores a vector after VecGetArray() has been called.
1810: Logically Collective on Vec
1812: Input Parameters:
1813: + x - the vector
1814: - a - location of pointer to array obtained from VecGetArray()
1816: Level: beginner
1818: Notes:
1819: For regular PETSc vectors this routine does not involve any copies. For
1820: any special vectors that do not store local vector data in a contiguous
1821: array, this routine will copy the data back into the underlying
1822: vector data structure from the array obtained with VecGetArray().
1824: This routine actually zeros out the a pointer. This is to prevent accidental
1825: us of the array after it has been restored. If you pass null for a it will
1826: not zero the array pointer a.
1828: Fortran Note:
1829: This routine is used differently from Fortran 77
1830: $ Vec x
1831: $ PetscScalar x_array(1)
1832: $ PetscOffset i_x
1833: $ PetscErrorCode ierr
1834: $ call VecGetArray(x,x_array,i_x,ierr)
1835: $
1836: $ Access first local entry in vector with
1837: $ value = x_array(i_x + 1)
1838: $
1839: $ ...... other code
1840: $ call VecRestoreArray(x,x_array,i_x,ierr)
1842: See the Fortran chapter of the users manual and
1843: petsc/src/snes/examples/tutorials/ex5f.F for details.
1844: For Fortran 90 see VecRestoreArrayF90()
1846: .seealso: VecGetArray(), VecRestoreArrayRead(), VecRestoreArrays(), VecRestoreArrayF90(), VecRestoreArrayReadF90(), VecPlaceArray(), VecRestoreArray2d()
1847: @*/
1848: PetscErrorCode VecRestoreArray(Vec x,PetscScalar **a)
1849: {
1854: if (x->petscnative) {
1855: #if defined(PETSC_HAVE_CUSP)
1856: x->valid_GPU_array = PETSC_CUSP_CPU;
1857: #endif
1858: #if defined(PETSC_HAVE_VIENNACL)
1859: x->valid_GPU_array = PETSC_VIENNACL_CPU;
1860: #endif
1861: } else {
1862: (*x->ops->restorearray)(x,a);
1863: }
1864: if (a) *a = NULL;
1865: PetscObjectStateIncrease((PetscObject)x);
1866: return(0);
1867: }
1871: /*@C
1872: VecRestoreArrayRead - Restore array obtained with VecGetArrayRead()
1874: Not Collective
1876: Input Parameters:
1877: + vec - the vector
1878: - array - the array
1880: Level: beginner
1882: .seealso: VecGetArray(), VecRestoreArray()
1883: @*/
1884: PetscErrorCode VecRestoreArrayRead(Vec x,const PetscScalar **a)
1885: {
1890: if (x->petscnative) {
1891: #if defined(PETSC_HAVE_VIENNACL)
1892: x->valid_GPU_array = PETSC_VIENNACL_CPU;
1893: #endif
1894: } else if (x->ops->restorearrayread) {
1895: (*x->ops->restorearrayread)(x,a);
1896: } else {
1897: (*x->ops->restorearray)(x,(PetscScalar**)a);
1898: }
1899: if (a) *a = NULL;
1900: return(0);
1901: }
1905: /*@
1906: VecPlaceArray - Allows one to replace the array in a vector with an
1907: array provided by the user. This is useful to avoid copying an array
1908: into a vector.
1910: Not Collective
1912: Input Parameters:
1913: + vec - the vector
1914: - array - the array
1916: Notes:
1917: You can return to the original array with a call to VecResetArray()
1919: Level: developer
1921: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecResetArray()
1923: @*/
1924: PetscErrorCode VecPlaceArray(Vec vec,const PetscScalar array[])
1925: {
1932: if (vec->ops->placearray) {
1933: (*vec->ops->placearray)(vec,array);
1934: } else SETERRQ(PetscObjectComm((PetscObject)vec),PETSC_ERR_SUP,"Cannot place array in this type of vector");
1935: PetscObjectStateIncrease((PetscObject)vec);
1936: return(0);
1937: }
1941: /*@C
1942: VecReplaceArray - Allows one to replace the array in a vector with an
1943: array provided by the user. This is useful to avoid copying an array
1944: into a vector.
1946: Not Collective
1948: Input Parameters:
1949: + vec - the vector
1950: - array - the array
1952: Notes:
1953: This permanently replaces the array and frees the memory associated
1954: with the old array.
1956: The memory passed in MUST be obtained with PetscMalloc() and CANNOT be
1957: freed by the user. It will be freed when the vector is destroy.
1959: Not supported from Fortran
1961: Level: developer
1963: .seealso: VecGetArray(), VecRestoreArray(), VecPlaceArray(), VecResetArray()
1965: @*/
1966: PetscErrorCode VecReplaceArray(Vec vec,const PetscScalar array[])
1967: {
1973: if (vec->ops->replacearray) {
1974: (*vec->ops->replacearray)(vec,array);
1975: } else SETERRQ(PetscObjectComm((PetscObject)vec),PETSC_ERR_SUP,"Cannot replace array in this type of vector");
1976: PetscObjectStateIncrease((PetscObject)vec);
1977: return(0);
1978: }
1980: /*MC
1981: VecDuplicateVecsF90 - Creates several vectors of the same type as an existing vector
1982: and makes them accessible via a Fortran90 pointer.
1984: Synopsis:
1985: VecDuplicateVecsF90(Vec x,PetscInt n,{Vec, pointer :: y(:)},integer ierr)
1987: Collective on Vec
1989: Input Parameters:
1990: + x - a vector to mimic
1991: - n - the number of vectors to obtain
1993: Output Parameters:
1994: + y - Fortran90 pointer to the array of vectors
1995: - ierr - error code
1997: Example of Usage:
1998: .vb
1999: Vec x
2000: Vec, pointer :: y(:)
2001: ....
2002: call VecDuplicateVecsF90(x,2,y,ierr)
2003: call VecSet(y(2),alpha,ierr)
2004: call VecSet(y(2),alpha,ierr)
2005: ....
2006: call VecDestroyVecsF90(2,y,ierr)
2007: .ve
2009: Notes:
2010: Not yet supported for all F90 compilers
2012: Use VecDestroyVecsF90() to free the space.
2014: Level: beginner
2016: .seealso: VecDestroyVecsF90(), VecDuplicateVecs()
2018: M*/
2020: /*MC
2021: VecRestoreArrayF90 - Restores a vector to a usable state after a call to
2022: VecGetArrayF90().
2024: Synopsis:
2025: VecRestoreArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2027: Logically Collective on Vec
2029: Input Parameters:
2030: + x - vector
2031: - xx_v - the Fortran90 pointer to the array
2033: Output Parameter:
2034: . ierr - error code
2036: Example of Usage:
2037: .vb
2038: PetscScalar, pointer :: xx_v(:)
2039: ....
2040: call VecGetArrayF90(x,xx_v,ierr)
2041: xx_v(3) = a
2042: call VecRestoreArrayF90(x,xx_v,ierr)
2043: .ve
2045: Level: beginner
2047: .seealso: VecGetArrayF90(), VecGetArray(), VecRestoreArray(), UsingFortran, VecRestoreArrayReadF90()
2049: M*/
2051: /*MC
2052: VecDestroyVecsF90 - Frees a block of vectors obtained with VecDuplicateVecsF90().
2054: Synopsis:
2055: VecDestroyVecsF90(PetscInt n,{Vec, pointer :: x(:)},PetscErrorCode ierr)
2057: Collective on Vec
2059: Input Parameters:
2060: + n - the number of vectors previously obtained
2061: - x - pointer to array of vector pointers
2063: Output Parameter:
2064: . ierr - error code
2066: Notes:
2067: Not yet supported for all F90 compilers
2069: Level: beginner
2071: .seealso: VecDestroyVecs(), VecDuplicateVecsF90()
2073: M*/
2075: /*MC
2076: VecGetArrayF90 - Accesses a vector array from Fortran90. For default PETSc
2077: vectors, VecGetArrayF90() returns a pointer to the local data array. Otherwise,
2078: this routine is implementation dependent. You MUST call VecRestoreArrayF90()
2079: when you no longer need access to the array.
2081: Synopsis:
2082: VecGetArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2084: Logically Collective on Vec
2086: Input Parameter:
2087: . x - vector
2089: Output Parameters:
2090: + xx_v - the Fortran90 pointer to the array
2091: - ierr - error code
2093: Example of Usage:
2094: .vb
2095: PetscScalar, pointer :: xx_v(:)
2096: ....
2097: call VecGetArrayF90(x,xx_v,ierr)
2098: xx_v(3) = a
2099: call VecRestoreArrayF90(x,xx_v,ierr)
2100: .ve
2102: If you ONLY intend to read entries from the array and not change any entries you should use VecGetArrayReadF90().
2104: Level: beginner
2106: .seealso: VecRestoreArrayF90(), VecGetArray(), VecRestoreArray(), VecGetArrayReadF90(), UsingFortran
2108: M*/
2110: /*MC
2111: VecGetArrayReadF90 - Accesses a read only array from Fortran90. For default PETSc
2112: vectors, VecGetArrayF90() returns a pointer to the local data array. Otherwise,
2113: this routine is implementation dependent. You MUST call VecRestoreArrayReadF90()
2114: when you no longer need access to the array.
2116: Synopsis:
2117: VecGetArrayReadF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2119: Logically Collective on Vec
2121: Input Parameter:
2122: . x - vector
2124: Output Parameters:
2125: + xx_v - the Fortran90 pointer to the array
2126: - ierr - error code
2128: Example of Usage:
2129: .vb
2130: PetscScalar, pointer :: xx_v(:)
2131: ....
2132: call VecGetArrayReadF90(x,xx_v,ierr)
2133: a = xx_v(3)
2134: call VecRestoreArrayReadF90(x,xx_v,ierr)
2135: .ve
2137: If you intend to write entries into the array you must use VecGetArrayF90().
2139: Level: beginner
2141: .seealso: VecRestoreArrayReadF90(), VecGetArray(), VecRestoreArray(), VecGetArrayRead(), VecRestoreArrayRead(), VecGetArrayF90(), UsingFortran
2143: M*/
2145: /*MC
2146: VecRestoreArrayReadF90 - Restores a readonly vector to a usable state after a call to
2147: VecGetArrayReadF90().
2149: Synopsis:
2150: VecRestoreArrayReadF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2152: Logically Collective on Vec
2154: Input Parameters:
2155: + x - vector
2156: - xx_v - the Fortran90 pointer to the array
2158: Output Parameter:
2159: . ierr - error code
2161: Example of Usage:
2162: .vb
2163: PetscScalar, pointer :: xx_v(:)
2164: ....
2165: call VecGetArrayReadF90(x,xx_v,ierr)
2166: a = xx_v(3)
2167: call VecRestoreArrayReadF90(x,xx_v,ierr)
2168: .ve
2170: Level: beginner
2172: .seealso: VecGetArrayReadF90(), VecGetArray(), VecRestoreArray(), VecGetArrayRead(), VecRestoreArrayRead(),UsingFortran, VecRestoreArrayF90()
2174: M*/
2178: /*@C
2179: VecGetArray2d - Returns a pointer to a 2d contiguous array that contains this
2180: processor's portion of the vector data. You MUST call VecRestoreArray2d()
2181: when you no longer need access to the array.
2183: Logically Collective
2185: Input Parameter:
2186: + x - the vector
2187: . m - first dimension of two dimensional array
2188: . n - second dimension of two dimensional array
2189: . mstart - first index you will use in first coordinate direction (often 0)
2190: - nstart - first index in the second coordinate direction (often 0)
2192: Output Parameter:
2193: . a - location to put pointer to the array
2195: Level: developer
2197: Notes:
2198: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
2199: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2200: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2201: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
2203: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2205: Concepts: vector^accessing local values as 2d array
2207: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2208: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2209: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2210: @*/
2211: PetscErrorCode VecGetArray2d(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2212: {
2214: PetscInt i,N;
2215: PetscScalar *aa;
2221: VecGetLocalSize(x,&N);
2222: 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);
2223: VecGetArray(x,&aa);
2225: PetscMalloc1(m,a);
2226: for (i=0; i<m; i++) (*a)[i] = aa + i*n - nstart;
2227: *a -= mstart;
2228: return(0);
2229: }
2233: /*@C
2234: VecRestoreArray2d - Restores a vector after VecGetArray2d() has been called.
2236: Logically Collective
2238: Input Parameters:
2239: + x - the vector
2240: . m - first dimension of two dimensional array
2241: . n - second dimension of the two dimensional array
2242: . mstart - first index you will use in first coordinate direction (often 0)
2243: . nstart - first index in the second coordinate direction (often 0)
2244: - a - location of pointer to array obtained from VecGetArray2d()
2246: Level: developer
2248: Notes:
2249: For regular PETSc vectors this routine does not involve any copies. For
2250: any special vectors that do not store local vector data in a contiguous
2251: array, this routine will copy the data back into the underlying
2252: vector data structure from the array obtained with VecGetArray().
2254: This routine actually zeros out the a pointer.
2256: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2257: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2258: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2259: @*/
2260: PetscErrorCode VecRestoreArray2d(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2261: {
2263: void *dummy;
2269: dummy = (void*)(*a + mstart);
2270: PetscFree(dummy);
2271: VecRestoreArray(x,NULL);
2272: return(0);
2273: }
2277: /*@C
2278: VecGetArray1d - Returns a pointer to a 1d contiguous array that contains this
2279: processor's portion of the vector data. You MUST call VecRestoreArray1d()
2280: when you no longer need access to the array.
2282: Logically Collective
2284: Input Parameter:
2285: + x - the vector
2286: . m - first dimension of two dimensional array
2287: - mstart - first index you will use in first coordinate direction (often 0)
2289: Output Parameter:
2290: . a - location to put pointer to the array
2292: Level: developer
2294: Notes:
2295: For a vector obtained from DMCreateLocalVector() mstart are likely
2296: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2297: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
2299: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2301: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2302: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2303: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2304: @*/
2305: PetscErrorCode VecGetArray1d(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2306: {
2308: PetscInt N;
2314: VecGetLocalSize(x,&N);
2315: if (m != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Local array size %D does not match 1d array dimensions %D",N,m);
2316: VecGetArray(x,a);
2317: *a -= mstart;
2318: return(0);
2319: }
2323: /*@C
2324: VecRestoreArray1d - Restores a vector after VecGetArray1d() has been called.
2326: Logically Collective
2328: Input Parameters:
2329: + x - the vector
2330: . m - first dimension of two dimensional array
2331: . mstart - first index you will use in first coordinate direction (often 0)
2332: - a - location of pointer to array obtained from VecGetArray21()
2334: Level: developer
2336: Notes:
2337: For regular PETSc vectors this routine does not involve any copies. For
2338: any special vectors that do not store local vector data in a contiguous
2339: array, this routine will copy the data back into the underlying
2340: vector data structure from the array obtained with VecGetArray1d().
2342: This routine actually zeros out the a pointer.
2344: Concepts: vector^accessing local values as 1d array
2346: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2347: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2348: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
2349: @*/
2350: PetscErrorCode VecRestoreArray1d(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2351: {
2357: VecRestoreArray(x,NULL);
2358: return(0);
2359: }
2364: /*@C
2365: VecGetArray3d - Returns a pointer to a 3d contiguous array that contains this
2366: processor's portion of the vector data. You MUST call VecRestoreArray3d()
2367: when you no longer need access to the array.
2369: Logically Collective
2371: Input Parameter:
2372: + x - the vector
2373: . m - first dimension of three dimensional array
2374: . n - second dimension of three dimensional array
2375: . p - third dimension of three dimensional array
2376: . mstart - first index you will use in first coordinate direction (often 0)
2377: . nstart - first index in the second coordinate direction (often 0)
2378: - pstart - first index in the third coordinate direction (often 0)
2380: Output Parameter:
2381: . a - location to put pointer to the array
2383: Level: developer
2385: Notes:
2386: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
2387: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2388: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2389: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
2391: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2393: Concepts: vector^accessing local values as 3d array
2395: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2396: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2397: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2398: @*/
2399: PetscErrorCode VecGetArray3d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
2400: {
2402: PetscInt i,N,j;
2403: PetscScalar *aa,**b;
2409: VecGetLocalSize(x,&N);
2410: 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);
2411: VecGetArray(x,&aa);
2413: PetscMalloc1(m*sizeof(PetscScalar**)+m*n,a);
2414: b = (PetscScalar**)((*a) + m);
2415: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
2416: for (i=0; i<m; i++)
2417: for (j=0; j<n; j++)
2418: b[i*n+j] = aa + i*n*p + j*p - pstart;
2420: *a -= mstart;
2421: return(0);
2422: }
2426: /*@C
2427: VecRestoreArray3d - Restores a vector after VecGetArray3d() has been called.
2429: Logically Collective
2431: Input Parameters:
2432: + x - the vector
2433: . m - first dimension of three dimensional array
2434: . n - second dimension of the three dimensional array
2435: . p - third dimension of the three dimensional array
2436: . mstart - first index you will use in first coordinate direction (often 0)
2437: . nstart - first index in the second coordinate direction (often 0)
2438: . pstart - first index in the third coordinate direction (often 0)
2439: - a - location of pointer to array obtained from VecGetArray3d()
2441: Level: developer
2443: Notes:
2444: For regular PETSc vectors this routine does not involve any copies. For
2445: any special vectors that do not store local vector data in a contiguous
2446: array, this routine will copy the data back into the underlying
2447: vector data structure from the array obtained with VecGetArray().
2449: This routine actually zeros out the a pointer.
2451: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2452: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2453: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
2454: @*/
2455: PetscErrorCode VecRestoreArray3d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
2456: {
2458: void *dummy;
2464: dummy = (void*)(*a + mstart);
2465: PetscFree(dummy);
2466: VecRestoreArray(x,NULL);
2467: return(0);
2468: }
2472: /*@C
2473: VecGetArray4d - Returns a pointer to a 4d contiguous array that contains this
2474: processor's portion of the vector data. You MUST call VecRestoreArray4d()
2475: when you no longer need access to the array.
2477: Logically Collective
2479: Input Parameter:
2480: + x - the vector
2481: . m - first dimension of four dimensional array
2482: . n - second dimension of four dimensional array
2483: . p - third dimension of four dimensional array
2484: . q - fourth dimension of four dimensional array
2485: . mstart - first index you will use in first coordinate direction (often 0)
2486: . nstart - first index in the second coordinate direction (often 0)
2487: . pstart - first index in the third coordinate direction (often 0)
2488: - qstart - first index in the fourth coordinate direction (often 0)
2490: Output Parameter:
2491: . a - location to put pointer to the array
2493: Level: beginner
2495: Notes:
2496: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
2497: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2498: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2499: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
2501: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2503: Concepts: vector^accessing local values as 3d array
2505: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2506: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2507: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2508: @*/
2509: PetscErrorCode VecGetArray4d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
2510: {
2512: PetscInt i,N,j,k;
2513: PetscScalar *aa,***b,**c;
2519: VecGetLocalSize(x,&N);
2520: 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);
2521: VecGetArray(x,&aa);
2523: PetscMalloc1(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p,a);
2524: b = (PetscScalar***)((*a) + m);
2525: c = (PetscScalar**)(b + m*n);
2526: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
2527: for (i=0; i<m; i++)
2528: for (j=0; j<n; j++)
2529: b[i*n+j] = c + i*n*p + j*p - pstart;
2530: for (i=0; i<m; i++)
2531: for (j=0; j<n; j++)
2532: for (k=0; k<p; k++)
2533: c[i*n*p+j*p+k] = aa + i*n*p*q + j*p*q + k*q - qstart;
2534: *a -= mstart;
2535: return(0);
2536: }
2540: /*@C
2541: VecRestoreArray4d - Restores a vector after VecGetArray3d() has been called.
2543: Logically Collective
2545: Input Parameters:
2546: + x - the vector
2547: . m - first dimension of four dimensional array
2548: . n - second dimension of the four dimensional array
2549: . p - third dimension of the four dimensional array
2550: . q - fourth dimension of the four dimensional array
2551: . mstart - first index you will use in first coordinate direction (often 0)
2552: . nstart - first index in the second coordinate direction (often 0)
2553: . pstart - first index in the third coordinate direction (often 0)
2554: . qstart - first index in the fourth coordinate direction (often 0)
2555: - a - location of pointer to array obtained from VecGetArray4d()
2557: Level: beginner
2559: Notes:
2560: For regular PETSc vectors this routine does not involve any copies. For
2561: any special vectors that do not store local vector data in a contiguous
2562: array, this routine will copy the data back into the underlying
2563: vector data structure from the array obtained with VecGetArray().
2565: This routine actually zeros out the a pointer.
2567: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2568: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2569: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
2570: @*/
2571: PetscErrorCode VecRestoreArray4d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
2572: {
2574: void *dummy;
2580: dummy = (void*)(*a + mstart);
2581: PetscFree(dummy);
2582: VecRestoreArray(x,NULL);
2583: return(0);
2584: }
2588: /*@C
2589: VecGetArray2dRead - Returns a pointer to a 2d contiguous array that contains this
2590: processor's portion of the vector data. You MUST call VecRestoreArray2dRead()
2591: when you no longer need access to the array.
2593: Logically Collective
2595: Input Parameter:
2596: + x - the vector
2597: . m - first dimension of two dimensional array
2598: . n - second dimension of two dimensional array
2599: . mstart - first index you will use in first coordinate direction (often 0)
2600: - nstart - first index in the second coordinate direction (often 0)
2602: Output Parameter:
2603: . a - location to put pointer to the array
2605: Level: developer
2607: Notes:
2608: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
2609: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2610: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2611: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
2613: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2615: Concepts: vector^accessing local values as 2d array
2617: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2618: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2619: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2620: @*/
2621: PetscErrorCode VecGetArray2dRead(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2622: {
2623: PetscErrorCode ierr;
2624: PetscInt i,N;
2625: const PetscScalar *aa;
2631: VecGetLocalSize(x,&N);
2632: 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);
2633: VecGetArrayRead(x,&aa);
2635: PetscMalloc1(m,a);
2636: for (i=0; i<m; i++) (*a)[i] = (PetscScalar*) aa + i*n - nstart;
2637: *a -= mstart;
2638: return(0);
2639: }
2643: /*@C
2644: VecRestoreArray2dRead - Restores a vector after VecGetArray2dRead() has been called.
2646: Logically Collective
2648: Input Parameters:
2649: + x - the vector
2650: . m - first dimension of two dimensional array
2651: . n - second dimension of the two dimensional array
2652: . mstart - first index you will use in first coordinate direction (often 0)
2653: . nstart - first index in the second coordinate direction (often 0)
2654: - a - location of pointer to array obtained from VecGetArray2d()
2656: Level: developer
2658: Notes:
2659: For regular PETSc vectors this routine does not involve any copies. For
2660: any special vectors that do not store local vector data in a contiguous
2661: array, this routine will copy the data back into the underlying
2662: vector data structure from the array obtained with VecGetArray().
2664: This routine actually zeros out the a pointer.
2666: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2667: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2668: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2669: @*/
2670: PetscErrorCode VecRestoreArray2dRead(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2671: {
2673: void *dummy;
2679: dummy = (void*)(*a + mstart);
2680: PetscFree(dummy);
2681: VecRestoreArrayRead(x,NULL);
2682: return(0);
2683: }
2687: /*@C
2688: VecGetArray1dRead - Returns a pointer to a 1d contiguous array that contains this
2689: processor's portion of the vector data. You MUST call VecRestoreArray1dRead()
2690: when you no longer need access to the array.
2692: Logically Collective
2694: Input Parameter:
2695: + x - the vector
2696: . m - first dimension of two dimensional array
2697: - mstart - first index you will use in first coordinate direction (often 0)
2699: Output Parameter:
2700: . a - location to put pointer to the array
2702: Level: developer
2704: Notes:
2705: For a vector obtained from DMCreateLocalVector() mstart are likely
2706: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2707: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
2709: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2711: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2712: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2713: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2714: @*/
2715: PetscErrorCode VecGetArray1dRead(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2716: {
2718: PetscInt N;
2724: VecGetLocalSize(x,&N);
2725: if (m != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Local array size %D does not match 1d array dimensions %D",N,m);
2726: VecGetArrayRead(x,(const PetscScalar**)a);
2727: *a -= mstart;
2728: return(0);
2729: }
2733: /*@C
2734: VecRestoreArray1dRead - Restores a vector after VecGetArray1dRead() has been called.
2736: Logically Collective
2738: Input Parameters:
2739: + x - the vector
2740: . m - first dimension of two dimensional array
2741: . mstart - first index you will use in first coordinate direction (often 0)
2742: - a - location of pointer to array obtained from VecGetArray21()
2744: Level: developer
2746: Notes:
2747: For regular PETSc vectors this routine does not involve any copies. For
2748: any special vectors that do not store local vector data in a contiguous
2749: array, this routine will copy the data back into the underlying
2750: vector data structure from the array obtained with VecGetArray1dRead().
2752: This routine actually zeros out the a pointer.
2754: Concepts: vector^accessing local values as 1d array
2756: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2757: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2758: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
2759: @*/
2760: PetscErrorCode VecRestoreArray1dRead(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2761: {
2767: VecRestoreArrayRead(x,NULL);
2768: return(0);
2769: }
2774: /*@C
2775: VecGetArray3dRead - Returns a pointer to a 3d contiguous array that contains this
2776: processor's portion of the vector data. You MUST call VecRestoreArray3dRead()
2777: when you no longer need access to the array.
2779: Logically Collective
2781: Input Parameter:
2782: + x - the vector
2783: . m - first dimension of three dimensional array
2784: . n - second dimension of three dimensional array
2785: . p - third dimension of three dimensional array
2786: . mstart - first index you will use in first coordinate direction (often 0)
2787: . nstart - first index in the second coordinate direction (often 0)
2788: - pstart - first index in the third coordinate direction (often 0)
2790: Output Parameter:
2791: . a - location to put pointer to the array
2793: Level: developer
2795: Notes:
2796: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
2797: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2798: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2799: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3dRead().
2801: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2803: Concepts: vector^accessing local values as 3d array
2805: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2806: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2807: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2808: @*/
2809: PetscErrorCode VecGetArray3dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
2810: {
2811: PetscErrorCode ierr;
2812: PetscInt i,N,j;
2813: const PetscScalar *aa;
2814: PetscScalar **b;
2820: VecGetLocalSize(x,&N);
2821: 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);
2822: VecGetArrayRead(x,&aa);
2824: PetscMalloc1(m*sizeof(PetscScalar**)+m*n,a);
2825: b = (PetscScalar**)((*a) + m);
2826: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
2827: for (i=0; i<m; i++)
2828: for (j=0; j<n; j++)
2829: b[i*n+j] = (PetscScalar *)aa + i*n*p + j*p - pstart;
2831: *a -= mstart;
2832: return(0);
2833: }
2837: /*@C
2838: VecRestoreArray3dRead - Restores a vector after VecGetArray3dRead() has been called.
2840: Logically Collective
2842: Input Parameters:
2843: + x - the vector
2844: . m - first dimension of three dimensional array
2845: . n - second dimension of the three dimensional array
2846: . p - third dimension of the three dimensional array
2847: . mstart - first index you will use in first coordinate direction (often 0)
2848: . nstart - first index in the second coordinate direction (often 0)
2849: . pstart - first index in the third coordinate direction (often 0)
2850: - a - location of pointer to array obtained from VecGetArray3dRead()
2852: Level: developer
2854: Notes:
2855: For regular PETSc vectors this routine does not involve any copies. For
2856: any special vectors that do not store local vector data in a contiguous
2857: array, this routine will copy the data back into the underlying
2858: vector data structure from the array obtained with VecGetArray().
2860: This routine actually zeros out the a pointer.
2862: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2863: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2864: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
2865: @*/
2866: PetscErrorCode VecRestoreArray3dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
2867: {
2869: void *dummy;
2875: dummy = (void*)(*a + mstart);
2876: PetscFree(dummy);
2877: VecRestoreArrayRead(x,NULL);
2878: return(0);
2879: }
2883: /*@C
2884: VecGetArray4dRead - Returns a pointer to a 4d contiguous array that contains this
2885: processor's portion of the vector data. You MUST call VecRestoreArray4dRead()
2886: when you no longer need access to the array.
2888: Logically Collective
2890: Input Parameter:
2891: + x - the vector
2892: . m - first dimension of four dimensional array
2893: . n - second dimension of four dimensional array
2894: . p - third dimension of four dimensional array
2895: . q - fourth dimension of four dimensional array
2896: . mstart - first index you will use in first coordinate direction (often 0)
2897: . nstart - first index in the second coordinate direction (often 0)
2898: . pstart - first index in the third coordinate direction (often 0)
2899: - qstart - first index in the fourth coordinate direction (often 0)
2901: Output Parameter:
2902: . a - location to put pointer to the array
2904: Level: beginner
2906: Notes:
2907: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
2908: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2909: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2910: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
2912: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2914: Concepts: vector^accessing local values as 3d array
2916: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2917: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2918: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2919: @*/
2920: PetscErrorCode VecGetArray4dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
2921: {
2922: PetscErrorCode ierr;
2923: PetscInt i,N,j,k;
2924: const PetscScalar *aa;
2925: PetscScalar ***b,**c;
2931: VecGetLocalSize(x,&N);
2932: 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);
2933: VecGetArrayRead(x,&aa);
2935: PetscMalloc1(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p,a);
2936: b = (PetscScalar***)((*a) + m);
2937: c = (PetscScalar**)(b + m*n);
2938: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
2939: for (i=0; i<m; i++)
2940: for (j=0; j<n; j++)
2941: b[i*n+j] = c + i*n*p + j*p - pstart;
2942: for (i=0; i<m; i++)
2943: for (j=0; j<n; j++)
2944: for (k=0; k<p; k++)
2945: c[i*n*p+j*p+k] = (PetscScalar*) aa + i*n*p*q + j*p*q + k*q - qstart;
2946: *a -= mstart;
2947: return(0);
2948: }
2952: /*@C
2953: VecRestoreArray4dRead - Restores a vector after VecGetArray3d() has been called.
2955: Logically Collective
2957: Input Parameters:
2958: + x - the vector
2959: . m - first dimension of four dimensional array
2960: . n - second dimension of the four dimensional array
2961: . p - third dimension of the four dimensional array
2962: . q - fourth dimension of the four dimensional array
2963: . mstart - first index you will use in first coordinate direction (often 0)
2964: . nstart - first index in the second coordinate direction (often 0)
2965: . pstart - first index in the third coordinate direction (often 0)
2966: . qstart - first index in the fourth coordinate direction (often 0)
2967: - a - location of pointer to array obtained from VecGetArray4dRead()
2969: Level: beginner
2971: Notes:
2972: For regular PETSc vectors this routine does not involve any copies. For
2973: any special vectors that do not store local vector data in a contiguous
2974: array, this routine will copy the data back into the underlying
2975: vector data structure from the array obtained with VecGetArray().
2977: This routine actually zeros out the a pointer.
2979: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2980: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2981: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
2982: @*/
2983: PetscErrorCode VecRestoreArray4dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
2984: {
2986: void *dummy;
2992: dummy = (void*)(*a + mstart);
2993: PetscFree(dummy);
2994: VecRestoreArrayRead(x,NULL);
2995: return(0);
2996: }
2998: #if defined(PETSC_USE_DEBUG)
3002: /*@
3003: VecLockGet - Gets the current lock status of a vector
3005: Logically Collective on Vec
3007: Input Parameter:
3008: . x - the vector
3010: Output Parameter:
3011: . state - greater than zero indicates the vector is still locked
3013: Level: beginner
3015: Concepts: vector^accessing local values
3017: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockPush(), VecLockGet()
3018: @*/
3019: PetscErrorCode VecLockGet(Vec x,PetscInt *state)
3020: {
3023: *state = x->lock;
3024: return(0);
3025: }
3029: /*@
3030: VecLockPush - Lock a vector from writing
3032: Logically Collective on Vec
3034: Input Parameter:
3035: . x - the vector
3037: Notes: If this is set then calls to VecGetArray() or VecSetValues() or any other routines that change the vectors values will fail.
3039: Call VecLockPop() to remove the latest lock
3041: Level: beginner
3043: Concepts: vector^accessing local values
3045: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockPop(), VecLockGet()
3046: @*/
3047: PetscErrorCode VecLockPush(Vec x)
3048: {
3051: x->lock++;
3052: return(0);
3053: }
3057: /*@
3058: VecLockPop - Unlock a vector from writing
3060: Logically Collective on Vec
3062: Input Parameter:
3063: . x - the vector
3065: Level: beginner
3067: Concepts: vector^accessing local values
3069: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockPush(), VecLockGet()
3070: @*/
3071: PetscErrorCode VecLockPop(Vec x)
3072: {
3075: x->lock--;
3076: if (x->lock < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Vector has been unlocked too many times");
3077: return(0);
3078: }
3080: #endif