Actual source code: rvector.c
petsc-3.9.4 2018-09-11
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>
7: static PetscInt VecGetSubVectorSavedStateId = -1;
9: PETSC_EXTERN PetscErrorCode VecValidValues(Vec vec,PetscInt argnum,PetscBool begin)
10: {
11: #if defined(PETSC_USE_DEBUG)
12: PetscErrorCode ierr;
13: PetscInt n,i;
14: const PetscScalar *x;
17: #if defined(PETSC_HAVE_VECCUDA) || defined(PETSC_HAVE_VIENNACL)
18: if ((vec->petscnative || vec->ops->getarray) && (vec->valid_GPU_array == PETSC_OFFLOAD_CPU || vec->valid_GPU_array == PETSC_OFFLOAD_BOTH)) {
19: #else
20: if (vec->petscnative || vec->ops->getarray) {
21: #endif
22: VecGetLocalSize(vec,&n);
23: VecGetArrayRead(vec,&x);
24: for (i=0; i<n; i++) {
25: if (begin) {
26: 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);
27: } else {
28: 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);
29: }
30: }
31: VecRestoreArrayRead(vec,&x);
32: }
33: return(0);
34: #else
35: return 0;
36: #endif
37: }
39: /*@
40: VecMaxPointwiseDivide - Computes the maximum of the componentwise division max = max_i abs(x_i/y_i).
42: Logically Collective on Vec
44: Input Parameters:
45: . x, y - the vectors
47: Output Parameter:
48: . max - the result
50: Level: advanced
52: Notes: x and y may be the same vector
53: if a particular y_i is zero, it is treated as 1 in the above formula
55: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMax(), VecPointwiseMin(), VecPointwiseMaxAbs()
56: @*/
57: PetscErrorCode VecMaxPointwiseDivide(Vec x,Vec y,PetscReal *max)
58: {
68: VecCheckSameSize(x,1,y,2);
69: (*x->ops->maxpointwisedivide)(x,y,max);
70: return(0);
71: }
73: /*@
74: VecDot - Computes the vector dot product.
76: Collective on Vec
78: Input Parameters:
79: . x, y - the vectors
81: Output Parameter:
82: . val - the dot product
84: Performance Issues:
85: $ per-processor memory bandwidth
86: $ interprocessor latency
87: $ work load inbalance that causes certain processes to arrive much earlier than others
89: Notes for Users of Complex Numbers:
90: For complex vectors, VecDot() computes
91: $ val = (x,y) = y^H x,
92: where y^H denotes the conjugate transpose of y. Note that this corresponds to the usual "mathematicians" complex
93: inner product where the SECOND argument gets the complex conjugate. Since the BLASdot() complex conjugates the first
94: first argument we call the BLASdot() with the arguments reversed.
96: Use VecTDot() for the indefinite form
97: $ val = (x,y) = y^T x,
98: where y^T denotes the transpose of y.
100: Level: intermediate
102: Concepts: inner product
103: Concepts: vector^inner product
105: .seealso: VecMDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecDotRealPart()
106: @*/
107: PetscErrorCode VecDot(Vec x,Vec y,PetscScalar *val)
108: {
118: VecCheckSameSize(x,1,y,2);
120: PetscLogEventBarrierBegin(VEC_DotBarrier,x,y,0,0,PetscObjectComm((PetscObject)x));
121: (*x->ops->dot)(x,y,val);
122: PetscLogEventBarrierEnd(VEC_DotBarrier,x,y,0,0,PetscObjectComm((PetscObject)x));
123: return(0);
124: }
126: /*@
127: VecDotRealPart - Computes the real part of the vector dot product.
129: Collective on Vec
131: Input Parameters:
132: . x, y - the vectors
134: Output Parameter:
135: . val - the real part of the dot product;
137: Performance Issues:
138: $ per-processor memory bandwidth
139: $ interprocessor latency
140: $ work load inbalance that causes certain processes to arrive much earlier than others
142: Notes for Users of Complex Numbers:
143: See VecDot() for more details on the definition of the dot product for complex numbers
145: For real numbers this returns the same value as VecDot()
147: 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
148: the space R^{2n} (that is a vector of 2n components with the real or imaginary part of the complex numbers for components)
150: Developer Note: This is not currently optimized to compute only the real part of the dot product.
152: Level: intermediate
154: Concepts: inner product
155: Concepts: vector^inner product
157: .seealso: VecMDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecDot(), VecDotNorm2()
158: @*/
159: PetscErrorCode VecDotRealPart(Vec x,Vec y,PetscReal *val)
160: {
162: PetscScalar fdot;
165: VecDot(x,y,&fdot);
166: *val = PetscRealPart(fdot);
167: return(0);
168: }
170: /*@
171: VecNorm - Computes the vector norm.
173: Collective on Vec
175: Input Parameters:
176: + x - the vector
177: - type - one of NORM_1, NORM_2, NORM_INFINITY. Also available
178: NORM_1_AND_2, which computes both norms and stores them
179: in a two element array.
181: Output Parameter:
182: . val - the norm
184: Notes:
185: $ NORM_1 denotes sum_i |x_i|
186: $ NORM_2 denotes sqrt(sum_i |x_i|^2)
187: $ NORM_INFINITY denotes max_i |x_i|
189: For complex numbers NORM_1 will return the traditional 1 norm of the 2 norm of the complex numbers; that is the 1
190: norm of the absolutely values of the complex entries. In PETSc 3.6 and earlier releases it returned the 1 norm of
191: the 1 norm of the complex entries (what is returned by the BLAS routine asum()). Both are valid norms but most
192: people expect the former.
194: Level: intermediate
196: Performance Issues:
197: $ per-processor memory bandwidth
198: $ interprocessor latency
199: $ work load inbalance that causes certain processes to arrive much earlier than others
201: Concepts: norm
202: Concepts: vector^norm
204: .seealso: VecDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecNormAvailable(),
205: VecNormBegin(), VecNormEnd()
207: @*/
208: PetscErrorCode VecNorm(Vec x,NormType type,PetscReal *val)
209: {
210: PetscBool flg;
218: /*
219: * Cached data?
220: */
221: if (type!=NORM_1_AND_2) {
222: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[type],*val,flg);
223: if (flg) return(0);
224: }
225: PetscLogEventBarrierBegin(VEC_NormBarrier,x,0,0,0,PetscObjectComm((PetscObject)x));
226: (*x->ops->norm)(x,type,val);
227: PetscLogEventBarrierEnd(VEC_NormBarrier,x,0,0,0,PetscObjectComm((PetscObject)x));
229: if (type!=NORM_1_AND_2) {
230: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[type],*val);
231: }
232: return(0);
233: }
235: /*@
236: VecNormAvailable - Returns the vector norm if it is already known.
238: Not Collective
240: Input Parameters:
241: + x - the vector
242: - type - one of NORM_1, NORM_2, NORM_INFINITY. Also available
243: NORM_1_AND_2, which computes both norms and stores them
244: in a two element array.
246: Output Parameter:
247: + available - PETSC_TRUE if the val returned is valid
248: - val - the norm
250: Notes:
251: $ NORM_1 denotes sum_i |x_i|
252: $ NORM_2 denotes sqrt(sum_i (x_i)^2)
253: $ NORM_INFINITY denotes max_i |x_i|
255: Level: intermediate
257: Performance Issues:
258: $ per-processor memory bandwidth
259: $ interprocessor latency
260: $ work load inbalance that causes certain processes to arrive much earlier than others
262: Compile Option:
263: PETSC_HAVE_SLOW_BLAS_NORM2 will cause a C (loop unrolled) version of the norm to be used, rather
264: than the BLAS. This should probably only be used when one is using the FORTRAN BLAS routines
265: (as opposed to vendor provided) because the FORTRAN BLAS NRM2() routine is very slow.
267: Concepts: norm
268: Concepts: vector^norm
270: .seealso: VecDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecNorm()
271: VecNormBegin(), VecNormEnd()
273: @*/
274: PetscErrorCode VecNormAvailable(Vec x,NormType type,PetscBool *available,PetscReal *val)
275: {
283: *available = PETSC_FALSE;
284: if (type!=NORM_1_AND_2) {
285: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[type],*val,*available);
286: }
287: return(0);
288: }
290: /*@
291: VecNormalize - Normalizes a vector by 2-norm.
293: Collective on Vec
295: Input Parameters:
296: + x - the vector
298: Output Parameter:
299: . x - the normalized vector
300: - val - the vector norm before normalization
302: Level: intermediate
304: Concepts: vector^normalizing
305: Concepts: normalizing^vector
307: @*/
308: PetscErrorCode VecNormalize(Vec x,PetscReal *val)
309: {
311: PetscReal norm;
316: PetscLogEventBegin(VEC_Normalize,x,0,0,0);
317: VecNorm(x,NORM_2,&norm);
318: if (norm == 0.0) {
319: PetscInfo(x,"Vector of zero norm can not be normalized; Returning only the zero norm\n");
320: } else if (norm != 1.0) {
321: PetscScalar tmp = 1.0/norm;
322: VecScale(x,tmp);
323: }
324: if (val) *val = norm;
325: PetscLogEventEnd(VEC_Normalize,x,0,0,0);
326: return(0);
327: }
329: /*@C
330: VecMax - Determines the maximum vector component and its location.
332: Collective on Vec
334: Input Parameter:
335: . x - the vector
337: Output Parameters:
338: + val - the maximum component
339: - p - the location of val (pass NULL if you don't want this)
341: Notes:
342: Returns the value PETSC_MIN_REAL and p = -1 if the vector is of length 0.
344: Returns the smallest index with the maximum value
345: Level: intermediate
347: Concepts: maximum^of vector
348: Concepts: vector^maximum value
350: .seealso: VecNorm(), VecMin()
351: @*/
352: PetscErrorCode VecMax(Vec x,PetscInt *p,PetscReal *val)
353: {
360: PetscLogEventBegin(VEC_Max,x,0,0,0);
361: (*x->ops->max)(x,p,val);
362: PetscLogEventEnd(VEC_Max,x,0,0,0);
363: return(0);
364: }
366: /*@C
367: VecMin - Determines the minimum vector component and its location.
369: Collective on Vec
371: Input Parameters:
372: . x - the vector
374: Output Parameter:
375: + val - the minimum component
376: - p - the location of val (pass NULL if you don't want this location)
378: Level: intermediate
380: Notes:
381: Returns the value PETSC_MAX_REAL and p = -1 if the vector is of length 0.
383: This returns the smallest index with the minumum value
385: Concepts: minimum^of vector
386: Concepts: vector^minimum entry
388: .seealso: VecMax()
389: @*/
390: PetscErrorCode VecMin(Vec x,PetscInt *p,PetscReal *val)
391: {
398: PetscLogEventBegin(VEC_Min,x,0,0,0);
399: (*x->ops->min)(x,p,val);
400: PetscLogEventEnd(VEC_Min,x,0,0,0);
401: return(0);
402: }
404: /*@
405: VecTDot - Computes an indefinite vector dot product. That is, this
406: routine does NOT use the complex conjugate.
408: Collective on Vec
410: Input Parameters:
411: . x, y - the vectors
413: Output Parameter:
414: . val - the dot product
416: Notes for Users of Complex Numbers:
417: For complex vectors, VecTDot() computes the indefinite form
418: $ val = (x,y) = y^T x,
419: where y^T denotes the transpose of y.
421: Use VecDot() for the inner product
422: $ val = (x,y) = y^H x,
423: where y^H denotes the conjugate transpose of y.
425: Level: intermediate
427: Concepts: inner product^non-Hermitian
428: Concepts: vector^inner product
429: Concepts: non-Hermitian inner product
431: .seealso: VecDot(), VecMTDot()
432: @*/
433: PetscErrorCode VecTDot(Vec x,Vec y,PetscScalar *val)
434: {
444: VecCheckSameSize(x,1,y,2);
446: PetscLogEventBegin(VEC_TDot,x,y,0,0);
447: (*x->ops->tdot)(x,y,val);
448: PetscLogEventEnd(VEC_TDot,x,y,0,0);
449: return(0);
450: }
452: /*@
453: VecScale - Scales a vector.
455: Not collective on Vec
457: Input Parameters:
458: + x - the vector
459: - alpha - the scalar
461: Output Parameter:
462: . x - the scaled vector
464: Note:
465: For a vector with n components, VecScale() computes
466: $ x[i] = alpha * x[i], for i=1,...,n.
468: Level: intermediate
470: Concepts: vector^scaling
471: Concepts: scaling^vector
473: @*/
474: PetscErrorCode VecScale(Vec x, PetscScalar alpha)
475: {
476: PetscReal norms[4] = {0.0,0.0,0.0, 0.0};
477: PetscBool flgs[4];
479: PetscInt i;
484: if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
485: PetscLogEventBegin(VEC_Scale,x,0,0,0);
486: if (alpha != (PetscScalar)1.0) {
487: /* get current stashed norms */
488: for (i=0; i<4; i++) {
489: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[i],norms[i],flgs[i]);
490: }
491: (*x->ops->scale)(x,alpha);
492: PetscObjectStateIncrease((PetscObject)x);
493: /* put the scaled stashed norms back into the Vec */
494: for (i=0; i<4; i++) {
495: if (flgs[i]) {
496: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[i],PetscAbsScalar(alpha)*norms[i]);
497: }
498: }
499: }
500: PetscLogEventEnd(VEC_Scale,x,0,0,0);
501: return(0);
502: }
504: /*@
505: VecSet - Sets all components of a vector to a single scalar value.
507: Logically Collective on Vec
509: Input Parameters:
510: + x - the vector
511: - alpha - the scalar
513: Output Parameter:
514: . x - the vector
516: Note:
517: For a vector of dimension n, VecSet() computes
518: $ x[i] = alpha, for i=1,...,n,
519: so that all vector entries then equal the identical
520: scalar value, alpha. Use the more general routine
521: VecSetValues() to set different vector entries.
523: You CANNOT call this after you have called VecSetValues() but before you call
524: VecAssemblyBegin/End().
526: Level: beginner
528: .seealso VecSetValues(), VecSetValuesBlocked(), VecSetRandom()
530: Concepts: vector^setting to constant
532: @*/
533: PetscErrorCode VecSet(Vec x,PetscScalar alpha)
534: {
535: PetscReal val;
541: 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()");
543: VecLocked(x,1);
545: PetscLogEventBegin(VEC_Set,x,0,0,0);
546: (*x->ops->set)(x,alpha);
547: PetscLogEventEnd(VEC_Set,x,0,0,0);
548: PetscObjectStateIncrease((PetscObject)x);
550: /* norms can be simply set (if |alpha|*N not too large) */
551: val = PetscAbsScalar(alpha);
552: if (x->map->N == 0) {
553: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_1],0.0l);
554: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],0.0);
555: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_2],0.0);
556: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_FROBENIUS],0.0);
557: } else if (val > PETSC_MAX_REAL/x->map->N) {
558: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],val);
559: } else {
560: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_1],x->map->N * val);
561: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],val);
562: val = PetscSqrtReal((PetscReal)x->map->N) * val;
563: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_2],val);
564: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_FROBENIUS],val);
565: }
566: return(0);
567: }
570: /*@
571: VecAXPY - Computes y = alpha x + y.
573: Logically Collective on Vec
575: Input Parameters:
576: + alpha - the scalar
577: - x, y - the vectors
579: Output Parameter:
580: . y - output vector
582: Level: intermediate
584: Notes: x and y MUST be different vectors
586: Concepts: vector^BLAS
587: Concepts: BLAS
589: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY()
590: @*/
591: PetscErrorCode VecAXPY(Vec y,PetscScalar alpha,Vec x)
592: {
601: VecCheckSameSize(x,1,y,3);
602: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y cannot be the same vector");
604: VecLocked(y,1);
606: VecLockPush(x);
607: PetscLogEventBegin(VEC_AXPY,x,y,0,0);
608: (*y->ops->axpy)(y,alpha,x);
609: PetscLogEventEnd(VEC_AXPY,x,y,0,0);
610: VecLockPop(x);
611: PetscObjectStateIncrease((PetscObject)y);
612: return(0);
613: }
615: /*@
616: VecAXPBY - Computes y = alpha x + beta y.
618: Logically Collective on Vec
620: Input Parameters:
621: + alpha,beta - the scalars
622: - x, y - the vectors
624: Output Parameter:
625: . y - output vector
627: Level: intermediate
629: Notes: x and y MUST be different vectors
631: Concepts: BLAS
632: Concepts: vector^BLAS
634: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPY()
635: @*/
636: PetscErrorCode VecAXPBY(Vec y,PetscScalar alpha,PetscScalar beta,Vec x)
637: {
646: VecCheckSameSize(x,1,y,4);
647: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y cannot be the same vector");
651: PetscLogEventBegin(VEC_AXPY,x,y,0,0);
652: (*y->ops->axpby)(y,alpha,beta,x);
653: PetscLogEventEnd(VEC_AXPY,x,y,0,0);
654: PetscObjectStateIncrease((PetscObject)y);
655: return(0);
656: }
658: /*@
659: VecAXPBYPCZ - Computes z = alpha x + beta y + gamma z
661: Logically Collective on Vec
663: Input Parameters:
664: + alpha,beta, gamma - the scalars
665: - x, y, z - the vectors
667: Output Parameter:
668: . z - output vector
670: Level: intermediate
672: Notes: x, y and z must be different vectors
674: Developer Note: alpha = 1 or gamma = 1 or gamma = 0.0 are handled as special cases
676: Concepts: BLAS
677: Concepts: vector^BLAS
679: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPY()
680: @*/
681: PetscErrorCode VecAXPBYPCZ(Vec z,PetscScalar alpha,PetscScalar beta,PetscScalar gamma,Vec x,Vec y)
682: {
694: VecCheckSameSize(x,1,y,5);
695: VecCheckSameSize(x,1,z,6);
696: if (x == y || x == z) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x, y, and z must be different vectors");
697: if (y == z) SETERRQ(PetscObjectComm((PetscObject)y),PETSC_ERR_ARG_IDN,"x, y, and z must be different vectors");
702: PetscLogEventBegin(VEC_AXPBYPCZ,x,y,z,0);
703: (*y->ops->axpbypcz)(z,alpha,beta,gamma,x,y);
704: PetscLogEventEnd(VEC_AXPBYPCZ,x,y,z,0);
705: PetscObjectStateIncrease((PetscObject)z);
706: return(0);
707: }
709: /*@
710: VecAYPX - Computes y = x + alpha y.
712: Logically Collective on Vec
714: Input Parameters:
715: + alpha - the scalar
716: - x, y - the vectors
718: Output Parameter:
719: . y - output vector
721: Level: intermediate
723: Notes: x and y MUST be different vectors
725: Concepts: vector^BLAS
726: Concepts: BLAS
728: .seealso: VecAXPY(), VecWAXPY()
729: @*/
730: PetscErrorCode VecAYPX(Vec y,PetscScalar alpha,Vec x)
731: {
740: VecCheckSameSize(x,1,y,3);
741: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y must be different vectors");
744: PetscLogEventBegin(VEC_AYPX,x,y,0,0);
745: (*y->ops->aypx)(y,alpha,x);
746: PetscLogEventEnd(VEC_AYPX,x,y,0,0);
747: PetscObjectStateIncrease((PetscObject)y);
748: return(0);
749: }
752: /*@
753: VecWAXPY - Computes w = alpha x + y.
755: Logically Collective on Vec
757: Input Parameters:
758: + alpha - the scalar
759: - x, y - the vectors
761: Output Parameter:
762: . w - the result
764: Level: intermediate
766: Notes: w cannot be either x or y, but x and y can be the same
768: Concepts: vector^BLAS
769: Concepts: BLAS
771: .seealso: VecAXPY(), VecAYPX(), VecAXPBY()
772: @*/
773: PetscErrorCode VecWAXPY(Vec w,PetscScalar alpha,Vec x,Vec y)
774: {
786: VecCheckSameSize(x,3,y,4);
787: VecCheckSameSize(x,3,w,1);
788: if (w == y) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Result vector w cannot be same as input vector y, suggest VecAXPY()");
789: if (w == x) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Result vector w cannot be same as input vector x, suggest VecAYPX()");
792: PetscLogEventBegin(VEC_WAXPY,x,y,w,0);
793: (*w->ops->waxpy)(w,alpha,x,y);
794: PetscLogEventEnd(VEC_WAXPY,x,y,w,0);
795: PetscObjectStateIncrease((PetscObject)w);
796: return(0);
797: }
800: /*@C
801: VecSetValues - Inserts or adds values into certain locations of a vector.
803: Not Collective
805: Input Parameters:
806: + x - vector to insert in
807: . ni - number of elements to add
808: . ix - indices where to add
809: . y - array of values
810: - iora - either INSERT_VALUES or ADD_VALUES, where
811: ADD_VALUES adds values to any existing entries, and
812: INSERT_VALUES replaces existing entries with new values
814: Notes:
815: VecSetValues() sets x[ix[i]] = y[i], for i=0,...,ni-1.
817: Calls to VecSetValues() with the INSERT_VALUES and ADD_VALUES
818: options cannot be mixed without intervening calls to the assembly
819: routines.
821: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
822: MUST be called after all calls to VecSetValues() have been completed.
824: VecSetValues() uses 0-based indices in Fortran as well as in C.
826: If you call VecSetOption(x, VEC_IGNORE_NEGATIVE_INDICES,PETSC_TRUE),
827: negative indices may be passed in ix. These rows are
828: simply ignored. This allows easily inserting element load matrices
829: with homogeneous Dirchlet boundary conditions that you don't want represented
830: in the vector.
832: Level: beginner
834: Concepts: vector^setting values
836: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValuesLocal(),
837: VecSetValue(), VecSetValuesBlocked(), InsertMode, INSERT_VALUES, ADD_VALUES, VecGetValues()
838: @*/
839: PetscErrorCode VecSetValues(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
840: {
845: if (!ni) return(0);
849: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
850: (*x->ops->setvalues)(x,ni,ix,y,iora);
851: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
852: PetscObjectStateIncrease((PetscObject)x);
853: return(0);
854: }
856: /*@
857: VecGetValues - Gets values from certain locations of a vector. Currently
858: can only get values on the same processor
860: Not Collective
862: Input Parameters:
863: + x - vector to get values from
864: . ni - number of elements to get
865: - ix - indices where to get them from (in global 1d numbering)
867: Output Parameter:
868: . y - array of values
870: Notes:
871: The user provides the allocated array y; it is NOT allocated in this routine
873: VecGetValues() gets y[i] = x[ix[i]], for i=0,...,ni-1.
875: VecAssemblyBegin() and VecAssemblyEnd() MUST be called before calling this
877: VecGetValues() uses 0-based indices in Fortran as well as in C.
879: If you call VecSetOption(x, VEC_IGNORE_NEGATIVE_INDICES,PETSC_TRUE),
880: negative indices may be passed in ix. These rows are
881: simply ignored.
883: Level: beginner
885: Concepts: vector^getting values
887: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues()
888: @*/
889: PetscErrorCode VecGetValues(Vec x,PetscInt ni,const PetscInt ix[],PetscScalar y[])
890: {
895: if (!ni) return(0);
899: (*x->ops->getvalues)(x,ni,ix,y);
900: return(0);
901: }
903: /*@C
904: VecSetValuesBlocked - Inserts or adds blocks of values into certain locations of a vector.
906: Not Collective
908: Input Parameters:
909: + x - vector to insert in
910: . ni - number of blocks to add
911: . ix - indices where to add in block count, rather than element count
912: . y - array of values
913: - iora - either INSERT_VALUES or ADD_VALUES, where
914: ADD_VALUES adds values to any existing entries, and
915: INSERT_VALUES replaces existing entries with new values
917: Notes:
918: VecSetValuesBlocked() sets x[bs*ix[i]+j] = y[bs*i+j],
919: for j=0,...,bs-1, for i=0,...,ni-1. where bs was set with VecSetBlockSize().
921: Calls to VecSetValuesBlocked() with the INSERT_VALUES and ADD_VALUES
922: options cannot be mixed without intervening calls to the assembly
923: routines.
925: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
926: MUST be called after all calls to VecSetValuesBlocked() have been completed.
928: VecSetValuesBlocked() uses 0-based indices in Fortran as well as in C.
930: Negative indices may be passed in ix, these rows are
931: simply ignored. This allows easily inserting element load matrices
932: with homogeneous Dirchlet boundary conditions that you don't want represented
933: in the vector.
935: Level: intermediate
937: Concepts: vector^setting values blocked
939: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValuesBlockedLocal(),
940: VecSetValues()
941: @*/
942: PetscErrorCode VecSetValuesBlocked(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
943: {
951: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
952: (*x->ops->setvaluesblocked)(x,ni,ix,y,iora);
953: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
954: PetscObjectStateIncrease((PetscObject)x);
955: return(0);
956: }
959: /*@C
960: VecSetValuesLocal - Inserts or adds values into certain locations of a vector,
961: using a local ordering of the nodes.
963: Not Collective
965: Input Parameters:
966: + x - vector to insert in
967: . ni - number of elements to add
968: . ix - indices where to add
969: . y - array of values
970: - iora - either INSERT_VALUES or ADD_VALUES, where
971: ADD_VALUES adds values to any existing entries, and
972: INSERT_VALUES replaces existing entries with new values
974: Level: intermediate
976: Notes:
977: VecSetValuesLocal() sets x[ix[i]] = y[i], for i=0,...,ni-1.
979: Calls to VecSetValues() with the INSERT_VALUES and ADD_VALUES
980: options cannot be mixed without intervening calls to the assembly
981: routines.
983: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
984: MUST be called after all calls to VecSetValuesLocal() have been completed.
986: VecSetValuesLocal() uses 0-based indices in Fortran as well as in C.
988: Concepts: vector^setting values with local numbering
990: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetLocalToGlobalMapping(),
991: VecSetValuesBlockedLocal()
992: @*/
993: PetscErrorCode VecSetValuesLocal(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
994: {
996: PetscInt lixp[128],*lix = lixp;
1000: if (!ni) return(0);
1005: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
1006: if (!x->ops->setvalueslocal) {
1007: if (!x->map->mapping) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Local to global never set with VecSetLocalToGlobalMapping()");
1008: if (ni > 128) {
1009: PetscMalloc1(ni,&lix);
1010: }
1011: ISLocalToGlobalMappingApply(x->map->mapping,ni,(PetscInt*)ix,lix);
1012: (*x->ops->setvalues)(x,ni,lix,y,iora);
1013: if (ni > 128) {
1014: PetscFree(lix);
1015: }
1016: } else {
1017: (*x->ops->setvalueslocal)(x,ni,ix,y,iora);
1018: }
1019: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
1020: PetscObjectStateIncrease((PetscObject)x);
1021: return(0);
1022: }
1024: /*@
1025: VecSetValuesBlockedLocal - Inserts or adds values into certain locations of a vector,
1026: using a local ordering of the nodes.
1028: Not Collective
1030: Input Parameters:
1031: + x - vector to insert in
1032: . ni - number of blocks to add
1033: . ix - indices where to add in block count, not element count
1034: . y - array of values
1035: - iora - either INSERT_VALUES or ADD_VALUES, where
1036: ADD_VALUES adds values to any existing entries, and
1037: INSERT_VALUES replaces existing entries with new values
1039: Level: intermediate
1041: Notes:
1042: VecSetValuesBlockedLocal() sets x[bs*ix[i]+j] = y[bs*i+j],
1043: for j=0,..bs-1, for i=0,...,ni-1, where bs has been set with VecSetBlockSize().
1045: Calls to VecSetValuesBlockedLocal() with the INSERT_VALUES and ADD_VALUES
1046: options cannot be mixed without intervening calls to the assembly
1047: routines.
1049: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
1050: MUST be called after all calls to VecSetValuesBlockedLocal() have been completed.
1052: VecSetValuesBlockedLocal() uses 0-based indices in Fortran as well as in C.
1055: Concepts: vector^setting values blocked with local numbering
1057: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetValuesBlocked(),
1058: VecSetLocalToGlobalMapping()
1059: @*/
1060: PetscErrorCode VecSetValuesBlockedLocal(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
1061: {
1063: PetscInt lixp[128],*lix = lixp;
1070: if (ni > 128) {
1071: PetscMalloc1(ni,&lix);
1072: }
1074: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
1075: ISLocalToGlobalMappingApplyBlock(x->map->mapping,ni,(PetscInt*)ix,lix);
1076: (*x->ops->setvaluesblocked)(x,ni,lix,y,iora);
1077: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
1078: if (ni > 128) {
1079: PetscFree(lix);
1080: }
1081: PetscObjectStateIncrease((PetscObject)x);
1082: return(0);
1083: }
1085: /*@
1086: VecMTDot - Computes indefinite vector multiple dot products.
1087: That is, it does NOT use the complex conjugate.
1089: Collective on Vec
1091: Input Parameters:
1092: + x - one vector
1093: . nv - number of vectors
1094: - y - array of vectors. Note that vectors are pointers
1096: Output Parameter:
1097: . val - array of the dot products
1099: Notes for Users of Complex Numbers:
1100: For complex vectors, VecMTDot() computes the indefinite form
1101: $ val = (x,y) = y^T x,
1102: where y^T denotes the transpose of y.
1104: Use VecMDot() for the inner product
1105: $ val = (x,y) = y^H x,
1106: where y^H denotes the conjugate transpose of y.
1108: Level: intermediate
1110: Concepts: inner product^multiple
1111: Concepts: vector^multiple inner products
1113: .seealso: VecMDot(), VecTDot()
1114: @*/
1115: PetscErrorCode VecMTDot(Vec x,PetscInt nv,const Vec y[],PetscScalar val[])
1116: {
1127: VecCheckSameSize(x,1,*y,3);
1129: PetscLogEventBegin(VEC_MTDot,x,*y,0,0);
1130: (*x->ops->mtdot)(x,nv,y,val);
1131: PetscLogEventEnd(VEC_MTDot,x,*y,0,0);
1132: return(0);
1133: }
1135: /*@
1136: VecMDot - Computes vector multiple dot products.
1138: Collective on Vec
1140: Input Parameters:
1141: + x - one vector
1142: . nv - number of vectors
1143: - y - array of vectors.
1145: Output Parameter:
1146: . val - array of the dot products (does not allocate the array)
1148: Notes for Users of Complex Numbers:
1149: For complex vectors, VecMDot() computes
1150: $ val = (x,y) = y^H x,
1151: where y^H denotes the conjugate transpose of y.
1153: Use VecMTDot() for the indefinite form
1154: $ val = (x,y) = y^T x,
1155: where y^T denotes the transpose of y.
1157: Level: intermediate
1159: Concepts: inner product^multiple
1160: Concepts: vector^multiple inner products
1162: .seealso: VecMTDot(), VecDot()
1163: @*/
1164: PetscErrorCode VecMDot(Vec x,PetscInt nv,const Vec y[],PetscScalar val[])
1165: {
1170: if (!nv) return(0);
1171: if (nv < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Number of vectors (given %D) cannot be negative",nv);
1178: VecCheckSameSize(x,1,*y,3);
1180: PetscLogEventBarrierBegin(VEC_MDotBarrier,x,*y,0,0,PetscObjectComm((PetscObject)x));
1181: (*x->ops->mdot)(x,nv,y,val);
1182: PetscLogEventBarrierEnd(VEC_MDotBarrier,x,*y,0,0,PetscObjectComm((PetscObject)x));
1183: return(0);
1184: }
1186: /*@
1187: VecMAXPY - Computes y = y + sum alpha[j] x[j]
1189: Logically Collective on Vec
1191: Input Parameters:
1192: + nv - number of scalars and x-vectors
1193: . alpha - array of scalars
1194: . y - one vector
1195: - x - array of vectors
1197: Level: intermediate
1199: Notes: y cannot be any of the x vectors
1201: Concepts: BLAS
1203: .seealso: VecAXPY(), VecWAXPY(), VecAYPX()
1204: @*/
1205: PetscErrorCode VecMAXPY(Vec y,PetscInt nv,const PetscScalar alpha[],Vec x[])
1206: {
1208: PetscInt i;
1212: if (!nv) return(0);
1213: if (nv < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Number of vectors (given %D) cannot be negative",nv);
1220: VecCheckSameSize(y,1,*x,4);
1223: PetscLogEventBegin(VEC_MAXPY,*x,y,0,0);
1224: (*y->ops->maxpy)(y,nv,alpha,x);
1225: PetscLogEventEnd(VEC_MAXPY,*x,y,0,0);
1226: PetscObjectStateIncrease((PetscObject)y);
1227: return(0);
1228: }
1230: /*@
1231: VecGetSubVector - Gets a vector representing part of another vector
1233: Collective on IS (and Vec if nonlocal entries are needed)
1235: Input Arguments:
1236: + X - vector from which to extract a subvector
1237: - is - index set representing portion of X to extract
1239: Output Arguments:
1240: . Y - subvector corresponding to is
1242: Level: advanced
1244: Notes:
1245: The subvector Y should be returned with VecRestoreSubVector().
1247: This function may return a subvector without making a copy, therefore it is not safe to use the original vector while
1248: modifying the subvector. Other non-overlapping subvectors can still be obtained from X using this function.
1250: .seealso: MatCreateSubMatrix()
1251: @*/
1252: PetscErrorCode VecGetSubVector(Vec X,IS is,Vec *Y)
1253: {
1254: PetscErrorCode ierr;
1255: Vec Z;
1261: if (X->ops->getsubvector) {
1262: (*X->ops->getsubvector)(X,is,&Z);
1263: } else { /* Default implementation currently does no caching */
1264: PetscInt gstart,gend,start;
1265: PetscBool contiguous,gcontiguous;
1266: VecGetOwnershipRange(X,&gstart,&gend);
1267: ISContiguousLocal(is,gstart,gend,&start,&contiguous);
1268: MPIU_Allreduce(&contiguous,&gcontiguous,1,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)is));
1269: if (gcontiguous) { /* We can do a no-copy implementation */
1270: PetscInt n,N,bs;
1271: PetscMPIInt size;
1272: PetscInt state;
1274: ISGetLocalSize(is,&n);
1275: VecGetBlockSize(X,&bs);
1276: if (n%bs || bs == 1 || !n) bs = -1; /* Do not decide block size if we do not have to */
1277: MPI_Comm_size(PetscObjectComm((PetscObject)X),&size);
1278: VecLockGet(X,&state);
1279: if (state) {
1280: const PetscScalar *x;
1281: VecGetArrayRead(X,&x);
1282: if (size == 1) {
1283: VecCreateSeqWithArray(PetscObjectComm((PetscObject)X),bs,n,(PetscScalar*)x+start,&Z);
1284: } else {
1285: ISGetSize(is,&N);
1286: VecCreateMPIWithArray(PetscObjectComm((PetscObject)X),bs,n,N,(PetscScalar*)x+start,&Z);
1287: }
1288: VecRestoreArrayRead(X,&x);
1289: VecLockPush(Z);
1290: } else {
1291: PetscScalar *x;
1292: VecGetArray(X,&x);
1293: if (size == 1) {
1294: VecCreateSeqWithArray(PetscObjectComm((PetscObject)X),bs,n,x+start,&Z);
1295: } else {
1296: ISGetSize(is,&N);
1297: VecCreateMPIWithArray(PetscObjectComm((PetscObject)X),bs,n,N,x+start,&Z);
1298: }
1299: VecRestoreArray(X,&x);
1300: }
1301: } else { /* Have to create a scatter and do a copy */
1302: VecScatter scatter;
1303: PetscInt n,N;
1304: ISGetLocalSize(is,&n);
1305: ISGetSize(is,&N);
1306: VecCreate(PetscObjectComm((PetscObject)is),&Z);
1307: VecSetSizes(Z,n,N);
1308: VecSetType(Z,((PetscObject)X)->type_name);
1309: VecScatterCreate(X,is,Z,NULL,&scatter);
1310: VecScatterBegin(scatter,X,Z,INSERT_VALUES,SCATTER_FORWARD);
1311: VecScatterEnd(scatter,X,Z,INSERT_VALUES,SCATTER_FORWARD);
1312: PetscObjectCompose((PetscObject)Z,"VecGetSubVector_Scatter",(PetscObject)scatter);
1313: VecScatterDestroy(&scatter);
1314: }
1315: }
1316: /* Record the state when the subvector was gotten so we know whether its values need to be put back */
1317: if (VecGetSubVectorSavedStateId < 0) {PetscObjectComposedDataRegister(&VecGetSubVectorSavedStateId);}
1318: PetscObjectComposedDataSetInt((PetscObject)Z,VecGetSubVectorSavedStateId,1);
1319: *Y = Z;
1320: return(0);
1321: }
1323: /*@
1324: VecRestoreSubVector - Restores a subvector extracted using VecGetSubVector()
1326: Collective on IS (and Vec if nonlocal entries need to be written)
1328: Input Arguments:
1329: + X - vector from which subvector was obtained
1330: . is - index set representing the subset of X
1331: - Y - subvector being restored
1333: Level: advanced
1335: .seealso: VecGetSubVector()
1336: @*/
1337: PetscErrorCode VecRestoreSubVector(Vec X,IS is,Vec *Y)
1338: {
1346: if (X->ops->restoresubvector) {
1347: (*X->ops->restoresubvector)(X,is,Y);
1348: } else {
1349: PETSC_UNUSED PetscObjectState dummystate = 0;
1350: PetscBool valid;
1351: PetscObjectComposedDataGetInt((PetscObject)*Y,VecGetSubVectorSavedStateId,dummystate,valid);
1352: if (!valid) {
1353: VecScatter scatter;
1355: PetscObjectQuery((PetscObject)*Y,"VecGetSubVector_Scatter",(PetscObject*)&scatter);
1356: if (scatter) {
1357: VecScatterBegin(scatter,*Y,X,INSERT_VALUES,SCATTER_REVERSE);
1358: VecScatterEnd(scatter,*Y,X,INSERT_VALUES,SCATTER_REVERSE);
1359: }
1360: }
1361: VecDestroy(Y);
1362: }
1363: return(0);
1364: }
1366: /*@
1367: VecGetLocalVectorRead - Maps the local portion of a vector into a
1368: vector. You must call VecRestoreLocalVectorRead() when the local
1369: vector is no longer needed.
1371: Not collective.
1373: Input parameter:
1374: . v - The vector for which the local vector is desired.
1376: Output parameter:
1377: . w - Upon exit this contains the local vector.
1379: Level: beginner
1381: Notes:
1382: This function is similar to VecGetArrayRead() which maps the local
1383: portion into a raw pointer. VecGetLocalVectorRead() is usually
1384: almost as efficient as VecGetArrayRead() but in certain circumstances
1385: VecGetLocalVectorRead() can be much more efficient than
1386: VecGetArrayRead(). This is because the construction of a contiguous
1387: array representing the vector data required by VecGetArrayRead() can
1388: be an expensive operation for certain vector types. For example, for
1389: GPU vectors VecGetArrayRead() requires that the data between device
1390: and host is synchronized.
1392: Unlike VecGetLocalVector(), this routine is not collective and
1393: preserves cached information.
1395: .seealso: VecRestoreLocalVectorRead(), VecGetLocalVector(), VecGetArrayRead(), VecGetArray()
1396: @*/
1397: PetscErrorCode VecGetLocalVectorRead(Vec v,Vec w)
1398: {
1400: PetscScalar *a;
1405: VecCheckSameLocalSize(v,1,w,2);
1406: if (v->ops->getlocalvectorread) {
1407: (*v->ops->getlocalvectorread)(v,w);
1408: } else {
1409: VecGetArrayRead(v,(const PetscScalar**)&a);
1410: VecPlaceArray(w,a);
1411: }
1412: return(0);
1413: }
1415: /*@
1416: VecRestoreLocalVectorRead - Unmaps the local portion of a vector
1417: previously mapped into a vector using VecGetLocalVectorRead().
1419: Not collective.
1421: Input parameter:
1422: . v - The local portion of this vector was previously mapped into w using VecGetLocalVectorRead().
1423: . w - The vector into which the local portion of v was mapped.
1425: Level: beginner
1427: .seealso: VecGetLocalVectorRead(), VecGetLocalVector(), VecGetArrayRead(), VecGetArray()
1428: @*/
1429: PetscErrorCode VecRestoreLocalVectorRead(Vec v,Vec w)
1430: {
1432: PetscScalar *a;
1437: if (v->ops->restorelocalvectorread) {
1438: (*v->ops->restorelocalvectorread)(v,w);
1439: } else {
1440: VecGetArrayRead(w,(const PetscScalar**)&a);
1441: VecRestoreArrayRead(v,(const PetscScalar**)&a);
1442: VecResetArray(w);
1443: }
1444: return(0);
1445: }
1447: /*@
1448: VecGetLocalVector - Maps the local portion of a vector into a
1449: vector.
1451: Collective on v, not collective on w.
1453: Input parameter:
1454: . v - The vector for which the local vector is desired.
1456: Output parameter:
1457: . w - Upon exit this contains the local vector.
1459: Level: beginner
1461: Notes:
1462: This function is similar to VecGetArray() which maps the local
1463: portion into a raw pointer. VecGetLocalVector() is usually about as
1464: efficient as VecGetArray() but in certain circumstances
1465: VecGetLocalVector() can be much more efficient than VecGetArray().
1466: This is because the construction of a contiguous array representing
1467: the vector data required by VecGetArray() can be an expensive
1468: operation for certain vector types. For example, for GPU vectors
1469: VecGetArray() requires that the data between device and host is
1470: synchronized.
1472: .seealso: VecRestoreLocalVector(), VecGetLocalVectorRead(), VecGetArrayRead(), VecGetArray()
1473: @*/
1474: PetscErrorCode VecGetLocalVector(Vec v,Vec w)
1475: {
1477: PetscScalar *a;
1482: VecCheckSameLocalSize(v,1,w,2);
1483: if (v->ops->getlocalvector) {
1484: (*v->ops->getlocalvector)(v,w);
1485: } else {
1486: VecGetArray(v,&a);
1487: VecPlaceArray(w,a);
1488: }
1489: return(0);
1490: }
1492: /*@
1493: VecRestoreLocalVector - Unmaps the local portion of a vector
1494: previously mapped into a vector using VecGetLocalVector().
1496: Logically collective.
1498: Input parameter:
1499: . v - The local portion of this vector was previously mapped into w using VecGetLocalVector().
1500: . w - The vector into which the local portion of v was mapped.
1502: Level: beginner
1504: .seealso: VecGetLocalVector(), VecGetLocalVectorRead(), VecRestoreLocalVectorRead(), LocalVectorRead(), VecGetArrayRead(), VecGetArray()
1505: @*/
1506: PetscErrorCode VecRestoreLocalVector(Vec v,Vec w)
1507: {
1509: PetscScalar *a;
1514: if (v->ops->restorelocalvector) {
1515: (*v->ops->restorelocalvector)(v,w);
1516: } else {
1517: VecGetArray(w,&a);
1518: VecRestoreArray(v,&a);
1519: VecResetArray(w);
1520: }
1521: return(0);
1522: }
1524: /*@C
1525: VecGetArray - Returns a pointer to a contiguous array that contains this
1526: processor's portion of the vector data. For the standard PETSc
1527: vectors, VecGetArray() returns a pointer to the local data array and
1528: does not use any copies. If the underlying vector data is not stored
1529: in a contiguous array this routine will copy the data to a contiguous
1530: array and return a pointer to that. You MUST call VecRestoreArray()
1531: when you no longer need access to the array.
1533: Logically Collective on Vec
1535: Input Parameter:
1536: . x - the vector
1538: Output Parameter:
1539: . a - location to put pointer to the array
1541: Fortran Note:
1542: This routine is used differently from Fortran 77
1543: $ Vec x
1544: $ PetscScalar x_array(1)
1545: $ PetscOffset i_x
1546: $ PetscErrorCode ierr
1547: $ call VecGetArray(x,x_array,i_x,ierr)
1548: $
1549: $ Access first local entry in vector with
1550: $ value = x_array(i_x + 1)
1551: $
1552: $ ...... other code
1553: $ call VecRestoreArray(x,x_array,i_x,ierr)
1554: For Fortran 90 see VecGetArrayF90()
1556: See the Fortran chapter of the users manual and
1557: petsc/src/snes/examples/tutorials/ex5f.F for details.
1559: Level: beginner
1561: Concepts: vector^accessing local values
1563: .seealso: VecRestoreArray(), VecGetArrayRead(), VecGetArrays(), VecGetArrayF90(), VecGetArrayReadF90(), VecPlaceArray(), VecGetArray2d(),
1564: VecGetArrayPair(), VecRestoreArrayPair()
1565: @*/
1566: PetscErrorCode VecGetArray(Vec x,PetscScalar **a)
1567: {
1569: #if defined(PETSC_HAVE_VIENNACL)
1570: PetscBool is_viennacltype = PETSC_FALSE;
1571: #endif
1575: VecLocked(x,1);
1576: if (x->petscnative) {
1577: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_VECCUDA)
1578: if (x->valid_GPU_array == PETSC_OFFLOAD_GPU) {
1579: #if defined(PETSC_HAVE_VIENNACL)
1580: PetscObjectTypeCompareAny((PetscObject)x,&is_viennacltype,VECSEQVIENNACL,VECMPIVIENNACL,VECVIENNACL,"");
1581: if (is_viennacltype) {
1582: VecViennaCLCopyFromGPU(x);
1583: } else
1584: #endif
1585: {
1586: #if defined(PETSC_HAVE_VECCUDA)
1587: VecCUDACopyFromGPU(x);
1588: #endif
1589: }
1590: } else if (x->valid_GPU_array == PETSC_OFFLOAD_UNALLOCATED) {
1591: #if defined(PETSC_HAVE_VIENNACL)
1592: PetscObjectTypeCompareAny((PetscObject)x,&is_viennacltype,VECSEQVIENNACL,VECMPIVIENNACL,VECVIENNACL,"");
1593: if (is_viennacltype) {
1594: VecViennaCLAllocateCheckHost(x);
1595: } else
1596: #endif
1597: {
1598: #if defined(PETSC_HAVE_VECCUDA)
1599: VecCUDAAllocateCheckHost(x);
1600: #endif
1601: }
1602: }
1603: #endif
1604: *a = *((PetscScalar**)x->data);
1605: } else {
1606: (*x->ops->getarray)(x,a);
1607: }
1608: return(0);
1609: }
1611: /*@C
1612: VecGetArrayRead - Get read-only pointer to contiguous array containing this processor's portion of the vector data.
1614: Not Collective
1616: Input Parameters:
1617: . x - the vector
1619: Output Parameter:
1620: . a - the array
1622: Level: beginner
1624: Notes:
1625: The array must be returned using a matching call to VecRestoreArrayRead().
1627: Unlike VecGetArray(), this routine is not collective and preserves cached information like vector norms.
1629: Standard PETSc vectors use contiguous storage so that this routine does not perform a copy. Other vector
1630: implementations may require a copy, but must such implementations should cache the contiguous representation so that
1631: only one copy is performed when this routine is called multiple times in sequence.
1633: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair()
1634: @*/
1635: PetscErrorCode VecGetArrayRead(Vec x,const PetscScalar **a)
1636: {
1638: #if defined(PETSC_HAVE_VIENNACL)
1639: PetscBool is_viennacltype = PETSC_FALSE;
1640: #endif
1644: if (x->petscnative) {
1645: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_VECCUDA)
1646: if (x->valid_GPU_array == PETSC_OFFLOAD_GPU) {
1647: #if defined(PETSC_HAVE_VIENNACL)
1648: PetscObjectTypeCompareAny((PetscObject)x,&is_viennacltype,VECSEQVIENNACL,VECMPIVIENNACL,VECVIENNACL,"");
1649: if (is_viennacltype) {
1650: VecViennaCLCopyFromGPU(x);
1651: } else
1652: #endif
1653: {
1654: #if defined(PETSC_HAVE_VECCUDA)
1655: VecCUDACopyFromGPU(x);
1656: #endif
1657: }
1658: }
1659: #endif
1660: *a = *((PetscScalar **)x->data);
1661: } else if (x->ops->getarrayread) {
1662: (*x->ops->getarrayread)(x,a);
1663: } else {
1664: (*x->ops->getarray)(x,(PetscScalar**)a);
1665: }
1666: return(0);
1667: }
1669: /*@C
1670: VecGetArrays - Returns a pointer to the arrays in a set of vectors
1671: that were created by a call to VecDuplicateVecs(). You MUST call
1672: VecRestoreArrays() when you no longer need access to the array.
1674: Logically Collective on Vec
1676: Input Parameter:
1677: + x - the vectors
1678: - n - the number of vectors
1680: Output Parameter:
1681: . a - location to put pointer to the array
1683: Fortran Note:
1684: This routine is not supported in Fortran.
1686: Level: intermediate
1688: .seealso: VecGetArray(), VecRestoreArrays()
1689: @*/
1690: PetscErrorCode VecGetArrays(const Vec x[],PetscInt n,PetscScalar **a[])
1691: {
1693: PetscInt i;
1694: PetscScalar **q;
1700: if (n <= 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must get at least one array n = %D",n);
1701: PetscMalloc1(n,&q);
1702: for (i=0; i<n; ++i) {
1703: VecGetArray(x[i],&q[i]);
1704: }
1705: *a = q;
1706: return(0);
1707: }
1709: /*@C
1710: VecRestoreArrays - Restores a group of vectors after VecGetArrays()
1711: has been called.
1713: Logically Collective on Vec
1715: Input Parameters:
1716: + x - the vector
1717: . n - the number of vectors
1718: - a - location of pointer to arrays obtained from VecGetArrays()
1720: Notes:
1721: For regular PETSc vectors this routine does not involve any copies. For
1722: any special vectors that do not store local vector data in a contiguous
1723: array, this routine will copy the data back into the underlying
1724: vector data structure from the arrays obtained with VecGetArrays().
1726: Fortran Note:
1727: This routine is not supported in Fortran.
1729: Level: intermediate
1731: .seealso: VecGetArrays(), VecRestoreArray()
1732: @*/
1733: PetscErrorCode VecRestoreArrays(const Vec x[],PetscInt n,PetscScalar **a[])
1734: {
1736: PetscInt i;
1737: PetscScalar **q = *a;
1744: for (i=0; i<n; ++i) {
1745: VecRestoreArray(x[i],&q[i]);
1746: }
1747: PetscFree(q);
1748: return(0);
1749: }
1751: /*@C
1752: VecRestoreArray - Restores a vector after VecGetArray() has been called.
1754: Logically Collective on Vec
1756: Input Parameters:
1757: + x - the vector
1758: - a - location of pointer to array obtained from VecGetArray()
1760: Level: beginner
1762: Notes:
1763: For regular PETSc vectors this routine does not involve any copies. For
1764: any special vectors that do not store local vector data in a contiguous
1765: array, this routine will copy the data back into the underlying
1766: vector data structure from the array obtained with VecGetArray().
1768: This routine actually zeros out the a pointer. This is to prevent accidental
1769: us of the array after it has been restored. If you pass null for a it will
1770: not zero the array pointer a.
1772: Fortran Note:
1773: This routine is used differently from Fortran 77
1774: $ Vec x
1775: $ PetscScalar x_array(1)
1776: $ PetscOffset i_x
1777: $ PetscErrorCode ierr
1778: $ call VecGetArray(x,x_array,i_x,ierr)
1779: $
1780: $ Access first local entry in vector with
1781: $ value = x_array(i_x + 1)
1782: $
1783: $ ...... other code
1784: $ call VecRestoreArray(x,x_array,i_x,ierr)
1786: See the Fortran chapter of the users manual and
1787: petsc/src/snes/examples/tutorials/ex5f.F for details.
1788: For Fortran 90 see VecRestoreArrayF90()
1790: .seealso: VecGetArray(), VecRestoreArrayRead(), VecRestoreArrays(), VecRestoreArrayF90(), VecRestoreArrayReadF90(), VecPlaceArray(), VecRestoreArray2d(),
1791: VecGetArrayPair(), VecRestoreArrayPair()
1792: @*/
1793: PetscErrorCode VecRestoreArray(Vec x,PetscScalar **a)
1794: {
1799: if (x->petscnative) {
1800: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_VECCUDA)
1801: x->valid_GPU_array = PETSC_OFFLOAD_CPU;
1802: #endif
1803: } else {
1804: (*x->ops->restorearray)(x,a);
1805: }
1806: if (a) *a = NULL;
1807: PetscObjectStateIncrease((PetscObject)x);
1808: return(0);
1809: }
1811: /*@C
1812: VecRestoreArrayRead - Restore array obtained with VecGetArrayRead()
1814: Not Collective
1816: Input Parameters:
1817: + vec - the vector
1818: - array - the array
1820: Level: beginner
1822: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair()
1823: @*/
1824: PetscErrorCode VecRestoreArrayRead(Vec x,const PetscScalar **a)
1825: {
1830: if (x->petscnative) {
1831: /* nothing */
1832: } else if (x->ops->restorearrayread) {
1833: (*x->ops->restorearrayread)(x,a);
1834: } else {
1835: (*x->ops->restorearray)(x,(PetscScalar**)a);
1836: }
1837: if (a) *a = NULL;
1838: return(0);
1839: }
1841: /*@
1842: VecPlaceArray - Allows one to replace the array in a vector with an
1843: array provided by the user. This is useful to avoid copying an array
1844: into a vector.
1846: Not Collective
1848: Input Parameters:
1849: + vec - the vector
1850: - array - the array
1852: Notes:
1853: You can return to the original array with a call to VecResetArray()
1855: Level: developer
1857: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecResetArray()
1859: @*/
1860: PetscErrorCode VecPlaceArray(Vec vec,const PetscScalar array[])
1861: {
1868: if (vec->ops->placearray) {
1869: (*vec->ops->placearray)(vec,array);
1870: } else SETERRQ(PetscObjectComm((PetscObject)vec),PETSC_ERR_SUP,"Cannot place array in this type of vector");
1871: PetscObjectStateIncrease((PetscObject)vec);
1872: return(0);
1873: }
1875: /*@C
1876: VecReplaceArray - Allows one to replace the array in a vector with an
1877: array provided by the user. This is useful to avoid copying an array
1878: into a vector.
1880: Not Collective
1882: Input Parameters:
1883: + vec - the vector
1884: - array - the array
1886: Notes:
1887: This permanently replaces the array and frees the memory associated
1888: with the old array.
1890: The memory passed in MUST be obtained with PetscMalloc() and CANNOT be
1891: freed by the user. It will be freed when the vector is destroy.
1893: Not supported from Fortran
1895: Level: developer
1897: .seealso: VecGetArray(), VecRestoreArray(), VecPlaceArray(), VecResetArray()
1899: @*/
1900: PetscErrorCode VecReplaceArray(Vec vec,const PetscScalar array[])
1901: {
1907: if (vec->ops->replacearray) {
1908: (*vec->ops->replacearray)(vec,array);
1909: } else SETERRQ(PetscObjectComm((PetscObject)vec),PETSC_ERR_SUP,"Cannot replace array in this type of vector");
1910: PetscObjectStateIncrease((PetscObject)vec);
1911: return(0);
1912: }
1914: /*MC
1915: VecDuplicateVecsF90 - Creates several vectors of the same type as an existing vector
1916: and makes them accessible via a Fortran90 pointer.
1918: Synopsis:
1919: VecDuplicateVecsF90(Vec x,PetscInt n,{Vec, pointer :: y(:)},integer ierr)
1921: Collective on Vec
1923: Input Parameters:
1924: + x - a vector to mimic
1925: - n - the number of vectors to obtain
1927: Output Parameters:
1928: + y - Fortran90 pointer to the array of vectors
1929: - ierr - error code
1931: Example of Usage:
1932: .vb
1933: #include <petsc/finclude/petscvec.h>
1934: use petscvec
1936: Vec x
1937: Vec, pointer :: y(:)
1938: ....
1939: call VecDuplicateVecsF90(x,2,y,ierr)
1940: call VecSet(y(2),alpha,ierr)
1941: call VecSet(y(2),alpha,ierr)
1942: ....
1943: call VecDestroyVecsF90(2,y,ierr)
1944: .ve
1946: Notes:
1947: Not yet supported for all F90 compilers
1949: Use VecDestroyVecsF90() to free the space.
1951: Level: beginner
1953: .seealso: VecDestroyVecsF90(), VecDuplicateVecs()
1955: M*/
1957: /*MC
1958: VecRestoreArrayF90 - Restores a vector to a usable state after a call to
1959: VecGetArrayF90().
1961: Synopsis:
1962: VecRestoreArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
1964: Logically Collective on Vec
1966: Input Parameters:
1967: + x - vector
1968: - xx_v - the Fortran90 pointer to the array
1970: Output Parameter:
1971: . ierr - error code
1973: Example of Usage:
1974: .vb
1975: #include <petsc/finclude/petscvec.h>
1976: use petscvec
1978: PetscScalar, pointer :: xx_v(:)
1979: ....
1980: call VecGetArrayF90(x,xx_v,ierr)
1981: xx_v(3) = a
1982: call VecRestoreArrayF90(x,xx_v,ierr)
1983: .ve
1985: Level: beginner
1987: .seealso: VecGetArrayF90(), VecGetArray(), VecRestoreArray(), UsingFortran, VecRestoreArrayReadF90()
1989: M*/
1991: /*MC
1992: VecDestroyVecsF90 - Frees a block of vectors obtained with VecDuplicateVecsF90().
1994: Synopsis:
1995: VecDestroyVecsF90(PetscInt n,{Vec, pointer :: x(:)},PetscErrorCode ierr)
1997: Collective on Vec
1999: Input Parameters:
2000: + n - the number of vectors previously obtained
2001: - x - pointer to array of vector pointers
2003: Output Parameter:
2004: . ierr - error code
2006: Notes:
2007: Not yet supported for all F90 compilers
2009: Level: beginner
2011: .seealso: VecDestroyVecs(), VecDuplicateVecsF90()
2013: M*/
2015: /*MC
2016: VecGetArrayF90 - Accesses a vector array from Fortran90. For default PETSc
2017: vectors, VecGetArrayF90() returns a pointer to the local data array. Otherwise,
2018: this routine is implementation dependent. You MUST call VecRestoreArrayF90()
2019: when you no longer need access to the array.
2021: Synopsis:
2022: VecGetArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2024: Logically Collective on Vec
2026: Input Parameter:
2027: . x - vector
2029: Output Parameters:
2030: + xx_v - the Fortran90 pointer to the array
2031: - ierr - error code
2033: Example of Usage:
2034: .vb
2035: #include <petsc/finclude/petscvec.h>
2036: use petscvec
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: If you ONLY intend to read entries from the array and not change any entries you should use VecGetArrayReadF90().
2047: Level: beginner
2049: .seealso: VecRestoreArrayF90(), VecGetArray(), VecRestoreArray(), VecGetArrayReadF90(), UsingFortran
2051: M*/
2053: /*MC
2054: VecGetArrayReadF90 - Accesses a read only array from Fortran90. For default PETSc
2055: vectors, VecGetArrayF90() returns a pointer to the local data array. Otherwise,
2056: this routine is implementation dependent. You MUST call VecRestoreArrayReadF90()
2057: when you no longer need access to the array.
2059: Synopsis:
2060: VecGetArrayReadF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2062: Logically Collective on Vec
2064: Input Parameter:
2065: . x - vector
2067: Output Parameters:
2068: + xx_v - the Fortran90 pointer to the array
2069: - ierr - error code
2071: Example of Usage:
2072: .vb
2073: #include <petsc/finclude/petscvec.h>
2074: use petscvec
2076: PetscScalar, pointer :: xx_v(:)
2077: ....
2078: call VecGetArrayReadF90(x,xx_v,ierr)
2079: a = xx_v(3)
2080: call VecRestoreArrayReadF90(x,xx_v,ierr)
2081: .ve
2083: If you intend to write entries into the array you must use VecGetArrayF90().
2085: Level: beginner
2087: .seealso: VecRestoreArrayReadF90(), VecGetArray(), VecRestoreArray(), VecGetArrayRead(), VecRestoreArrayRead(), VecGetArrayF90(), UsingFortran
2089: M*/
2091: /*MC
2092: VecRestoreArrayReadF90 - Restores a readonly vector to a usable state after a call to
2093: VecGetArrayReadF90().
2095: Synopsis:
2096: VecRestoreArrayReadF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2098: Logically Collective on Vec
2100: Input Parameters:
2101: + x - vector
2102: - xx_v - the Fortran90 pointer to the array
2104: Output Parameter:
2105: . ierr - error code
2107: Example of Usage:
2108: .vb
2109: #include <petsc/finclude/petscvec.h>
2110: use petscvec
2112: PetscScalar, pointer :: xx_v(:)
2113: ....
2114: call VecGetArrayReadF90(x,xx_v,ierr)
2115: a = xx_v(3)
2116: call VecRestoreArrayReadF90(x,xx_v,ierr)
2117: .ve
2119: Level: beginner
2121: .seealso: VecGetArrayReadF90(), VecGetArray(), VecRestoreArray(), VecGetArrayRead(), VecRestoreArrayRead(),UsingFortran, VecRestoreArrayF90()
2123: M*/
2125: /*@C
2126: VecGetArray2d - Returns a pointer to a 2d contiguous array that contains this
2127: processor's portion of the vector data. You MUST call VecRestoreArray2d()
2128: when you no longer need access to the array.
2130: Logically Collective
2132: Input Parameter:
2133: + x - the vector
2134: . m - first dimension of two dimensional array
2135: . n - second dimension of two dimensional array
2136: . mstart - first index you will use in first coordinate direction (often 0)
2137: - nstart - first index in the second coordinate direction (often 0)
2139: Output Parameter:
2140: . a - location to put pointer to the array
2142: Level: developer
2144: Notes:
2145: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
2146: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2147: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2148: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
2150: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2152: Concepts: vector^accessing local values as 2d array
2154: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2155: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2156: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2157: @*/
2158: PetscErrorCode VecGetArray2d(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2159: {
2161: PetscInt i,N;
2162: PetscScalar *aa;
2168: VecGetLocalSize(x,&N);
2169: 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);
2170: VecGetArray(x,&aa);
2172: PetscMalloc1(m,a);
2173: for (i=0; i<m; i++) (*a)[i] = aa + i*n - nstart;
2174: *a -= mstart;
2175: return(0);
2176: }
2178: /*@C
2179: VecRestoreArray2d - Restores a vector after VecGetArray2d() has been called.
2181: Logically Collective
2183: Input Parameters:
2184: + x - the vector
2185: . m - first dimension of two dimensional array
2186: . n - second dimension of the two dimensional array
2187: . mstart - first index you will use in first coordinate direction (often 0)
2188: . nstart - first index in the second coordinate direction (often 0)
2189: - a - location of pointer to array obtained from VecGetArray2d()
2191: Level: developer
2193: Notes:
2194: For regular PETSc vectors this routine does not involve any copies. For
2195: any special vectors that do not store local vector data in a contiguous
2196: array, this routine will copy the data back into the underlying
2197: vector data structure from the array obtained with VecGetArray().
2199: This routine actually zeros out the a pointer.
2201: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2202: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2203: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2204: @*/
2205: PetscErrorCode VecRestoreArray2d(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2206: {
2208: void *dummy;
2214: dummy = (void*)(*a + mstart);
2215: PetscFree(dummy);
2216: VecRestoreArray(x,NULL);
2217: return(0);
2218: }
2220: /*@C
2221: VecGetArray1d - Returns a pointer to a 1d contiguous array that contains this
2222: processor's portion of the vector data. You MUST call VecRestoreArray1d()
2223: when you no longer need access to the array.
2225: Logically Collective
2227: Input Parameter:
2228: + x - the vector
2229: . m - first dimension of two dimensional array
2230: - mstart - first index you will use in first coordinate direction (often 0)
2232: Output Parameter:
2233: . a - location to put pointer to the array
2235: Level: developer
2237: Notes:
2238: For a vector obtained from DMCreateLocalVector() mstart are likely
2239: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2240: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
2242: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2244: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2245: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2246: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2247: @*/
2248: PetscErrorCode VecGetArray1d(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2249: {
2251: PetscInt N;
2257: VecGetLocalSize(x,&N);
2258: if (m != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Local array size %D does not match 1d array dimensions %D",N,m);
2259: VecGetArray(x,a);
2260: *a -= mstart;
2261: return(0);
2262: }
2264: /*@C
2265: VecRestoreArray1d - Restores a vector after VecGetArray1d() has been called.
2267: Logically Collective
2269: Input Parameters:
2270: + x - the vector
2271: . m - first dimension of two dimensional array
2272: . mstart - first index you will use in first coordinate direction (often 0)
2273: - a - location of pointer to array obtained from VecGetArray21()
2275: Level: developer
2277: Notes:
2278: For regular PETSc vectors this routine does not involve any copies. For
2279: any special vectors that do not store local vector data in a contiguous
2280: array, this routine will copy the data back into the underlying
2281: vector data structure from the array obtained with VecGetArray1d().
2283: This routine actually zeros out the a pointer.
2285: Concepts: vector^accessing local values as 1d array
2287: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2288: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2289: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
2290: @*/
2291: PetscErrorCode VecRestoreArray1d(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2292: {
2298: VecRestoreArray(x,NULL);
2299: return(0);
2300: }
2303: /*@C
2304: VecGetArray3d - Returns a pointer to a 3d contiguous array that contains this
2305: processor's portion of the vector data. You MUST call VecRestoreArray3d()
2306: when you no longer need access to the array.
2308: Logically Collective
2310: Input Parameter:
2311: + x - the vector
2312: . m - first dimension of three dimensional array
2313: . n - second dimension of three dimensional array
2314: . p - third dimension of three dimensional array
2315: . mstart - first index you will use in first coordinate direction (often 0)
2316: . nstart - first index in the second coordinate direction (often 0)
2317: - pstart - first index in the third coordinate direction (often 0)
2319: Output Parameter:
2320: . a - location to put pointer to the array
2322: Level: developer
2324: Notes:
2325: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
2326: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2327: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2328: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
2330: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2332: Concepts: vector^accessing local values as 3d array
2334: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2335: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2336: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2337: @*/
2338: PetscErrorCode VecGetArray3d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
2339: {
2341: PetscInt i,N,j;
2342: PetscScalar *aa,**b;
2348: VecGetLocalSize(x,&N);
2349: 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);
2350: VecGetArray(x,&aa);
2352: PetscMalloc1(m*sizeof(PetscScalar**)+m*n,a);
2353: b = (PetscScalar**)((*a) + m);
2354: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
2355: for (i=0; i<m; i++)
2356: for (j=0; j<n; j++)
2357: b[i*n+j] = aa + i*n*p + j*p - pstart;
2359: *a -= mstart;
2360: return(0);
2361: }
2363: /*@C
2364: VecRestoreArray3d - Restores a vector after VecGetArray3d() has been called.
2366: Logically Collective
2368: Input Parameters:
2369: + x - the vector
2370: . m - first dimension of three dimensional array
2371: . n - second dimension of the three dimensional array
2372: . p - third dimension of the three dimensional array
2373: . mstart - first index you will use in first coordinate direction (often 0)
2374: . nstart - first index in the second coordinate direction (often 0)
2375: . pstart - first index in the third coordinate direction (often 0)
2376: - a - location of pointer to array obtained from VecGetArray3d()
2378: Level: developer
2380: Notes:
2381: For regular PETSc vectors this routine does not involve any copies. For
2382: any special vectors that do not store local vector data in a contiguous
2383: array, this routine will copy the data back into the underlying
2384: vector data structure from the array obtained with VecGetArray().
2386: This routine actually zeros out the a pointer.
2388: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2389: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2390: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
2391: @*/
2392: PetscErrorCode VecRestoreArray3d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
2393: {
2395: void *dummy;
2401: dummy = (void*)(*a + mstart);
2402: PetscFree(dummy);
2403: VecRestoreArray(x,NULL);
2404: return(0);
2405: }
2407: /*@C
2408: VecGetArray4d - Returns a pointer to a 4d contiguous array that contains this
2409: processor's portion of the vector data. You MUST call VecRestoreArray4d()
2410: when you no longer need access to the array.
2412: Logically Collective
2414: Input Parameter:
2415: + x - the vector
2416: . m - first dimension of four dimensional array
2417: . n - second dimension of four dimensional array
2418: . p - third dimension of four dimensional array
2419: . q - fourth dimension of four dimensional array
2420: . mstart - first index you will use in first coordinate direction (often 0)
2421: . nstart - first index in the second coordinate direction (often 0)
2422: . pstart - first index in the third coordinate direction (often 0)
2423: - qstart - first index in the fourth coordinate direction (often 0)
2425: Output Parameter:
2426: . a - location to put pointer to the array
2428: Level: beginner
2430: Notes:
2431: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
2432: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2433: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2434: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
2436: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2438: Concepts: vector^accessing local values as 3d array
2440: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2441: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2442: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2443: @*/
2444: PetscErrorCode VecGetArray4d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
2445: {
2447: PetscInt i,N,j,k;
2448: PetscScalar *aa,***b,**c;
2454: VecGetLocalSize(x,&N);
2455: 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);
2456: VecGetArray(x,&aa);
2458: PetscMalloc1(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p,a);
2459: b = (PetscScalar***)((*a) + m);
2460: c = (PetscScalar**)(b + m*n);
2461: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
2462: for (i=0; i<m; i++)
2463: for (j=0; j<n; j++)
2464: b[i*n+j] = c + i*n*p + j*p - pstart;
2465: for (i=0; i<m; i++)
2466: for (j=0; j<n; j++)
2467: for (k=0; k<p; k++)
2468: c[i*n*p+j*p+k] = aa + i*n*p*q + j*p*q + k*q - qstart;
2469: *a -= mstart;
2470: return(0);
2471: }
2473: /*@C
2474: VecRestoreArray4d - Restores a vector after VecGetArray3d() has been called.
2476: Logically Collective
2478: Input Parameters:
2479: + x - the vector
2480: . m - first dimension of four dimensional array
2481: . n - second dimension of the four dimensional array
2482: . p - third dimension of the four dimensional array
2483: . q - fourth dimension of the four dimensional array
2484: . mstart - first index you will use in first coordinate direction (often 0)
2485: . nstart - first index in the second coordinate direction (often 0)
2486: . pstart - first index in the third coordinate direction (often 0)
2487: . qstart - first index in the fourth coordinate direction (often 0)
2488: - a - location of pointer to array obtained from VecGetArray4d()
2490: Level: beginner
2492: Notes:
2493: For regular PETSc vectors this routine does not involve any copies. For
2494: any special vectors that do not store local vector data in a contiguous
2495: array, this routine will copy the data back into the underlying
2496: vector data structure from the array obtained with VecGetArray().
2498: This routine actually zeros out the a pointer.
2500: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2501: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2502: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
2503: @*/
2504: PetscErrorCode VecRestoreArray4d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
2505: {
2507: void *dummy;
2513: dummy = (void*)(*a + mstart);
2514: PetscFree(dummy);
2515: VecRestoreArray(x,NULL);
2516: return(0);
2517: }
2519: /*@C
2520: VecGetArray2dRead - Returns a pointer to a 2d contiguous array that contains this
2521: processor's portion of the vector data. You MUST call VecRestoreArray2dRead()
2522: when you no longer need access to the array.
2524: Logically Collective
2526: Input Parameter:
2527: + x - the vector
2528: . m - first dimension of two dimensional array
2529: . n - second dimension of two dimensional array
2530: . mstart - first index you will use in first coordinate direction (often 0)
2531: - nstart - first index in the second coordinate direction (often 0)
2533: Output Parameter:
2534: . a - location to put pointer to the array
2536: Level: developer
2538: Notes:
2539: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
2540: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2541: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2542: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
2544: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2546: Concepts: vector^accessing local values as 2d array
2548: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2549: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2550: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2551: @*/
2552: PetscErrorCode VecGetArray2dRead(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2553: {
2554: PetscErrorCode ierr;
2555: PetscInt i,N;
2556: const PetscScalar *aa;
2562: VecGetLocalSize(x,&N);
2563: 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);
2564: VecGetArrayRead(x,&aa);
2566: PetscMalloc1(m,a);
2567: for (i=0; i<m; i++) (*a)[i] = (PetscScalar*) aa + i*n - nstart;
2568: *a -= mstart;
2569: return(0);
2570: }
2572: /*@C
2573: VecRestoreArray2dRead - Restores a vector after VecGetArray2dRead() has been called.
2575: Logically Collective
2577: Input Parameters:
2578: + x - the vector
2579: . m - first dimension of two dimensional array
2580: . n - second dimension of the two dimensional array
2581: . mstart - first index you will use in first coordinate direction (often 0)
2582: . nstart - first index in the second coordinate direction (often 0)
2583: - a - location of pointer to array obtained from VecGetArray2d()
2585: Level: developer
2587: Notes:
2588: For regular PETSc vectors this routine does not involve any copies. For
2589: any special vectors that do not store local vector data in a contiguous
2590: array, this routine will copy the data back into the underlying
2591: vector data structure from the array obtained with VecGetArray().
2593: This routine actually zeros out the a pointer.
2595: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2596: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2597: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2598: @*/
2599: PetscErrorCode VecRestoreArray2dRead(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2600: {
2602: void *dummy;
2608: dummy = (void*)(*a + mstart);
2609: PetscFree(dummy);
2610: VecRestoreArrayRead(x,NULL);
2611: return(0);
2612: }
2614: /*@C
2615: VecGetArray1dRead - Returns a pointer to a 1d contiguous array that contains this
2616: processor's portion of the vector data. You MUST call VecRestoreArray1dRead()
2617: when you no longer need access to the array.
2619: Logically Collective
2621: Input Parameter:
2622: + x - the vector
2623: . m - first dimension of two dimensional array
2624: - mstart - first index you will use in first coordinate direction (often 0)
2626: Output Parameter:
2627: . a - location to put pointer to the array
2629: Level: developer
2631: Notes:
2632: For a vector obtained from DMCreateLocalVector() mstart are likely
2633: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2634: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
2636: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2638: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2639: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2640: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2641: @*/
2642: PetscErrorCode VecGetArray1dRead(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2643: {
2645: PetscInt N;
2651: VecGetLocalSize(x,&N);
2652: if (m != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Local array size %D does not match 1d array dimensions %D",N,m);
2653: VecGetArrayRead(x,(const PetscScalar**)a);
2654: *a -= mstart;
2655: return(0);
2656: }
2658: /*@C
2659: VecRestoreArray1dRead - Restores a vector after VecGetArray1dRead() has been called.
2661: Logically Collective
2663: Input Parameters:
2664: + x - the vector
2665: . m - first dimension of two dimensional array
2666: . mstart - first index you will use in first coordinate direction (often 0)
2667: - a - location of pointer to array obtained from VecGetArray21()
2669: Level: developer
2671: Notes:
2672: For regular PETSc vectors this routine does not involve any copies. For
2673: any special vectors that do not store local vector data in a contiguous
2674: array, this routine will copy the data back into the underlying
2675: vector data structure from the array obtained with VecGetArray1dRead().
2677: This routine actually zeros out the a pointer.
2679: Concepts: vector^accessing local values as 1d array
2681: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2682: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2683: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
2684: @*/
2685: PetscErrorCode VecRestoreArray1dRead(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2686: {
2692: VecRestoreArrayRead(x,NULL);
2693: return(0);
2694: }
2697: /*@C
2698: VecGetArray3dRead - Returns a pointer to a 3d contiguous array that contains this
2699: processor's portion of the vector data. You MUST call VecRestoreArray3dRead()
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 three dimensional array
2707: . n - second dimension of three dimensional array
2708: . p - third dimension of three dimensional array
2709: . mstart - first index you will use in first coordinate direction (often 0)
2710: . nstart - first index in the second coordinate direction (often 0)
2711: - pstart - first index in the third coordinate direction (often 0)
2713: Output Parameter:
2714: . a - location to put pointer to the array
2716: Level: developer
2718: Notes:
2719: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
2720: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2721: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2722: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3dRead().
2724: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2726: Concepts: vector^accessing local values as 3d array
2728: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2729: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2730: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2731: @*/
2732: PetscErrorCode VecGetArray3dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
2733: {
2734: PetscErrorCode ierr;
2735: PetscInt i,N,j;
2736: const PetscScalar *aa;
2737: PetscScalar **b;
2743: VecGetLocalSize(x,&N);
2744: 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);
2745: VecGetArrayRead(x,&aa);
2747: PetscMalloc1(m*sizeof(PetscScalar**)+m*n,a);
2748: b = (PetscScalar**)((*a) + m);
2749: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
2750: for (i=0; i<m; i++)
2751: for (j=0; j<n; j++)
2752: b[i*n+j] = (PetscScalar *)aa + i*n*p + j*p - pstart;
2754: *a -= mstart;
2755: return(0);
2756: }
2758: /*@C
2759: VecRestoreArray3dRead - Restores a vector after VecGetArray3dRead() has been called.
2761: Logically Collective
2763: Input Parameters:
2764: + x - the vector
2765: . m - first dimension of three dimensional array
2766: . n - second dimension of the three dimensional array
2767: . p - third dimension of the three dimensional array
2768: . mstart - first index you will use in first coordinate direction (often 0)
2769: . nstart - first index in the second coordinate direction (often 0)
2770: . pstart - first index in the third coordinate direction (often 0)
2771: - a - location of pointer to array obtained from VecGetArray3dRead()
2773: Level: developer
2775: Notes:
2776: For regular PETSc vectors this routine does not involve any copies. For
2777: any special vectors that do not store local vector data in a contiguous
2778: array, this routine will copy the data back into the underlying
2779: vector data structure from the array obtained with VecGetArray().
2781: This routine actually zeros out the a pointer.
2783: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2784: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2785: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
2786: @*/
2787: PetscErrorCode VecRestoreArray3dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
2788: {
2790: void *dummy;
2796: dummy = (void*)(*a + mstart);
2797: PetscFree(dummy);
2798: VecRestoreArrayRead(x,NULL);
2799: return(0);
2800: }
2802: /*@C
2803: VecGetArray4dRead - Returns a pointer to a 4d contiguous array that contains this
2804: processor's portion of the vector data. You MUST call VecRestoreArray4dRead()
2805: when you no longer need access to the array.
2807: Logically Collective
2809: Input Parameter:
2810: + x - the vector
2811: . m - first dimension of four dimensional array
2812: . n - second dimension of four dimensional array
2813: . p - third dimension of four dimensional array
2814: . q - fourth dimension of four dimensional array
2815: . mstart - first index you will use in first coordinate direction (often 0)
2816: . nstart - first index in the second coordinate direction (often 0)
2817: . pstart - first index in the third coordinate direction (often 0)
2818: - qstart - first index in the fourth coordinate direction (often 0)
2820: Output Parameter:
2821: . a - location to put pointer to the array
2823: Level: beginner
2825: Notes:
2826: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
2827: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2828: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2829: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
2831: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2833: Concepts: vector^accessing local values as 3d array
2835: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2836: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2837: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2838: @*/
2839: PetscErrorCode VecGetArray4dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
2840: {
2841: PetscErrorCode ierr;
2842: PetscInt i,N,j,k;
2843: const PetscScalar *aa;
2844: PetscScalar ***b,**c;
2850: VecGetLocalSize(x,&N);
2851: 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);
2852: VecGetArrayRead(x,&aa);
2854: PetscMalloc1(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p,a);
2855: b = (PetscScalar***)((*a) + m);
2856: c = (PetscScalar**)(b + m*n);
2857: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
2858: for (i=0; i<m; i++)
2859: for (j=0; j<n; j++)
2860: b[i*n+j] = c + i*n*p + j*p - pstart;
2861: for (i=0; i<m; i++)
2862: for (j=0; j<n; j++)
2863: for (k=0; k<p; k++)
2864: c[i*n*p+j*p+k] = (PetscScalar*) aa + i*n*p*q + j*p*q + k*q - qstart;
2865: *a -= mstart;
2866: return(0);
2867: }
2869: /*@C
2870: VecRestoreArray4dRead - Restores a vector after VecGetArray3d() has been called.
2872: Logically Collective
2874: Input Parameters:
2875: + x - the vector
2876: . m - first dimension of four dimensional array
2877: . n - second dimension of the four dimensional array
2878: . p - third dimension of the four dimensional array
2879: . q - fourth dimension of the four dimensional array
2880: . mstart - first index you will use in first coordinate direction (often 0)
2881: . nstart - first index in the second coordinate direction (often 0)
2882: . pstart - first index in the third coordinate direction (often 0)
2883: . qstart - first index in the fourth coordinate direction (often 0)
2884: - a - location of pointer to array obtained from VecGetArray4dRead()
2886: Level: beginner
2888: Notes:
2889: For regular PETSc vectors this routine does not involve any copies. For
2890: any special vectors that do not store local vector data in a contiguous
2891: array, this routine will copy the data back into the underlying
2892: vector data structure from the array obtained with VecGetArray().
2894: This routine actually zeros out the a pointer.
2896: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2897: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2898: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
2899: @*/
2900: PetscErrorCode VecRestoreArray4dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
2901: {
2903: void *dummy;
2909: dummy = (void*)(*a + mstart);
2910: PetscFree(dummy);
2911: VecRestoreArrayRead(x,NULL);
2912: return(0);
2913: }
2915: #if defined(PETSC_USE_DEBUG)
2917: /*@
2918: VecLockGet - Gets the current lock status of a vector
2920: Logically Collective on Vec
2922: Input Parameter:
2923: . x - the vector
2925: Output Parameter:
2926: . state - greater than zero indicates the vector is still locked
2928: Level: beginner
2930: Concepts: vector^accessing local values
2932: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockPush(), VecLockGet()
2933: @*/
2934: PetscErrorCode VecLockGet(Vec x,PetscInt *state)
2935: {
2938: *state = x->lock;
2939: return(0);
2940: }
2942: /*@
2943: VecLockPush - Lock a vector from writing
2945: Logically Collective on Vec
2947: Input Parameter:
2948: . x - the vector
2950: Notes: If this is set then calls to VecGetArray() or VecSetValues() or any other routines that change the vectors values will fail.
2952: Call VecLockPop() to remove the latest lock
2954: Level: beginner
2956: Concepts: vector^accessing local values
2958: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockPop(), VecLockGet()
2959: @*/
2960: PetscErrorCode VecLockPush(Vec x)
2961: {
2964: x->lock++;
2965: return(0);
2966: }
2968: /*@
2969: VecLockPop - Unlock a vector from writing
2971: Logically Collective on Vec
2973: Input Parameter:
2974: . x - the vector
2976: Level: beginner
2978: Concepts: vector^accessing local values
2980: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockPush(), VecLockGet()
2981: @*/
2982: PetscErrorCode VecLockPop(Vec x)
2983: {
2986: x->lock--;
2987: if (x->lock < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Vector has been unlocked too many times");
2988: return(0);
2989: }
2991: #endif