Actual source code: rvector.c
1: /*
2: Provides the interface functions for vector operations that have PetscScalar/PetscReal in the signature
3: These are the vector functions the user calls.
4: */
5: #include "petsc/private/sfimpl.h"
6: #include "petscsystypes.h"
7: #include <petsc/private/vecimpl.h>
8: #if defined(PETSC_HAVE_CUDA)
9: #include <../src/vec/vec/impls/dvecimpl.h>
10: #include <petsc/private/cudavecimpl.h>
11: #endif
12: #if defined(PETSC_HAVE_HIP)
13: #include <../src/vec/vec/impls/dvecimpl.h>
14: #include <petsc/private/hipvecimpl.h>
15: #endif
16: static PetscInt VecGetSubVectorSavedStateId = -1;
18: PETSC_EXTERN PetscErrorCode VecValidValues(Vec vec,PetscInt argnum,PetscBool begin)
19: {
20: #if defined(PETSC_USE_DEBUG)
21: PetscErrorCode ierr;
22: PetscInt n,i;
23: const PetscScalar *x;
26: #if defined(PETSC_HAVE_DEVICE)
27: if ((vec->petscnative || vec->ops->getarray) && (vec->offloadmask & PETSC_OFFLOAD_CPU)) {
28: #else
29: if (vec->petscnative || vec->ops->getarray) {
30: #endif
31: VecGetLocalSize(vec,&n);
32: VecGetArrayRead(vec,&x);
33: for (i=0; i<n; i++) {
34: if (begin) {
35: 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);
36: } else {
37: 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);
38: }
39: }
40: VecRestoreArrayRead(vec,&x);
41: }
42: #else
44: #endif
45: return(0);
46: }
48: /*@
49: VecMaxPointwiseDivide - Computes the maximum of the componentwise division max = max_i abs(x_i/y_i).
51: Logically Collective on Vec
53: Input Parameters:
54: . x, y - the vectors
56: Output Parameter:
57: . max - the result
59: Level: advanced
61: Notes:
62: x and y may be the same vector
63: if a particular y_i is zero, it is treated as 1 in the above formula
65: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMax(), VecPointwiseMin(), VecPointwiseMaxAbs()
66: @*/
67: PetscErrorCode VecMaxPointwiseDivide(Vec x,Vec y,PetscReal *max)
68: {
78: VecCheckSameSize(x,1,y,2);
79: (*x->ops->maxpointwisedivide)(x,y,max);
80: return(0);
81: }
83: /*@
84: VecDot - Computes the vector dot product.
86: Collective on Vec
88: Input Parameters:
89: . x, y - the vectors
91: Output Parameter:
92: . val - the dot product
94: Performance Issues:
95: $ per-processor memory bandwidth
96: $ interprocessor latency
97: $ work load inbalance that causes certain processes to arrive much earlier than others
99: Notes for Users of Complex Numbers:
100: For complex vectors, VecDot() computes
101: $ val = (x,y) = y^H x,
102: where y^H denotes the conjugate transpose of y. Note that this corresponds to the usual "mathematicians" complex
103: inner product where the SECOND argument gets the complex conjugate. Since the BLASdot() complex conjugates the first
104: first argument we call the BLASdot() with the arguments reversed.
106: Use VecTDot() for the indefinite form
107: $ val = (x,y) = y^T x,
108: where y^T denotes the transpose of y.
110: Level: intermediate
113: .seealso: VecMDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecDotRealPart()
114: @*/
115: PetscErrorCode VecDot(Vec x,Vec y,PetscScalar *val)
116: {
126: VecCheckSameSize(x,1,y,2);
128: PetscLogEventBegin(VEC_Dot,x,y,0,0);
129: (*x->ops->dot)(x,y,val);
130: PetscLogEventEnd(VEC_Dot,x,y,0,0);
131: return(0);
132: }
134: /*@
135: VecDotRealPart - Computes the real part of the vector dot product.
137: Collective on Vec
139: Input Parameters:
140: . x, y - the vectors
142: Output Parameter:
143: . val - the real part of the dot product;
145: Performance Issues:
146: $ per-processor memory bandwidth
147: $ interprocessor latency
148: $ work load inbalance that causes certain processes to arrive much earlier than others
150: Notes for Users of Complex Numbers:
151: See VecDot() for more details on the definition of the dot product for complex numbers
153: For real numbers this returns the same value as VecDot()
155: 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
156: the space R^{2n} (that is a vector of 2n components with the real or imaginary part of the complex numbers for components)
158: Developer Note: This is not currently optimized to compute only the real part of the dot product.
160: Level: intermediate
163: .seealso: VecMDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecDot(), VecDotNorm2()
164: @*/
165: PetscErrorCode VecDotRealPart(Vec x,Vec y,PetscReal *val)
166: {
168: PetscScalar fdot;
171: VecDot(x,y,&fdot);
172: *val = PetscRealPart(fdot);
173: return(0);
174: }
176: /*@
177: VecNorm - Computes the vector norm.
179: Collective on Vec
181: Input Parameters:
182: + x - the vector
183: - type - one of NORM_1, NORM_2, NORM_INFINITY. Also available
184: NORM_1_AND_2, which computes both norms and stores them
185: in a two element array.
187: Output Parameter:
188: . val - the norm
190: Notes:
191: $ NORM_1 denotes sum_i |x_i|
192: $ NORM_2 denotes sqrt(sum_i |x_i|^2)
193: $ NORM_INFINITY denotes max_i |x_i|
195: For complex numbers NORM_1 will return the traditional 1 norm of the 2 norm of the complex numbers; that is the 1
196: norm of the absolute values of the complex entries. In PETSc 3.6 and earlier releases it returned the 1 norm of
197: the 1 norm of the complex entries (what is returned by the BLAS routine asum()). Both are valid norms but most
198: people expect the former.
200: Level: intermediate
202: Performance Issues:
203: $ per-processor memory bandwidth
204: $ interprocessor latency
205: $ work load inbalance that causes certain processes to arrive much earlier than others
208: .seealso: VecDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecNormAvailable(),
209: VecNormBegin(), VecNormEnd()
211: @*/
213: PetscErrorCode VecNorm(Vec x,NormType type,PetscReal *val)
214: {
215: PetscBool flg;
223: /*
224: * Cached data?
225: */
226: if (type!=NORM_1_AND_2) {
227: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[type],*val,flg);
228: if (flg) return(0);
229: }
230: PetscLogEventBegin(VEC_Norm,x,0,0,0);
231: (*x->ops->norm)(x,type,val);
232: PetscLogEventEnd(VEC_Norm,x,0,0,0);
233: if (type!=NORM_1_AND_2) {
234: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[type],*val);
235: }
236: return(0);
237: }
239: /*@
240: VecNormAvailable - Returns the vector norm if it is already known.
242: Not Collective
244: Input Parameters:
245: + x - the vector
246: - type - one of NORM_1, NORM_2, NORM_INFINITY. Also available
247: NORM_1_AND_2, which computes both norms and stores them
248: in a two element array.
250: Output Parameter:
251: + available - PETSC_TRUE if the val returned is valid
252: - val - the norm
254: Notes:
255: $ NORM_1 denotes sum_i |x_i|
256: $ NORM_2 denotes sqrt(sum_i (x_i)^2)
257: $ NORM_INFINITY denotes max_i |x_i|
259: Level: intermediate
261: Performance Issues:
262: $ per-processor memory bandwidth
263: $ interprocessor latency
264: $ work load inbalance that causes certain processes to arrive much earlier than others
266: Compile Option:
267: PETSC_HAVE_SLOW_BLAS_NORM2 will cause a C (loop unrolled) version of the norm to be used, rather
268: than the BLAS. This should probably only be used when one is using the FORTRAN BLAS routines
269: (as opposed to vendor provided) because the FORTRAN BLAS NRM2() routine is very slow.
272: .seealso: VecDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecNorm()
273: VecNormBegin(), VecNormEnd()
275: @*/
276: PetscErrorCode VecNormAvailable(Vec x,NormType type,PetscBool *available,PetscReal *val)
277: {
285: *available = PETSC_FALSE;
286: if (type!=NORM_1_AND_2) {
287: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[type],*val,*available);
288: }
289: return(0);
290: }
292: /*@
293: VecNormalize - Normalizes a vector by 2-norm.
295: Collective on Vec
297: Input Parameters:
298: + x - the vector
300: Output Parameter:
301: . x - the normalized vector
302: - val - the vector norm before normalization
304: Level: intermediate
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 vector component with maximum real part and its location.
332: Collective on Vec
334: Input Parameter:
335: . x - the vector
337: Output Parameters:
338: + p - the location of val (pass NULL if you don't want this)
339: - val - the maximum component
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
348: .seealso: VecNorm(), VecMin()
349: @*/
350: PetscErrorCode VecMax(Vec x,PetscInt *p,PetscReal *val)
351: {
358: PetscLogEventBegin(VEC_Max,x,0,0,0);
359: (*x->ops->max)(x,p,val);
360: PetscLogEventEnd(VEC_Max,x,0,0,0);
361: return(0);
362: }
364: /*@C
365: VecMin - Determines the vector component with minimum real part and its location.
367: Collective on Vec
369: Input Parameters:
370: . x - the vector
372: Output Parameter:
373: + p - the location of val (pass NULL if you don't want this location)
374: - val - the minimum component
376: Level: intermediate
378: Notes:
379: Returns the value PETSC_MAX_REAL and p = -1 if the vector is of length 0.
381: This returns the smallest index with the minumum value
384: .seealso: VecMax()
385: @*/
386: PetscErrorCode VecMin(Vec x,PetscInt *p,PetscReal *val)
387: {
394: PetscLogEventBegin(VEC_Min,x,0,0,0);
395: (*x->ops->min)(x,p,val);
396: PetscLogEventEnd(VEC_Min,x,0,0,0);
397: return(0);
398: }
400: /*@
401: VecTDot - Computes an indefinite vector dot product. That is, this
402: routine does NOT use the complex conjugate.
404: Collective on Vec
406: Input Parameters:
407: . x, y - the vectors
409: Output Parameter:
410: . val - the dot product
412: Notes for Users of Complex Numbers:
413: For complex vectors, VecTDot() computes the indefinite form
414: $ val = (x,y) = y^T x,
415: where y^T denotes the transpose of y.
417: Use VecDot() for the inner product
418: $ val = (x,y) = y^H x,
419: where y^H denotes the conjugate transpose of y.
421: Level: intermediate
423: .seealso: VecDot(), VecMTDot()
424: @*/
425: PetscErrorCode VecTDot(Vec x,Vec y,PetscScalar *val)
426: {
436: VecCheckSameSize(x,1,y,2);
438: PetscLogEventBegin(VEC_TDot,x,y,0,0);
439: (*x->ops->tdot)(x,y,val);
440: PetscLogEventEnd(VEC_TDot,x,y,0,0);
441: return(0);
442: }
444: /*@
445: VecScale - Scales a vector.
447: Not collective on Vec
449: Input Parameters:
450: + x - the vector
451: - alpha - the scalar
453: Output Parameter:
454: . x - the scaled vector
456: Note:
457: For a vector with n components, VecScale() computes
458: $ x[i] = alpha * x[i], for i=1,...,n.
460: Level: intermediate
463: @*/
464: PetscErrorCode VecScale(Vec x, PetscScalar alpha)
465: {
466: PetscReal norms[4] = {0.0,0.0,0.0, 0.0};
467: PetscBool flgs[4];
469: PetscInt i;
474: if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
475: PetscLogEventBegin(VEC_Scale,x,0,0,0);
476: if (alpha != (PetscScalar)1.0) {
477: VecSetErrorIfLocked(x,1);
478: /* get current stashed norms */
479: for (i=0; i<4; i++) {
480: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[i],norms[i],flgs[i]);
481: }
482: (*x->ops->scale)(x,alpha);
483: PetscObjectStateIncrease((PetscObject)x);
484: /* put the scaled stashed norms back into the Vec */
485: for (i=0; i<4; i++) {
486: if (flgs[i]) {
487: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[i],PetscAbsScalar(alpha)*norms[i]);
488: }
489: }
490: }
491: PetscLogEventEnd(VEC_Scale,x,0,0,0);
492: return(0);
493: }
495: /*@
496: VecSet - Sets all components of a vector to a single scalar value.
498: Logically Collective on Vec
500: Input Parameters:
501: + x - the vector
502: - alpha - the scalar
504: Output Parameter:
505: . x - the vector
507: Note:
508: For a vector of dimension n, VecSet() computes
509: $ x[i] = alpha, for i=1,...,n,
510: so that all vector entries then equal the identical
511: scalar value, alpha. Use the more general routine
512: VecSetValues() to set different vector entries.
514: You CANNOT call this after you have called VecSetValues() but before you call
515: VecAssemblyBegin/End().
517: Level: beginner
519: .seealso VecSetValues(), VecSetValuesBlocked(), VecSetRandom()
521: @*/
522: PetscErrorCode VecSet(Vec x,PetscScalar alpha)
523: {
524: PetscReal val;
530: 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()");
532: VecSetErrorIfLocked(x,1);
534: PetscLogEventBegin(VEC_Set,x,0,0,0);
535: (*x->ops->set)(x,alpha);
536: PetscLogEventEnd(VEC_Set,x,0,0,0);
537: PetscObjectStateIncrease((PetscObject)x);
539: /* norms can be simply set (if |alpha|*N not too large) */
540: val = PetscAbsScalar(alpha);
541: if (x->map->N == 0) {
542: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_1],0.0l);
543: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],0.0);
544: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_2],0.0);
545: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_FROBENIUS],0.0);
546: } else if (val > PETSC_MAX_REAL/x->map->N) {
547: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],val);
548: } else {
549: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_1],x->map->N * val);
550: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],val);
551: val = PetscSqrtReal((PetscReal)x->map->N) * val;
552: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_2],val);
553: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_FROBENIUS],val);
554: }
555: return(0);
556: }
559: /*@
560: VecAXPY - Computes y = alpha x + y.
562: Logically Collective on Vec
564: Input Parameters:
565: + alpha - the scalar
566: - x, y - the vectors
568: Output Parameter:
569: . y - output vector
571: Level: intermediate
573: Notes:
574: x and y MUST be different vectors
575: This routine is optimized for alpha of 0.0, otherwise it calls the BLAS routine
577: $ VecAXPY(y,alpha,x) y = alpha x + y
578: $ VecAYPX(y,beta,x) y = x + beta y
579: $ VecAXPBY(y,alpha,beta,x) y = alpha x + beta y
580: $ VecWAXPY(w,alpha,x,y) w = alpha x + y
581: $ VecAXPBYPCZ(w,alpha,beta,gamma,x,y) z = alpha x + beta y + gamma z
582: $ VecMAXPY(y,nv,alpha[],x[]) y = sum alpha[i] x[i] + y
585: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPBYPCZ(), VecAXPBY()
586: @*/
587: PetscErrorCode VecAXPY(Vec y,PetscScalar alpha,Vec x)
588: {
597: VecCheckSameSize(x,3,y,1);
598: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y cannot be the same vector");
600: if (alpha == (PetscScalar)0.0) return(0);
601: VecSetErrorIfLocked(y,1);
603: VecLockReadPush(x);
604: PetscLogEventBegin(VEC_AXPY,x,y,0,0);
605: (*y->ops->axpy)(y,alpha,x);
606: PetscLogEventEnd(VEC_AXPY,x,y,0,0);
607: VecLockReadPop(x);
608: PetscObjectStateIncrease((PetscObject)y);
609: return(0);
610: }
612: /*@
613: VecAXPBY - Computes y = alpha x + beta y.
615: Logically Collective on Vec
617: Input Parameters:
618: + alpha,beta - the scalars
619: - x, y - the vectors
621: Output Parameter:
622: . y - output vector
624: Level: intermediate
626: Notes:
627: x and y MUST be different vectors
628: The implementation is optimized for alpha and/or beta values of 0.0 and 1.0
631: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPY(), VecAXPBYPCZ()
632: @*/
633: PetscErrorCode VecAXPBY(Vec y,PetscScalar alpha,PetscScalar beta,Vec x)
634: {
643: VecCheckSameSize(y,1,x,4);
644: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y cannot be the same vector");
647: if (alpha == (PetscScalar)0.0 && beta == (PetscScalar)1.0) return(0);
648: VecSetErrorIfLocked(y,1);
649: PetscLogEventBegin(VEC_AXPY,x,y,0,0);
650: (*y->ops->axpby)(y,alpha,beta,x);
651: PetscLogEventEnd(VEC_AXPY,x,y,0,0);
652: PetscObjectStateIncrease((PetscObject)y);
653: return(0);
654: }
656: /*@
657: VecAXPBYPCZ - Computes z = alpha x + beta y + gamma z
659: Logically Collective on Vec
661: Input Parameters:
662: + alpha,beta, gamma - the scalars
663: - x, y, z - the vectors
665: Output Parameter:
666: . z - output vector
668: Level: intermediate
670: Notes:
671: x, y and z must be different vectors
672: The implementation is optimized for alpha of 1.0 and gamma of 1.0 or 0.0
675: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPY(), VecAXPBY()
676: @*/
677: PetscErrorCode VecAXPBYPCZ(Vec z,PetscScalar alpha,PetscScalar beta,PetscScalar gamma,Vec x,Vec y)
678: {
690: VecCheckSameSize(x,1,y,5);
691: VecCheckSameSize(x,1,z,6);
692: if (x == y || x == z) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x, y, and z must be different vectors");
693: if (y == z) SETERRQ(PetscObjectComm((PetscObject)y),PETSC_ERR_ARG_IDN,"x, y, and z must be different vectors");
697: if (alpha == (PetscScalar)0.0 && beta == (PetscScalar)0.0 && gamma == (PetscScalar)1.0) return(0);
698: VecSetErrorIfLocked(z,1);
700: PetscLogEventBegin(VEC_AXPBYPCZ,x,y,z,0);
701: (*y->ops->axpbypcz)(z,alpha,beta,gamma,x,y);
702: PetscLogEventEnd(VEC_AXPBYPCZ,x,y,z,0);
703: PetscObjectStateIncrease((PetscObject)z);
704: return(0);
705: }
707: /*@
708: VecAYPX - Computes y = x + beta y.
710: Logically Collective on Vec
712: Input Parameters:
713: + beta - the scalar
714: - x, y - the vectors
716: Output Parameter:
717: . y - output vector
719: Level: intermediate
721: Notes:
722: x and y MUST be different vectors
723: The implementation is optimized for beta of -1.0, 0.0, and 1.0
726: .seealso: VecMAXPY(), VecWAXPY(), VecAXPY(), VecAXPBYPCZ(), VecAXPBY()
727: @*/
728: PetscErrorCode VecAYPX(Vec y,PetscScalar beta,Vec x)
729: {
738: VecCheckSameSize(x,1,y,3);
739: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y must be different vectors");
741: VecSetErrorIfLocked(y,1);
743: PetscLogEventBegin(VEC_AYPX,x,y,0,0);
744: (*y->ops->aypx)(y,beta,x);
745: PetscLogEventEnd(VEC_AYPX,x,y,0,0);
746: PetscObjectStateIncrease((PetscObject)y);
747: return(0);
748: }
751: /*@
752: VecWAXPY - Computes w = alpha x + y.
754: Logically Collective on Vec
756: Input Parameters:
757: + alpha - the scalar
758: - x, y - the vectors
760: Output Parameter:
761: . w - the result
763: Level: intermediate
765: Notes:
766: w cannot be either x or y, but x and y can be the same
767: The implementation is optimzed for alpha of -1.0, 0.0, and 1.0
770: .seealso: VecAXPY(), VecAYPX(), VecAXPBY(), VecMAXPY(), VecAXPBYPCZ()
771: @*/
772: PetscErrorCode VecWAXPY(Vec w,PetscScalar alpha,Vec x,Vec y)
773: {
785: VecCheckSameSize(x,3,y,4);
786: VecCheckSameSize(x,3,w,1);
787: if (w == y) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Result vector w cannot be same as input vector y, suggest VecAXPY()");
788: if (w == x) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Result vector w cannot be same as input vector x, suggest VecAYPX()");
790: VecSetErrorIfLocked(w,1);
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: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValuesLocal(),
835: VecSetValue(), VecSetValuesBlocked(), InsertMode, INSERT_VALUES, ADD_VALUES, VecGetValues()
836: @*/
837: PetscErrorCode VecSetValues(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
838: {
843: if (!ni) return(0);
848: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
849: (*x->ops->setvalues)(x,ni,ix,y,iora);
850: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
851: PetscObjectStateIncrease((PetscObject)x);
852: return(0);
853: }
855: /*@C
856: VecGetValues - Gets values from certain locations of a vector. Currently
857: can only get values on the same processor
859: Not Collective
861: Input Parameters:
862: + x - vector to get values from
863: . ni - number of elements to get
864: - ix - indices where to get them from (in global 1d numbering)
866: Output Parameter:
867: . y - array of values
869: Notes:
870: The user provides the allocated array y; it is NOT allocated in this routine
872: VecGetValues() gets y[i] = x[ix[i]], for i=0,...,ni-1.
874: VecAssemblyBegin() and VecAssemblyEnd() MUST be called before calling this
876: VecGetValues() uses 0-based indices in Fortran as well as in C.
878: If you call VecSetOption(x, VEC_IGNORE_NEGATIVE_INDICES,PETSC_TRUE),
879: negative indices may be passed in ix. These rows are
880: simply ignored.
882: Level: beginner
884: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues()
885: @*/
886: PetscErrorCode VecGetValues(Vec x,PetscInt ni,const PetscInt ix[],PetscScalar y[])
887: {
892: if (!ni) return(0);
896: (*x->ops->getvalues)(x,ni,ix,y);
897: return(0);
898: }
900: /*@C
901: VecSetValuesBlocked - Inserts or adds blocks of values into certain locations of a vector.
903: Not Collective
905: Input Parameters:
906: + x - vector to insert in
907: . ni - number of blocks to add
908: . ix - indices where to add in block count, rather than element count
909: . y - array of values
910: - iora - either INSERT_VALUES or ADD_VALUES, where
911: ADD_VALUES adds values to any existing entries, and
912: INSERT_VALUES replaces existing entries with new values
914: Notes:
915: VecSetValuesBlocked() sets x[bs*ix[i]+j] = y[bs*i+j],
916: for j=0,...,bs-1, for i=0,...,ni-1. where bs was set with VecSetBlockSize().
918: Calls to VecSetValuesBlocked() with the INSERT_VALUES and ADD_VALUES
919: options cannot be mixed without intervening calls to the assembly
920: routines.
922: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
923: MUST be called after all calls to VecSetValuesBlocked() have been completed.
925: VecSetValuesBlocked() uses 0-based indices in Fortran as well as in C.
927: Negative indices may be passed in ix, these rows are
928: simply ignored. This allows easily inserting element load matrices
929: with homogeneous Dirchlet boundary conditions that you don't want represented
930: in the vector.
932: Level: intermediate
934: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValuesBlockedLocal(),
935: VecSetValues()
936: @*/
937: PetscErrorCode VecSetValuesBlocked(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
938: {
943: if (!ni) return(0);
948: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
949: (*x->ops->setvaluesblocked)(x,ni,ix,y,iora);
950: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
951: PetscObjectStateIncrease((PetscObject)x);
952: return(0);
953: }
956: /*@C
957: VecSetValuesLocal - Inserts or adds values into certain locations of a vector,
958: using a local ordering of the nodes.
960: Not Collective
962: Input Parameters:
963: + x - vector to insert in
964: . ni - number of elements to add
965: . ix - indices where to add
966: . y - array of values
967: - iora - either INSERT_VALUES or ADD_VALUES, where
968: ADD_VALUES adds values to any existing entries, and
969: INSERT_VALUES replaces existing entries with new values
971: Level: intermediate
973: Notes:
974: VecSetValuesLocal() sets x[ix[i]] = y[i], for i=0,...,ni-1.
976: Calls to VecSetValues() with the INSERT_VALUES and ADD_VALUES
977: options cannot be mixed without intervening calls to the assembly
978: routines.
980: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
981: MUST be called after all calls to VecSetValuesLocal() have been completed.
983: VecSetValuesLocal() uses 0-based indices in Fortran as well as in C.
985: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetLocalToGlobalMapping(),
986: VecSetValuesBlockedLocal()
987: @*/
988: PetscErrorCode VecSetValuesLocal(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
989: {
991: PetscInt lixp[128],*lix = lixp;
995: if (!ni) return(0);
1000: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
1001: if (!x->ops->setvalueslocal) {
1002: if (!x->map->mapping) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Local to global never set with VecSetLocalToGlobalMapping()");
1003: if (ni > 128) {
1004: PetscMalloc1(ni,&lix);
1005: }
1006: ISLocalToGlobalMappingApply(x->map->mapping,ni,(PetscInt*)ix,lix);
1007: (*x->ops->setvalues)(x,ni,lix,y,iora);
1008: if (ni > 128) {
1009: PetscFree(lix);
1010: }
1011: } else {
1012: (*x->ops->setvalueslocal)(x,ni,ix,y,iora);
1013: }
1014: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
1015: PetscObjectStateIncrease((PetscObject)x);
1016: return(0);
1017: }
1019: /*@
1020: VecSetValuesBlockedLocal - Inserts or adds values into certain locations of a vector,
1021: using a local ordering of the nodes.
1023: Not Collective
1025: Input Parameters:
1026: + x - vector to insert in
1027: . ni - number of blocks to add
1028: . ix - indices where to add in block count, not element count
1029: . y - array of values
1030: - iora - either INSERT_VALUES or ADD_VALUES, where
1031: ADD_VALUES adds values to any existing entries, and
1032: INSERT_VALUES replaces existing entries with new values
1034: Level: intermediate
1036: Notes:
1037: VecSetValuesBlockedLocal() sets x[bs*ix[i]+j] = y[bs*i+j],
1038: for j=0,..bs-1, for i=0,...,ni-1, where bs has been set with VecSetBlockSize().
1040: Calls to VecSetValuesBlockedLocal() with the INSERT_VALUES and ADD_VALUES
1041: options cannot be mixed without intervening calls to the assembly
1042: routines.
1044: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
1045: MUST be called after all calls to VecSetValuesBlockedLocal() have been completed.
1047: VecSetValuesBlockedLocal() uses 0-based indices in Fortran as well as in C.
1050: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetValuesBlocked(),
1051: VecSetLocalToGlobalMapping()
1052: @*/
1053: PetscErrorCode VecSetValuesBlockedLocal(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
1054: {
1056: PetscInt lixp[128],*lix = lixp;
1060: if (!ni) return(0);
1064: if (ni > 128) {
1065: PetscMalloc1(ni,&lix);
1066: }
1068: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
1069: ISLocalToGlobalMappingApplyBlock(x->map->mapping,ni,(PetscInt*)ix,lix);
1070: (*x->ops->setvaluesblocked)(x,ni,lix,y,iora);
1071: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
1072: if (ni > 128) {
1073: PetscFree(lix);
1074: }
1075: PetscObjectStateIncrease((PetscObject)x);
1076: return(0);
1077: }
1079: /*@
1080: VecMTDot - Computes indefinite vector multiple dot products.
1081: That is, it does NOT use the complex conjugate.
1083: Collective on Vec
1085: Input Parameters:
1086: + x - one vector
1087: . nv - number of vectors
1088: - y - array of vectors. Note that vectors are pointers
1090: Output Parameter:
1091: . val - array of the dot products
1093: Notes for Users of Complex Numbers:
1094: For complex vectors, VecMTDot() computes the indefinite form
1095: $ val = (x,y) = y^T x,
1096: where y^T denotes the transpose of y.
1098: Use VecMDot() for the inner product
1099: $ val = (x,y) = y^H x,
1100: where y^H denotes the conjugate transpose of y.
1102: Level: intermediate
1105: .seealso: VecMDot(), VecTDot()
1106: @*/
1107: PetscErrorCode VecMTDot(Vec x,PetscInt nv,const Vec y[],PetscScalar val[])
1108: {
1114: if (!nv) return(0);
1121: VecCheckSameSize(x,1,*y,3);
1123: PetscLogEventBegin(VEC_MTDot,x,*y,0,0);
1124: (*x->ops->mtdot)(x,nv,y,val);
1125: PetscLogEventEnd(VEC_MTDot,x,*y,0,0);
1126: return(0);
1127: }
1129: /*@
1130: VecMDot - Computes vector multiple dot products.
1132: Collective on Vec
1134: Input Parameters:
1135: + x - one vector
1136: . nv - number of vectors
1137: - y - array of vectors.
1139: Output Parameter:
1140: . val - array of the dot products (does not allocate the array)
1142: Notes for Users of Complex Numbers:
1143: For complex vectors, VecMDot() computes
1144: $ val = (x,y) = y^H x,
1145: where y^H denotes the conjugate transpose of y.
1147: Use VecMTDot() for the indefinite form
1148: $ val = (x,y) = y^T x,
1149: where y^T denotes the transpose of y.
1151: Level: intermediate
1154: .seealso: VecMTDot(), VecDot()
1155: @*/
1156: PetscErrorCode VecMDot(Vec x,PetscInt nv,const Vec y[],PetscScalar val[])
1157: {
1163: if (!nv) return(0);
1164: if (nv < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Number of vectors (given %D) cannot be negative",nv);
1171: VecCheckSameSize(x,1,*y,3);
1173: PetscLogEventBegin(VEC_MDot,x,*y,0,0);
1174: (*x->ops->mdot)(x,nv,y,val);
1175: PetscLogEventEnd(VEC_MDot,x,*y,0,0);
1176: return(0);
1177: }
1179: /*@
1180: VecMAXPY - Computes y = y + sum alpha[i] x[i]
1182: Logically Collective on Vec
1184: Input Parameters:
1185: + nv - number of scalars and x-vectors
1186: . alpha - array of scalars
1187: . y - one vector
1188: - x - array of vectors
1190: Level: intermediate
1192: Notes:
1193: y cannot be any of the x vectors
1195: .seealso: VecAYPX(), VecWAXPY(), VecAXPY(), VecAXPBYPCZ(), VecAXPBY()
1196: @*/
1197: PetscErrorCode VecMAXPY(Vec y,PetscInt nv,const PetscScalar alpha[],Vec x[])
1198: {
1200: PetscInt i;
1201: PetscBool nonzero;
1206: if (!nv) return(0);
1207: if (nv < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Number of vectors (given %D) cannot be negative",nv);
1214: VecCheckSameSize(y,1,*x,4);
1216: for (i=0, nonzero = PETSC_FALSE; i<nv && !nonzero; i++) nonzero = (PetscBool)(nonzero || alpha[i] != (PetscScalar)0.0);
1217: if (!nonzero) return(0);
1218: VecSetErrorIfLocked(y,1);
1219: PetscLogEventBegin(VEC_MAXPY,*x,y,0,0);
1220: (*y->ops->maxpy)(y,nv,alpha,x);
1221: PetscLogEventEnd(VEC_MAXPY,*x,y,0,0);
1222: PetscObjectStateIncrease((PetscObject)y);
1223: return(0);
1224: }
1226: /*@
1227: VecConcatenate - Creates a new vector that is a vertical concatenation of all the given array of vectors
1228: in the order they appear in the array. The concatenated vector resides on the same
1229: communicator and is the same type as the source vectors.
1231: Collective on X
1233: Input Arguments:
1234: + nx - number of vectors to be concatenated
1235: - X - array containing the vectors to be concatenated in the order of concatenation
1237: Output Arguments:
1238: + Y - concatenated vector
1239: - x_is - array of index sets corresponding to the concatenated components of Y (NULL if not needed)
1241: Notes:
1242: Concatenation is similar to the functionality of a VecNest object; they both represent combination of
1243: different vector spaces. However, concatenated vectors do not store any information about their
1244: sub-vectors and own their own data. Consequently, this function provides index sets to enable the
1245: manipulation of data in the concatenated vector that corresponds to the original components at creation.
1247: This is a useful tool for outer loop algorithms, particularly constrained optimizers, where the solver
1248: has to operate on combined vector spaces and cannot utilize VecNest objects due to incompatibility with
1249: bound projections.
1251: Level: advanced
1253: .seealso: VECNEST, VECSCATTER, VecScatterCreate()
1254: @*/
1255: PetscErrorCode VecConcatenate(PetscInt nx, const Vec X[], Vec *Y, IS *x_is[])
1256: {
1257: MPI_Comm comm;
1258: VecType vec_type;
1259: Vec Ytmp, Xtmp;
1260: IS *is_tmp;
1261: PetscInt i, shift=0, Xnl, Xng, Xbegin;
1270: if ((*X)->ops->concatenate) {
1271: /* use the dedicated concatenation function if available */
1272: (*(*X)->ops->concatenate)(nx,X,Y,x_is);
1273: } else {
1274: /* loop over vectors and start creating IS */
1275: comm = PetscObjectComm((PetscObject)(*X));
1276: VecGetType(*X, &vec_type);
1277: PetscMalloc1(nx, &is_tmp);
1278: for (i=0; i<nx; i++) {
1279: VecGetSize(X[i], &Xng);
1280: VecGetLocalSize(X[i], &Xnl);
1281: VecGetOwnershipRange(X[i], &Xbegin, NULL);
1282: ISCreateStride(comm, Xnl, shift + Xbegin, 1, &is_tmp[i]);
1283: shift += Xng;
1284: }
1285: /* create the concatenated vector */
1286: VecCreate(comm, &Ytmp);
1287: VecSetType(Ytmp, vec_type);
1288: VecSetSizes(Ytmp, PETSC_DECIDE, shift);
1289: VecSetUp(Ytmp);
1290: /* copy data from X array to Y and return */
1291: for (i=0; i<nx; i++) {
1292: VecGetSubVector(Ytmp, is_tmp[i], &Xtmp);
1293: VecCopy(X[i], Xtmp);
1294: VecRestoreSubVector(Ytmp, is_tmp[i], &Xtmp);
1295: }
1296: *Y = Ytmp;
1297: if (x_is) {
1298: *x_is = is_tmp;
1299: } else {
1300: for (i=0; i<nx; i++) {
1301: ISDestroy(&is_tmp[i]);
1302: }
1303: PetscFree(is_tmp);
1304: }
1305: }
1306: return(0);
1307: }
1309: /*@
1310: VecGetSubVector - Gets a vector representing part of another vector
1312: Collective on X and IS
1314: Input Arguments:
1315: + X - vector from which to extract a subvector
1316: - is - index set representing portion of X to extract
1318: Output Arguments:
1319: . Y - subvector corresponding to is
1321: Level: advanced
1323: Notes:
1324: The subvector Y should be returned with VecRestoreSubVector().
1325: X and is must be defined on the same communicator
1327: This function may return a subvector without making a copy, therefore it is not safe to use the original vector while
1328: modifying the subvector. Other non-overlapping subvectors can still be obtained from X using this function.
1329: The resulting subvector inherits the block size from the IS if greater than one. Otherwise, the block size is guessed from the block size of the original vec.
1331: .seealso: MatCreateSubMatrix()
1332: @*/
1333: PetscErrorCode VecGetSubVector(Vec X,IS is,Vec *Y)
1334: {
1335: PetscErrorCode ierr;
1336: Vec Z;
1343: if (X->ops->getsubvector) {
1344: (*X->ops->getsubvector)(X,is,&Z);
1345: } else { /* Default implementation currently does no caching */
1346: PetscInt gstart,gend,start;
1347: PetscBool red[2] = { PETSC_TRUE, PETSC_TRUE };
1348: PetscInt n,N,ibs,vbs,bs = -1;
1350: ISGetLocalSize(is,&n);
1351: ISGetSize(is,&N);
1352: ISGetBlockSize(is,&ibs);
1353: VecGetBlockSize(X,&vbs);
1354: VecGetOwnershipRange(X,&gstart,&gend);
1355: ISContiguousLocal(is,gstart,gend,&start,&red[0]);
1356: /* block size is given by IS if ibs > 1; otherwise, check the vector */
1357: if (ibs > 1) {
1358: MPIU_Allreduce(MPI_IN_PLACE,red,1,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)is));
1359: bs = ibs;
1360: } else {
1361: if (n%vbs || vbs == 1) red[1] = PETSC_FALSE; /* this process invalidate the collectiveness of block size */
1362: MPIU_Allreduce(MPI_IN_PLACE,red,2,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)is));
1363: if (red[0] && red[1]) bs = vbs; /* all processes have a valid block size and the access will be contiguous */
1364: }
1365: if (red[0]) { /* We can do a no-copy implementation */
1366: const PetscScalar *x;
1367: PetscInt state = 0;
1368: PetscBool isstd,iscuda,iship;
1370: PetscObjectTypeCompareAny((PetscObject)X,&isstd,VECSEQ,VECMPI,VECSTANDARD,"");
1371: PetscObjectTypeCompareAny((PetscObject)X,&iscuda,VECSEQCUDA,VECMPICUDA,"");
1372: PetscObjectTypeCompareAny((PetscObject)X,&iship,VECSEQHIP,VECMPIHIP,"");
1373: if (iscuda) {
1374: #if defined(PETSC_HAVE_CUDA)
1375: const PetscScalar *x_d;
1376: PetscMPIInt size;
1377: PetscOffloadMask flg;
1379: VecCUDAGetArrays_Private(X,&x,&x_d,&flg);
1380: if (flg == PETSC_OFFLOAD_UNALLOCATED) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Not for PETSC_OFFLOAD_UNALLOCATED");
1381: if (n && !x && !x_d) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Missing vector data");
1382: if (x) x += start;
1383: if (x_d) x_d += start;
1384: MPI_Comm_size(PetscObjectComm((PetscObject)X),&size);
1385: if (size == 1) {
1386: VecCreateSeqCUDAWithArrays(PetscObjectComm((PetscObject)X),bs,n,x,x_d,&Z);
1387: } else {
1388: VecCreateMPICUDAWithArrays(PetscObjectComm((PetscObject)X),bs,n,N,x,x_d,&Z);
1389: }
1390: Z->offloadmask = flg;
1391: #endif
1392: } else if (iship) {
1393: #if defined(PETSC_HAVE_HIP)
1394: const PetscScalar *x_d;
1395: PetscMPIInt size;
1396: PetscOffloadMask flg;
1398: VecHIPGetArrays_Private(X,&x,&x_d,&flg);
1399: if (flg == PETSC_OFFLOAD_UNALLOCATED) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Not for PETSC_OFFLOAD_UNALLOCATED");
1400: if (n && !x && !x_d) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Missing vector data");
1401: if (x) x += start;
1402: if (x_d) x_d += start;
1403: MPI_Comm_size(PetscObjectComm((PetscObject)X),&size);
1404: if (size == 1) {
1405: VecCreateSeqHIPWithArrays(PetscObjectComm((PetscObject)X),bs,n,x,x_d,&Z);
1406: } else {
1407: VecCreateMPIHIPWithArrays(PetscObjectComm((PetscObject)X),bs,n,N,x,x_d,&Z);
1408: }
1409: Z->offloadmask = flg;
1410: #endif
1411: } else if (isstd) {
1412: PetscMPIInt size;
1414: MPI_Comm_size(PetscObjectComm((PetscObject)X),&size);
1415: VecGetArrayRead(X,&x);
1416: if (x) x += start;
1417: if (size == 1) {
1418: VecCreateSeqWithArray(PetscObjectComm((PetscObject)X),bs,n,x,&Z);
1419: } else {
1420: VecCreateMPIWithArray(PetscObjectComm((PetscObject)X),bs,n,N,x,&Z);
1421: }
1422: VecRestoreArrayRead(X,&x);
1423: } else { /* default implementation: use place array */
1424: VecGetArrayRead(X,&x);
1425: VecCreate(PetscObjectComm((PetscObject)X),&Z);
1426: VecSetType(Z,((PetscObject)X)->type_name);
1427: VecSetSizes(Z,n,N);
1428: VecSetBlockSize(Z,bs);
1429: VecPlaceArray(Z,x ? x+start : NULL);
1430: VecRestoreArrayRead(X,&x);
1431: }
1433: /* this is relevant only in debug mode */
1434: VecLockGet(X,&state);
1435: if (state) {
1436: VecLockReadPush(Z);
1437: }
1438: Z->ops->placearray = NULL;
1439: Z->ops->replacearray = NULL;
1440: } else { /* Have to create a scatter and do a copy */
1441: VecScatter scatter;
1443: VecCreate(PetscObjectComm((PetscObject)is),&Z);
1444: VecSetSizes(Z,n,N);
1445: VecSetBlockSize(Z,bs);
1446: VecSetType(Z,((PetscObject)X)->type_name);
1447: VecScatterCreate(X,is,Z,NULL,&scatter);
1448: VecScatterBegin(scatter,X,Z,INSERT_VALUES,SCATTER_FORWARD);
1449: VecScatterEnd(scatter,X,Z,INSERT_VALUES,SCATTER_FORWARD);
1450: PetscObjectCompose((PetscObject)Z,"VecGetSubVector_Scatter",(PetscObject)scatter);
1451: VecScatterDestroy(&scatter);
1452: }
1453: }
1454: /* Record the state when the subvector was gotten so we know whether its values need to be put back */
1455: if (VecGetSubVectorSavedStateId < 0) {PetscObjectComposedDataRegister(&VecGetSubVectorSavedStateId);}
1456: PetscObjectComposedDataSetInt((PetscObject)Z,VecGetSubVectorSavedStateId,1);
1457: *Y = Z;
1458: return(0);
1459: }
1461: /*@
1462: VecRestoreSubVector - Restores a subvector extracted using VecGetSubVector()
1464: Collective on IS
1466: Input Arguments:
1467: + X - vector from which subvector was obtained
1468: . is - index set representing the subset of X
1469: - Y - subvector being restored
1471: Level: advanced
1473: .seealso: VecGetSubVector()
1474: @*/
1475: PetscErrorCode VecRestoreSubVector(Vec X,IS is,Vec *Y)
1476: {
1485: if (X->ops->restoresubvector) {
1486: (*X->ops->restoresubvector)(X,is,Y);
1487: } else {
1488: PETSC_UNUSED PetscObjectState dummystate = 0;
1489: PetscBool valid;
1491: PetscObjectComposedDataGetInt((PetscObject)*Y,VecGetSubVectorSavedStateId,dummystate,valid);
1492: if (!valid) {
1493: VecScatter scatter;
1494: PetscInt state;
1496: VecLockGet(X,&state);
1497: if (state != 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Vec X is locked for read-only or read/write access");
1499: PetscObjectQuery((PetscObject)*Y,"VecGetSubVector_Scatter",(PetscObject*)&scatter);
1500: if (scatter) {
1501: VecScatterBegin(scatter,*Y,X,INSERT_VALUES,SCATTER_REVERSE);
1502: VecScatterEnd(scatter,*Y,X,INSERT_VALUES,SCATTER_REVERSE);
1503: } else {
1504: PetscBool iscuda,iship;
1505: PetscObjectTypeCompareAny((PetscObject)X,&iscuda,VECSEQCUDA,VECMPICUDA,"");
1506: PetscObjectTypeCompareAny((PetscObject)X,&iship,VECSEQHIP,VECMPIHIP,"");
1508: if (iscuda) {
1509: #if defined(PETSC_HAVE_CUDA)
1510: PetscOffloadMask ymask = (*Y)->offloadmask;
1512: /* The offloadmask of X dictates where to move memory
1513: If X GPU data is valid, then move Y data on GPU if needed
1514: Otherwise, move back to the CPU */
1515: switch (X->offloadmask) {
1516: case PETSC_OFFLOAD_BOTH:
1517: if (ymask == PETSC_OFFLOAD_CPU) {
1518: VecCUDAResetArray(*Y);
1519: } else if (ymask == PETSC_OFFLOAD_GPU) {
1520: X->offloadmask = PETSC_OFFLOAD_GPU;
1521: }
1522: break;
1523: case PETSC_OFFLOAD_GPU:
1524: if (ymask == PETSC_OFFLOAD_CPU) {
1525: VecCUDAResetArray(*Y);
1526: }
1527: break;
1528: case PETSC_OFFLOAD_CPU:
1529: if (ymask == PETSC_OFFLOAD_GPU) {
1530: VecResetArray(*Y);
1531: }
1532: break;
1533: case PETSC_OFFLOAD_UNALLOCATED:
1534: case PETSC_OFFLOAD_VECKOKKOS:
1535: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"This should not happen");
1536: break;
1537: }
1538: #endif
1539: } else if (iship) {
1540: #if defined(PETSC_HAVE_HIP)
1541: PetscOffloadMask ymask = (*Y)->offloadmask;
1543: /* The offloadmask of X dictates where to move memory
1544: If X GPU data is valid, then move Y data on GPU if needed
1545: Otherwise, move back to the CPU */
1546: switch (X->offloadmask) {
1547: case PETSC_OFFLOAD_BOTH:
1548: if (ymask == PETSC_OFFLOAD_CPU) {
1549: VecHIPResetArray(*Y);
1550: } else if (ymask == PETSC_OFFLOAD_GPU) {
1551: X->offloadmask = PETSC_OFFLOAD_GPU;
1552: }
1553: break;
1554: case PETSC_OFFLOAD_GPU:
1555: if (ymask == PETSC_OFFLOAD_CPU) {
1556: VecHIPResetArray(*Y);
1557: }
1558: break;
1559: case PETSC_OFFLOAD_CPU:
1560: if (ymask == PETSC_OFFLOAD_GPU) {
1561: VecResetArray(*Y);
1562: }
1563: break;
1564: case PETSC_OFFLOAD_UNALLOCATED:
1565: case PETSC_OFFLOAD_VECKOKKOS:
1566: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"This should not happen");
1567: break;
1568: }
1569: #endif
1570: } else {
1571: /* If OpenCL vecs updated the device memory, this triggers a copy on the CPU */
1572: VecResetArray(*Y);
1573: }
1574: PetscObjectStateIncrease((PetscObject)X);
1575: }
1576: }
1577: VecDestroy(Y);
1578: }
1579: return(0);
1580: }
1582: /*@
1583: VecGetLocalVectorRead - Maps the local portion of a vector into a
1584: vector. You must call VecRestoreLocalVectorRead() when the local
1585: vector is no longer needed.
1587: Not collective.
1589: Input parameter:
1590: . v - The vector for which the local vector is desired.
1592: Output parameter:
1593: . w - Upon exit this contains the local vector.
1595: Level: beginner
1597: Notes:
1598: This function is similar to VecGetArrayRead() which maps the local
1599: portion into a raw pointer. VecGetLocalVectorRead() is usually
1600: almost as efficient as VecGetArrayRead() but in certain circumstances
1601: VecGetLocalVectorRead() can be much more efficient than
1602: VecGetArrayRead(). This is because the construction of a contiguous
1603: array representing the vector data required by VecGetArrayRead() can
1604: be an expensive operation for certain vector types. For example, for
1605: GPU vectors VecGetArrayRead() requires that the data between device
1606: and host is synchronized.
1608: Unlike VecGetLocalVector(), this routine is not collective and
1609: preserves cached information.
1611: .seealso: VecRestoreLocalVectorRead(), VecGetLocalVector(), VecGetArrayRead(), VecGetArray()
1612: @*/
1613: PetscErrorCode VecGetLocalVectorRead(Vec v,Vec w)
1614: {
1616: PetscScalar *a;
1621: VecCheckSameLocalSize(v,1,w,2);
1622: if (v->ops->getlocalvectorread) {
1623: (*v->ops->getlocalvectorread)(v,w);
1624: } else {
1625: VecGetArrayRead(v,(const PetscScalar**)&a);
1626: VecPlaceArray(w,a);
1627: }
1628: PetscObjectStateIncrease((PetscObject)w);
1629: VecLockReadPush(v);
1630: VecLockReadPush(w);
1631: return(0);
1632: }
1634: /*@
1635: VecRestoreLocalVectorRead - Unmaps the local portion of a vector
1636: previously mapped into a vector using VecGetLocalVectorRead().
1638: Not collective.
1640: Input parameter:
1641: + v - The local portion of this vector was previously mapped into w using VecGetLocalVectorRead().
1642: - w - The vector into which the local portion of v was mapped.
1644: Level: beginner
1646: .seealso: VecGetLocalVectorRead(), VecGetLocalVector(), VecGetArrayRead(), VecGetArray()
1647: @*/
1648: PetscErrorCode VecRestoreLocalVectorRead(Vec v,Vec w)
1649: {
1651: PetscScalar *a;
1656: if (v->ops->restorelocalvectorread) {
1657: (*v->ops->restorelocalvectorread)(v,w);
1658: } else {
1659: VecGetArrayRead(w,(const PetscScalar**)&a);
1660: VecRestoreArrayRead(v,(const PetscScalar**)&a);
1661: VecResetArray(w);
1662: }
1663: VecLockReadPop(v);
1664: VecLockReadPop(w);
1665: PetscObjectStateIncrease((PetscObject)w);
1666: return(0);
1667: }
1669: /*@
1670: VecGetLocalVector - Maps the local portion of a vector into a
1671: vector.
1673: Collective on v, not collective on w.
1675: Input parameter:
1676: . v - The vector for which the local vector is desired.
1678: Output parameter:
1679: . w - Upon exit this contains the local vector.
1681: Level: beginner
1683: Notes:
1684: This function is similar to VecGetArray() which maps the local
1685: portion into a raw pointer. VecGetLocalVector() is usually about as
1686: efficient as VecGetArray() but in certain circumstances
1687: VecGetLocalVector() can be much more efficient than VecGetArray().
1688: This is because the construction of a contiguous array representing
1689: the vector data required by VecGetArray() can be an expensive
1690: operation for certain vector types. For example, for GPU vectors
1691: VecGetArray() requires that the data between device and host is
1692: synchronized.
1694: .seealso: VecRestoreLocalVector(), VecGetLocalVectorRead(), VecGetArrayRead(), VecGetArray()
1695: @*/
1696: PetscErrorCode VecGetLocalVector(Vec v,Vec w)
1697: {
1699: PetscScalar *a;
1704: VecCheckSameLocalSize(v,1,w,2);
1705: if (v->ops->getlocalvector) {
1706: (*v->ops->getlocalvector)(v,w);
1707: } else {
1708: VecGetArray(v,&a);
1709: VecPlaceArray(w,a);
1710: }
1711: PetscObjectStateIncrease((PetscObject)w);
1712: return(0);
1713: }
1715: /*@
1716: VecRestoreLocalVector - Unmaps the local portion of a vector
1717: previously mapped into a vector using VecGetLocalVector().
1719: Logically collective.
1721: Input parameter:
1722: + v - The local portion of this vector was previously mapped into w using VecGetLocalVector().
1723: - w - The vector into which the local portion of v was mapped.
1725: Level: beginner
1727: .seealso: VecGetLocalVector(), VecGetLocalVectorRead(), VecRestoreLocalVectorRead(), LocalVectorRead(), VecGetArrayRead(), VecGetArray()
1728: @*/
1729: PetscErrorCode VecRestoreLocalVector(Vec v,Vec w)
1730: {
1732: PetscScalar *a;
1737: if (v->ops->restorelocalvector) {
1738: (*v->ops->restorelocalvector)(v,w);
1739: } else {
1740: VecGetArray(w,&a);
1741: VecRestoreArray(v,&a);
1742: VecResetArray(w);
1743: }
1744: PetscObjectStateIncrease((PetscObject)w);
1745: PetscObjectStateIncrease((PetscObject)v);
1746: return(0);
1747: }
1749: /*@C
1750: VecGetArray - Returns a pointer to a contiguous array that contains this
1751: processor's portion of the vector data. For the standard PETSc
1752: vectors, VecGetArray() returns a pointer to the local data array and
1753: does not use any copies. If the underlying vector data is not stored
1754: in a contiguous array this routine will copy the data to a contiguous
1755: array and return a pointer to that. You MUST call VecRestoreArray()
1756: when you no longer need access to the array.
1758: Logically Collective on Vec
1760: Input Parameter:
1761: . x - the vector
1763: Output Parameter:
1764: . a - location to put pointer to the array
1766: Fortran Note:
1767: This routine is used differently from Fortran 77
1768: $ Vec x
1769: $ PetscScalar x_array(1)
1770: $ PetscOffset i_x
1771: $ PetscErrorCode ierr
1772: $ call VecGetArray(x,x_array,i_x,ierr)
1773: $
1774: $ Access first local entry in vector with
1775: $ value = x_array(i_x + 1)
1776: $
1777: $ ...... other code
1778: $ call VecRestoreArray(x,x_array,i_x,ierr)
1779: For Fortran 90 see VecGetArrayF90()
1781: See the Fortran chapter of the users manual and
1782: petsc/src/snes/tutorials/ex5f.F for details.
1784: Level: beginner
1786: .seealso: VecRestoreArray(), VecGetArrayRead(), VecGetArrays(), VecGetArrayF90(), VecGetArrayReadF90(), VecPlaceArray(), VecGetArray2d(),
1787: VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayWrite(), VecRestoreArrayWrite()
1788: @*/
1789: PetscErrorCode VecGetArray(Vec x,PetscScalar **a)
1790: {
1795: VecSetErrorIfLocked(x,1);
1796: if (x->ops->getarray) { /* The if-else order matters! VECNEST, VECCUDA etc should have ops->getarray while VECCUDA etc are petscnative */
1797: (*x->ops->getarray)(x,a);
1798: } else if (x->petscnative) { /* VECSTANDARD */
1799: *a = *((PetscScalar**)x->data);
1800: } else SETERRQ1(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Cannot get array for vector type \"%s\"",((PetscObject)x)->type_name);
1801: return(0);
1802: }
1804: /*@C
1805: VecRestoreArray - Restores a vector after VecGetArray() has been called.
1807: Logically Collective on Vec
1809: Input Parameters:
1810: + x - the vector
1811: - a - location of pointer to array obtained from VecGetArray()
1813: Level: beginner
1815: .seealso: VecGetArray(), VecRestoreArrayRead(), VecRestoreArrays(), VecRestoreArrayF90(), VecRestoreArrayReadF90(), VecPlaceArray(), VecRestoreArray2d(),
1816: VecGetArrayPair(), VecRestoreArrayPair()
1817: @*/
1818: PetscErrorCode VecRestoreArray(Vec x,PetscScalar **a)
1819: {
1824: if (x->ops->restorearray) { /* VECNEST, VECCUDA etc */
1825: (*x->ops->restorearray)(x,a);
1826: } else if (x->petscnative) { /* VECSTANDARD */
1827: /* nothing */
1828: } else SETERRQ1(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Cannot restore array for vector type \"%s\"",((PetscObject)x)->type_name);
1829: if (a) *a = NULL;
1830: PetscObjectStateIncrease((PetscObject)x);
1831: return(0);
1832: }
1833: /*@C
1834: VecGetArrayRead - Get read-only pointer to contiguous array containing this processor's portion of the vector data.
1836: Not Collective
1838: Input Parameter:
1839: . x - the vector
1841: Output Parameter:
1842: . a - the array
1844: Level: beginner
1846: Notes:
1847: The array must be returned using a matching call to VecRestoreArrayRead().
1849: Unlike VecGetArray(), this routine is not collective and preserves cached information like vector norms.
1851: Standard PETSc vectors use contiguous storage so that this routine does not perform a copy. Other vector
1852: implementations may require a copy, but must such implementations should cache the contiguous representation so that
1853: only one copy is performed when this routine is called multiple times in sequence.
1855: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair()
1856: @*/
1857: PetscErrorCode VecGetArrayRead(Vec x,const PetscScalar **a)
1858: {
1863: if (x->ops->getarray) { /* VECNEST, VECCUDA, VECKOKKOS etc */
1864: (*x->ops->getarray)(x,(PetscScalar**)a);
1865: } else if (x->petscnative) { /* VECSTANDARD */
1866: *a = *((PetscScalar**)x->data);
1867: } else SETERRQ1(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Cannot get array read for vector type \"%s\"",((PetscObject)x)->type_name);
1868: return(0);
1869: }
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(), VecGetArrayPair(), VecRestoreArrayPair()
1883: @*/
1884: PetscErrorCode VecRestoreArrayRead(Vec x,const PetscScalar **a)
1885: {
1890: if (x->petscnative) { /* VECSTANDARD, VECCUDA, VECKOKKOS etc */
1891: /* nothing */
1892: } else if (x->ops->restorearrayread) { /* VECNEST */
1893: (*x->ops->restorearrayread)(x,a);
1894: } else { /* No one? */
1895: (*x->ops->restorearray)(x,(PetscScalar**)a);
1896: }
1897: if (a) *a = NULL;
1898: return(0);
1899: }
1901: /*@C
1902: VecGetArrayWrite - Returns a pointer to a contiguous array that WILL contains this
1903: processor's portion of the vector data. The values in this array are NOT valid, the routine calling this
1904: routine is responsible for putting values into the array; any values it does not set will be invalid
1906: Logically Collective on Vec
1908: Input Parameter:
1909: . x - the vector
1911: Output Parameter:
1912: . a - location to put pointer to the array
1914: Level: intermediate
1916: This is for vectors associate with GPUs, the vector is not copied up before giving access. If you need correct
1917: values in the array use VecGetArray()
1919: Concepts: vector^accessing local values
1921: .seealso: VecRestoreArray(), VecGetArrayRead(), VecGetArrays(), VecGetArrayF90(), VecGetArrayReadF90(), VecPlaceArray(), VecGetArray2d(),
1922: VecGetArrayPair(), VecRestoreArrayPair(), VecGetArray(), VecRestoreArrayWrite()
1923: @*/
1924: PetscErrorCode VecGetArrayWrite(Vec x,PetscScalar **a)
1925: {
1930: VecSetErrorIfLocked(x,1);
1931: if (x->ops->getarraywrite) {
1932: (*x->ops->getarraywrite)(x,a);
1933: } else {
1934: VecGetArray(x,a);
1935: }
1936: return(0);
1937: }
1939: /*@C
1940: VecRestoreArrayWrite - Restores a vector after VecGetArrayWrite() has been called.
1942: Logically Collective on Vec
1944: Input Parameters:
1945: + x - the vector
1946: - a - location of pointer to array obtained from VecGetArray()
1948: Level: beginner
1950: .seealso: VecGetArray(), VecRestoreArrayRead(), VecRestoreArrays(), VecRestoreArrayF90(), VecRestoreArrayReadF90(), VecPlaceArray(), VecRestoreArray2d(),
1951: VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayWrite()
1952: @*/
1953: PetscErrorCode VecRestoreArrayWrite(Vec x,PetscScalar **a)
1954: {
1959: if (x->ops->restorearraywrite) {
1960: (*x->ops->restorearraywrite)(x,a);
1961: } else if (x->ops->restorearray) {
1962: (*x->ops->restorearray)(x,a);
1963: }
1964: if (a) *a = NULL;
1965: PetscObjectStateIncrease((PetscObject)x);
1966: return(0);
1967: }
1969: /*@C
1970: VecGetArrays - Returns a pointer to the arrays in a set of vectors
1971: that were created by a call to VecDuplicateVecs(). You MUST call
1972: VecRestoreArrays() when you no longer need access to the array.
1974: Logically Collective on Vec
1976: Input Parameters:
1977: + x - the vectors
1978: - n - the number of vectors
1980: Output Parameter:
1981: . a - location to put pointer to the array
1983: Fortran Note:
1984: This routine is not supported in Fortran.
1986: Level: intermediate
1988: .seealso: VecGetArray(), VecRestoreArrays()
1989: @*/
1990: PetscErrorCode VecGetArrays(const Vec x[],PetscInt n,PetscScalar **a[])
1991: {
1993: PetscInt i;
1994: PetscScalar **q;
2000: if (n <= 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must get at least one array n = %D",n);
2001: PetscMalloc1(n,&q);
2002: for (i=0; i<n; ++i) {
2003: VecGetArray(x[i],&q[i]);
2004: }
2005: *a = q;
2006: return(0);
2007: }
2009: /*@C
2010: VecRestoreArrays - Restores a group of vectors after VecGetArrays()
2011: has been called.
2013: Logically Collective on Vec
2015: Input Parameters:
2016: + x - the vector
2017: . n - the number of vectors
2018: - a - location of pointer to arrays obtained from VecGetArrays()
2020: Notes:
2021: For regular PETSc vectors this routine does not involve any copies. For
2022: any special vectors that do not store local vector data in a contiguous
2023: array, this routine will copy the data back into the underlying
2024: vector data structure from the arrays obtained with VecGetArrays().
2026: Fortran Note:
2027: This routine is not supported in Fortran.
2029: Level: intermediate
2031: .seealso: VecGetArrays(), VecRestoreArray()
2032: @*/
2033: PetscErrorCode VecRestoreArrays(const Vec x[],PetscInt n,PetscScalar **a[])
2034: {
2036: PetscInt i;
2037: PetscScalar **q = *a;
2044: for (i=0; i<n; ++i) {
2045: VecRestoreArray(x[i],&q[i]);
2046: }
2047: PetscFree(q);
2048: return(0);
2049: }
2051: /*@C
2052: VecGetArrayAndMemType - Like VecGetArray(), but if this is a device vector (e.g., VECCUDA) and the device has up-to-date data,
2053: the returned pointer will be a device pointer to the device memory that contains this processor's portion of the vector data.
2054: Otherwise, when this is a host vector (e.g., VECMPI), or a device vector with the host having newer data than device,
2055: it functions as VecGetArray() and returns a host pointer.
2057: Logically Collective on Vec
2059: Input Parameter:
2060: . x - the vector
2062: Output Parameters:
2063: + a - location to put pointer to the array
2064: - mtype - memory type of the array
2066: Level: beginner
2068: .seealso: VecRestoreArrayAndMemType(), VecRestoreArrayAndMemType(), VecRestoreArray(), VecGetArrayRead(), VecGetArrays(), VecGetArrayF90(), VecGetArrayReadF90(),
2069: VecPlaceArray(), VecGetArray2d(), VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayWrite(), VecRestoreArrayWrite()
2070: @*/
2071: PetscErrorCode VecGetArrayAndMemType(Vec x,PetscScalar **a,PetscMemType *mtype)
2072: {
2078: VecSetErrorIfLocked(x,1);
2079: if (x->ops->getarrayandmemtype) { /* VECCUDA, VECKOKKOS etc */
2080: (*x->ops->getarrayandmemtype)(x,a,mtype);
2081: } else { /* VECSTANDARD, VECNEST, VECVIENNACL */
2082: VecGetArray(x,a);
2083: if (mtype) *mtype = PETSC_MEMTYPE_HOST;
2084: }
2085: return(0);
2086: }
2088: /*@C
2089: VecRestoreArrayAndMemType - Restores a vector after VecGetArrayAndMemType() has been called.
2091: Logically Collective on Vec
2093: Input Parameters:
2094: + x - the vector
2095: - a - location of pointer to array obtained from VecGetArrayAndMemType()
2097: Level: beginner
2099: .seealso: VecGetArrayAndMemType(), VecGetArray(), VecRestoreArrayRead(), VecRestoreArrays(), VecRestoreArrayF90(), VecRestoreArrayReadF90(),
2100: VecPlaceArray(), VecRestoreArray2d(), VecGetArrayPair(), VecRestoreArrayPair()
2101: @*/
2102: PetscErrorCode VecRestoreArrayAndMemType(Vec x,PetscScalar **a)
2103: {
2109: if (x->ops->restorearrayandmemtype) { /* VECCUDA, VECKOKKOS etc */
2110: (*x->ops->restorearrayandmemtype)(x,a);
2111: } else if (x->ops->restorearray) { /* VECNEST, VECVIENNACL */
2112: (*x->ops->restorearray)(x,a);
2113: } /* VECSTANDARD does nothing */
2114: if (a) *a = NULL;
2115: PetscObjectStateIncrease((PetscObject)x);
2116: return(0);
2117: }
2119: /*@C
2120: VecGetArrayReadAndMemType - Like VecGetArrayRead(), but if this is a CUDA vector and it is currently offloaded to GPU,
2121: the returned pointer will be a GPU pointer to the GPU memory that contains this processor's portion of the
2122: vector data. Otherwise, it functions as VecGetArrayRead().
2124: Not Collective
2126: Input Parameter:
2127: . x - the vector
2129: Output Parameters:
2130: + a - the array
2131: - mtype - memory type of the array
2133: Level: beginner
2135: Notes:
2136: The array must be returned using a matching call to VecRestoreArrayReadAndMemType().
2139: .seealso: VecRestoreArrayReadAndMemType(), VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayAndMemType()
2140: @*/
2141: PetscErrorCode VecGetArrayReadAndMemType(Vec x,const PetscScalar **a,PetscMemType *mtype)
2142: {
2148: #if defined(PETSC_USE_DEBUG)
2149: if (x->ops->getarrayreadandmemtype) SETERRQ1(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Not expected vector type \"%s\" has ops->getarrayreadandmemtype",((PetscObject)x)->type_name);
2150: #endif
2152: if (x->ops->getarrayandmemtype) { /* VECCUDA, VECKOKKOS etc, though they are also petscnative */
2153: (*x->ops->getarrayandmemtype)(x,(PetscScalar**)a,mtype);
2154: } else if (x->ops->getarray) { /* VECNEST, VECVIENNACL */
2155: (*x->ops->getarray)(x,(PetscScalar**)a);
2156: if (mtype) *mtype = PETSC_MEMTYPE_HOST;
2157: } else if (x->petscnative) { /* VECSTANDARD */
2158: *a = *((PetscScalar**)x->data);
2159: if (mtype) *mtype = PETSC_MEMTYPE_HOST;
2160: } else SETERRQ1(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Cannot get array read in place for vector type \"%s\"",((PetscObject)x)->type_name);
2161: return(0);
2162: }
2164: /*@C
2165: VecRestoreArrayReadAndMemType - Restore array obtained with VecGetArrayReadAndMemType()
2167: Not Collective
2169: Input Parameters:
2170: + vec - the vector
2171: - array - the array
2173: Level: beginner
2175: .seealso: VecGetArrayReadAndMemType(), VecRestoreArrayAndMemType(), VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair()
2176: @*/
2177: PetscErrorCode VecRestoreArrayReadAndMemType(Vec x,const PetscScalar **a)
2178: {
2184: if (x->petscnative) { /* VECSTANDARD, VECCUDA, VECKOKKOS, VECVIENNACL etc */
2185: /* nothing */
2186: } else if (x->ops->restorearrayread) { /* VECNEST */
2187: (*x->ops->restorearrayread)(x,a);
2188: } else SETERRQ1(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Cannot restore array read in place for vector type \"%s\"",((PetscObject)x)->type_name);
2189: if (a) *a = NULL;
2190: return(0);
2191: }
2193: /*@
2194: VecPlaceArray - Allows one to replace the array in a vector with an
2195: array provided by the user. This is useful to avoid copying an array
2196: into a vector.
2198: Not Collective
2200: Input Parameters:
2201: + vec - the vector
2202: - array - the array
2204: Notes:
2205: You can return to the original array with a call to VecResetArray()
2207: Level: developer
2209: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecResetArray()
2211: @*/
2212: PetscErrorCode VecPlaceArray(Vec vec,const PetscScalar array[])
2213: {
2220: if (vec->ops->placearray) {
2221: (*vec->ops->placearray)(vec,array);
2222: } else SETERRQ(PetscObjectComm((PetscObject)vec),PETSC_ERR_SUP,"Cannot place array in this type of vector");
2223: PetscObjectStateIncrease((PetscObject)vec);
2224: return(0);
2225: }
2227: /*@C
2228: VecReplaceArray - Allows one to replace the array in a vector with an
2229: array provided by the user. This is useful to avoid copying an array
2230: into a vector.
2232: Not Collective
2234: Input Parameters:
2235: + vec - the vector
2236: - array - the array
2238: Notes:
2239: This permanently replaces the array and frees the memory associated
2240: with the old array.
2242: The memory passed in MUST be obtained with PetscMalloc() and CANNOT be
2243: freed by the user. It will be freed when the vector is destroyed.
2245: Not supported from Fortran
2247: Level: developer
2249: .seealso: VecGetArray(), VecRestoreArray(), VecPlaceArray(), VecResetArray()
2251: @*/
2252: PetscErrorCode VecReplaceArray(Vec vec,const PetscScalar array[])
2253: {
2259: if (vec->ops->replacearray) {
2260: (*vec->ops->replacearray)(vec,array);
2261: } else SETERRQ(PetscObjectComm((PetscObject)vec),PETSC_ERR_SUP,"Cannot replace array in this type of vector");
2262: PetscObjectStateIncrease((PetscObject)vec);
2263: return(0);
2264: }
2267: /*@C
2268: VecCUDAGetArray - Provides access to the CUDA buffer inside a vector.
2270: This function has semantics similar to VecGetArray(): the pointer
2271: returned by this function points to a consistent view of the vector
2272: data. This may involve a copy operation of data from the host to the
2273: device if the data on the device is out of date. If the device
2274: memory hasn't been allocated previously it will be allocated as part
2275: of this function call. VecCUDAGetArray() assumes that
2276: the user will modify the vector data. This is similar to
2277: intent(inout) in fortran.
2279: The CUDA device pointer has to be released by calling
2280: VecCUDARestoreArray(). Upon restoring the vector data
2281: the data on the host will be marked as out of date. A subsequent
2282: access of the host data will thus incur a data transfer from the
2283: device to the host.
2286: Input Parameter:
2287: . v - the vector
2289: Output Parameter:
2290: . a - the CUDA device pointer
2292: Fortran note:
2293: This function is not currently available from Fortran.
2295: Level: intermediate
2297: .seealso: VecCUDARestoreArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2298: @*/
2299: PETSC_EXTERN PetscErrorCode VecCUDAGetArray(Vec v, PetscScalar **a)
2300: {
2303: #if defined(PETSC_HAVE_CUDA)
2304: {
2306: VecCUDACopyToGPU(v);
2307: *a = ((Vec_CUDA*)v->spptr)->GPUarray;
2308: }
2309: #endif
2310: return(0);
2311: }
2313: /*@C
2314: VecCUDARestoreArray - Restore a CUDA device pointer previously acquired with VecCUDAGetArray().
2316: This marks the host data as out of date. Subsequent access to the
2317: vector data on the host side with for instance VecGetArray() incurs a
2318: data transfer.
2320: Input Parameter:
2321: + v - the vector
2322: - a - the CUDA device pointer. This pointer is invalid after
2323: VecCUDARestoreArray() returns.
2325: Fortran note:
2326: This function is not currently available from Fortran.
2328: Level: intermediate
2330: .seealso: VecCUDAGetArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2331: @*/
2332: PETSC_EXTERN PetscErrorCode VecCUDARestoreArray(Vec v, PetscScalar **a)
2333: {
2338: #if defined(PETSC_HAVE_CUDA)
2339: v->offloadmask = PETSC_OFFLOAD_GPU;
2340: #endif
2341: PetscObjectStateIncrease((PetscObject)v);
2342: return(0);
2343: }
2345: /*@C
2346: VecCUDAGetArrayRead - Provides read access to the CUDA buffer inside a vector.
2348: This function is analogous to VecGetArrayRead(): The pointer
2349: returned by this function points to a consistent view of the vector
2350: data. This may involve a copy operation of data from the host to the
2351: device if the data on the device is out of date. If the device
2352: memory hasn't been allocated previously it will be allocated as part
2353: of this function call. VecCUDAGetArrayRead() assumes that the
2354: user will not modify the vector data. This is analgogous to
2355: intent(in) in Fortran.
2357: The CUDA device pointer has to be released by calling
2358: VecCUDARestoreArrayRead(). If the data on the host side was
2359: previously up to date it will remain so, i.e. data on both the device
2360: and the host is up to date. Accessing data on the host side does not
2361: incur a device to host data transfer.
2363: Input Parameter:
2364: . v - the vector
2366: Output Parameter:
2367: . a - the CUDA pointer.
2369: Fortran note:
2370: This function is not currently available from Fortran.
2372: Level: intermediate
2374: .seealso: VecCUDARestoreArrayRead(), VecCUDAGetArray(), VecCUDAGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2375: @*/
2376: PETSC_EXTERN PetscErrorCode VecCUDAGetArrayRead(Vec v,const PetscScalar** a)
2377: {
2380: VecCUDAGetArray(v,(PetscScalar**)a);
2381: return(0);
2382: }
2384: /*@C
2385: VecCUDARestoreArrayRead - Restore a CUDA device pointer previously acquired with VecCUDAGetArrayRead().
2387: If the data on the host side was previously up to date it will remain
2388: so, i.e. data on both the device and the host is up to date.
2389: Accessing data on the host side e.g. with VecGetArray() does not
2390: incur a device to host data transfer.
2392: Input Parameter:
2393: + v - the vector
2394: - a - the CUDA device pointer. This pointer is invalid after
2395: VecCUDARestoreArrayRead() returns.
2397: Fortran note:
2398: This function is not currently available from Fortran.
2400: Level: intermediate
2402: .seealso: VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecCUDAGetArray(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2403: @*/
2404: PETSC_EXTERN PetscErrorCode VecCUDARestoreArrayRead(Vec v, const PetscScalar **a)
2405: {
2408: *a = NULL;
2409: return(0);
2410: }
2412: /*@C
2413: VecCUDAGetArrayWrite - Provides write access to the CUDA buffer inside a vector.
2415: The data pointed to by the device pointer is uninitialized. The user
2416: may not read from this data. Furthermore, the entire array needs to
2417: be filled by the user to obtain well-defined behaviour. The device
2418: memory will be allocated by this function if it hasn't been allocated
2419: previously. This is analogous to intent(out) in Fortran.
2421: The device pointer needs to be released with
2422: VecCUDARestoreArrayWrite(). When the pointer is released the
2423: host data of the vector is marked as out of data. Subsequent access
2424: of the host data with e.g. VecGetArray() incurs a device to host data
2425: transfer.
2428: Input Parameter:
2429: . v - the vector
2431: Output Parameter:
2432: . a - the CUDA pointer
2434: Fortran note:
2435: This function is not currently available from Fortran.
2437: Level: advanced
2439: .seealso: VecCUDARestoreArrayWrite(), VecCUDAGetArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2440: @*/
2441: PETSC_EXTERN PetscErrorCode VecCUDAGetArrayWrite(Vec v, PetscScalar **a)
2442: {
2445: #if defined(PETSC_HAVE_CUDA)
2446: {
2448: VecCUDAAllocateCheck(v);
2449: *a = ((Vec_CUDA*)v->spptr)->GPUarray;
2450: }
2451: #endif
2452: return(0);
2453: }
2455: /*@C
2456: VecCUDARestoreArrayWrite - Restore a CUDA device pointer previously acquired with VecCUDAGetArrayWrite().
2458: Data on the host will be marked as out of date. Subsequent access of
2459: the data on the host side e.g. with VecGetArray() will incur a device
2460: to host data transfer.
2462: Input Parameter:
2463: + v - the vector
2464: - a - the CUDA device pointer. This pointer is invalid after
2465: VecCUDARestoreArrayWrite() returns.
2467: Fortran note:
2468: This function is not currently available from Fortran.
2470: Level: intermediate
2472: .seealso: VecCUDAGetArrayWrite(), VecCUDAGetArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2473: @*/
2474: PETSC_EXTERN PetscErrorCode VecCUDARestoreArrayWrite(Vec v, PetscScalar **a)
2475: {
2480: #if defined(PETSC_HAVE_CUDA)
2481: v->offloadmask = PETSC_OFFLOAD_GPU;
2482: if (a) *a = NULL;
2483: #endif
2484: PetscObjectStateIncrease((PetscObject)v);
2485: return(0);
2486: }
2488: /*@C
2489: VecCUDAPlaceArray - Allows one to replace the GPU array in a vector with a
2490: GPU array provided by the user. This is useful to avoid copying an
2491: array into a vector.
2493: Not Collective
2495: Input Parameters:
2496: + vec - the vector
2497: - array - the GPU array
2499: Notes:
2500: You can return to the original GPU array with a call to VecCUDAResetArray()
2501: It is not possible to use VecCUDAPlaceArray() and VecPlaceArray() at the
2502: same time on the same vector.
2504: Level: developer
2506: .seealso: VecPlaceArray(), VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecResetArray(), VecCUDAResetArray(), VecCUDAReplaceArray()
2508: @*/
2509: PetscErrorCode VecCUDAPlaceArray(Vec vin,const PetscScalar a[])
2510: {
2515: #if defined(PETSC_HAVE_CUDA)
2516: VecCUDACopyToGPU(vin);
2517: if (((Vec_Seq*)vin->data)->unplacedarray) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"VecCUDAPlaceArray()/VecPlaceArray() was already called on this vector, without a call to VecCUDAResetArray()/VecResetArray()");
2518: ((Vec_Seq*)vin->data)->unplacedarray = (PetscScalar *) ((Vec_CUDA*)vin->spptr)->GPUarray; /* save previous GPU array so reset can bring it back */
2519: ((Vec_CUDA*)vin->spptr)->GPUarray = (PetscScalar*)a;
2520: vin->offloadmask = PETSC_OFFLOAD_GPU;
2521: #endif
2522: PetscObjectStateIncrease((PetscObject)vin);
2523: return(0);
2524: }
2526: /*@C
2527: VecCUDAReplaceArray - Allows one to replace the GPU array in a vector
2528: with a GPU array provided by the user. This is useful to avoid copying
2529: a GPU array into a vector.
2531: Not Collective
2533: Input Parameters:
2534: + vec - the vector
2535: - array - the GPU array
2537: Notes:
2538: This permanently replaces the GPU array and frees the memory associated
2539: with the old GPU array.
2541: The memory passed in CANNOT be freed by the user. It will be freed
2542: when the vector is destroyed.
2544: Not supported from Fortran
2546: Level: developer
2548: .seealso: VecGetArray(), VecRestoreArray(), VecPlaceArray(), VecResetArray(), VecCUDAResetArray(), VecCUDAPlaceArray(), VecReplaceArray()
2550: @*/
2551: PetscErrorCode VecCUDAReplaceArray(Vec vin,const PetscScalar a[])
2552: {
2553: #if defined(PETSC_HAVE_CUDA)
2554: cudaError_t err;
2555: #endif
2560: #if defined(PETSC_HAVE_CUDA)
2561: if (((Vec_CUDA*)vin->spptr)->GPUarray_allocated) {
2562: err = cudaFree(((Vec_CUDA*)vin->spptr)->GPUarray_allocated);CHKERRCUDA(err);
2563: }
2564: ((Vec_CUDA*)vin->spptr)->GPUarray_allocated = ((Vec_CUDA*)vin->spptr)->GPUarray = (PetscScalar*)a;
2565: vin->offloadmask = PETSC_OFFLOAD_GPU;
2566: #endif
2567: PetscObjectStateIncrease((PetscObject)vin);
2568: return(0);
2569: }
2571: /*@C
2572: VecCUDAResetArray - Resets a vector to use its default memory. Call this
2573: after the use of VecCUDAPlaceArray().
2575: Not Collective
2577: Input Parameters:
2578: . vec - the vector
2580: Level: developer
2582: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecPlaceArray(), VecResetArray(), VecCUDAPlaceArray(), VecCUDAReplaceArray()
2584: @*/
2585: PetscErrorCode VecCUDAResetArray(Vec vin)
2586: {
2591: #if defined(PETSC_HAVE_CUDA)
2592: VecCUDACopyToGPU(vin);
2593: ((Vec_CUDA*)vin->spptr)->GPUarray = (PetscScalar *) ((Vec_Seq*)vin->data)->unplacedarray;
2594: ((Vec_Seq*)vin->data)->unplacedarray = 0;
2595: vin->offloadmask = PETSC_OFFLOAD_GPU;
2596: #endif
2597: PetscObjectStateIncrease((PetscObject)vin);
2598: return(0);
2599: }
2601: /*@C
2602: VecHIPGetArray - Provides access to the HIP buffer inside a vector.
2604: This function has semantics similar to VecGetArray(): the pointer
2605: returned by this function points to a consistent view of the vector
2606: data. This may involve a copy operation of data from the host to the
2607: device if the data on the device is out of date. If the device
2608: memory hasn't been allocated previously it will be allocated as part
2609: of this function call. VecHIPGetArray() assumes that
2610: the user will modify the vector data. This is similar to
2611: intent(inout) in fortran.
2613: The HIP device pointer has to be released by calling
2614: VecHIPRestoreArray(). Upon restoring the vector data
2615: the data on the host will be marked as out of date. A subsequent
2616: access of the host data will thus incur a data transfer from the
2617: device to the host.
2620: Input Parameter:
2621: . v - the vector
2623: Output Parameter:
2624: . a - the HIP device pointer
2626: Fortran note:
2627: This function is not currently available from Fortran.
2629: Level: intermediate
2631: .seealso: VecHIPRestoreArray(), VecHIPGetArrayRead(), VecHIPGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2632: @*/
2633: PETSC_EXTERN PetscErrorCode VecHIPGetArray(Vec v, PetscScalar **a)
2634: {
2635: #if defined(PETSC_HAVE_HIP)
2637: #endif
2641: #if defined(PETSC_HAVE_HIP)
2642: *a = 0;
2643: VecHIPCopyToGPU(v);
2644: *a = ((Vec_HIP*)v->spptr)->GPUarray;
2645: #endif
2646: return(0);
2647: }
2649: /*@C
2650: VecHIPRestoreArray - Restore a HIP device pointer previously acquired with VecHIPGetArray().
2652: This marks the host data as out of date. Subsequent access to the
2653: vector data on the host side with for instance VecGetArray() incurs a
2654: data transfer.
2656: Input Parameter:
2657: + v - the vector
2658: - a - the HIP device pointer. This pointer is invalid after
2659: VecHIPRestoreArray() returns.
2661: Fortran note:
2662: This function is not currently available from Fortran.
2664: Level: intermediate
2666: .seealso: VecHIPGetArray(), VecHIPGetArrayRead(), VecHIPGetArrayWrite(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2667: @*/
2668: PETSC_EXTERN PetscErrorCode VecHIPRestoreArray(Vec v, PetscScalar **a)
2669: {
2674: #if defined(PETSC_HAVE_HIP)
2675: v->offloadmask = PETSC_OFFLOAD_GPU;
2676: #endif
2678: PetscObjectStateIncrease((PetscObject)v);
2679: return(0);
2680: }
2682: /*@C
2683: VecHIPGetArrayRead - Provides read access to the HIP buffer inside a vector.
2685: This function is analogous to VecGetArrayRead(): The pointer
2686: returned by this function points to a consistent view of the vector
2687: data. This may involve a copy operation of data from the host to the
2688: device if the data on the device is out of date. If the device
2689: memory hasn't been allocated previously it will be allocated as part
2690: of this function call. VecHIPGetArrayRead() assumes that the
2691: user will not modify the vector data. This is analgogous to
2692: intent(in) in Fortran.
2694: The HIP device pointer has to be released by calling
2695: VecHIPRestoreArrayRead(). If the data on the host side was
2696: previously up to date it will remain so, i.e. data on both the device
2697: and the host is up to date. Accessing data on the host side does not
2698: incur a device to host data transfer.
2700: Input Parameter:
2701: . v - the vector
2703: Output Parameter:
2704: . a - the HIP pointer.
2706: Fortran note:
2707: This function is not currently available from Fortran.
2709: Level: intermediate
2711: .seealso: VecHIPRestoreArrayRead(), VecHIPGetArray(), VecHIPGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2712: @*/
2713: PETSC_EXTERN PetscErrorCode VecHIPGetArrayRead(Vec v, const PetscScalar **a)
2714: {
2715: #if defined(PETSC_HAVE_HIP)
2717: #endif
2721: #if defined(PETSC_HAVE_HIP)
2722: *a = 0;
2723: VecHIPCopyToGPU(v);
2724: *a = ((Vec_HIP*)v->spptr)->GPUarray;
2725: #endif
2726: return(0);
2727: }
2729: /*@C
2730: VecHIPRestoreArrayRead - Restore a HIP device pointer previously acquired with VecHIPGetArrayRead().
2732: If the data on the host side was previously up to date it will remain
2733: so, i.e. data on both the device and the host is up to date.
2734: Accessing data on the host side e.g. with VecGetArray() does not
2735: incur a device to host data transfer.
2737: Input Parameter:
2738: + v - the vector
2739: - a - the HIP device pointer. This pointer is invalid after
2740: VecHIPRestoreArrayRead() returns.
2742: Fortran note:
2743: This function is not currently available from Fortran.
2745: Level: intermediate
2747: .seealso: VecHIPGetArrayRead(), VecHIPGetArrayWrite(), VecHIPGetArray(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2748: @*/
2749: PETSC_EXTERN PetscErrorCode VecHIPRestoreArrayRead(Vec v, const PetscScalar **a)
2750: {
2753: *a = NULL;
2754: return(0);
2755: }
2757: /*@C
2758: VecHIPGetArrayWrite - Provides write access to the HIP buffer inside a vector.
2760: The data pointed to by the device pointer is uninitialized. The user
2761: may not read from this data. Furthermore, the entire array needs to
2762: be filled by the user to obtain well-defined behaviour. The device
2763: memory will be allocated by this function if it hasn't been allocated
2764: previously. This is analogous to intent(out) in Fortran.
2766: The device pointer needs to be released with
2767: VecHIPRestoreArrayWrite(). When the pointer is released the
2768: host data of the vector is marked as out of data. Subsequent access
2769: of the host data with e.g. VecGetArray() incurs a device to host data
2770: transfer.
2773: Input Parameter:
2774: . v - the vector
2776: Output Parameter:
2777: . a - the HIP pointer
2779: Fortran note:
2780: This function is not currently available from Fortran.
2782: Level: advanced
2784: .seealso: VecHIPRestoreArrayWrite(), VecHIPGetArray(), VecHIPGetArrayRead(), VecHIPGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2785: @*/
2786: PETSC_EXTERN PetscErrorCode VecHIPGetArrayWrite(Vec v, PetscScalar **a)
2787: {
2788: #if defined(PETSC_HAVE_HIP)
2790: #endif
2794: #if defined(PETSC_HAVE_HIP)
2795: *a = 0;
2796: VecHIPAllocateCheck(v);
2797: *a = ((Vec_HIP*)v->spptr)->GPUarray;
2798: #endif
2799: return(0);
2800: }
2802: /*@C
2803: VecHIPRestoreArrayWrite - Restore a HIP device pointer previously acquired with VecHIPGetArrayWrite().
2805: Data on the host will be marked as out of date. Subsequent access of
2806: the data on the host side e.g. with VecGetArray() will incur a device
2807: to host data transfer.
2809: Input Parameter:
2810: + v - the vector
2811: - a - the HIP device pointer. This pointer is invalid after
2812: VecHIPRestoreArrayWrite() returns.
2814: Fortran note:
2815: This function is not currently available from Fortran.
2817: Level: intermediate
2819: .seealso: VecHIPGetArrayWrite(), VecHIPGetArray(), VecHIPGetArrayRead(), VecHIPGetArrayWrite(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2820: @*/
2821: PETSC_EXTERN PetscErrorCode VecHIPRestoreArrayWrite(Vec v, PetscScalar **a)
2822: {
2827: #if defined(PETSC_HAVE_HIP)
2828: v->offloadmask = PETSC_OFFLOAD_GPU;
2829: #endif
2831: PetscObjectStateIncrease((PetscObject)v);
2832: return(0);
2833: }
2835: /*@C
2836: VecHIPPlaceArray - Allows one to replace the GPU array in a vector with a
2837: GPU array provided by the user. This is useful to avoid copying an
2838: array into a vector.
2840: Not Collective
2842: Input Parameters:
2843: + vec - the vector
2844: - array - the GPU array
2846: Notes:
2847: You can return to the original GPU array with a call to VecHIPResetArray()
2848: It is not possible to use VecHIPPlaceArray() and VecPlaceArray() at the
2849: same time on the same vector.
2851: Level: developer
2853: .seealso: VecPlaceArray(), VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecResetArray(), VecHIPResetArray(), VecHIPReplaceArray()
2855: @*/
2856: PetscErrorCode VecHIPPlaceArray(Vec vin,const PetscScalar a[])
2857: {
2862: #if defined(PETSC_HAVE_HIP)
2863: VecHIPCopyToGPU(vin);
2864: if (((Vec_Seq*)vin->data)->unplacedarray) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"VecHIPPlaceArray()/VecPlaceArray() was already called on this vector, without a call to VecHIPResetArray()/VecResetArray()");
2865: ((Vec_Seq*)vin->data)->unplacedarray = (PetscScalar *) ((Vec_HIP*)vin->spptr)->GPUarray; /* save previous GPU array so reset can bring it back */
2866: ((Vec_HIP*)vin->spptr)->GPUarray = (PetscScalar*)a;
2867: vin->offloadmask = PETSC_OFFLOAD_GPU;
2868: #endif
2869: PetscObjectStateIncrease((PetscObject)vin);
2870: return(0);
2871: }
2873: /*@C
2874: VecHIPReplaceArray - Allows one to replace the GPU array in a vector
2875: with a GPU array provided by the user. This is useful to avoid copying
2876: a GPU array into a vector.
2878: Not Collective
2880: Input Parameters:
2881: + vec - the vector
2882: - array - the GPU array
2884: Notes:
2885: This permanently replaces the GPU array and frees the memory associated
2886: with the old GPU array.
2888: The memory passed in CANNOT be freed by the user. It will be freed
2889: when the vector is destroyed.
2891: Not supported from Fortran
2893: Level: developer
2895: .seealso: VecGetArray(), VecRestoreArray(), VecPlaceArray(), VecResetArray(), VecHIPResetArray(), VecHIPPlaceArray(), VecReplaceArray()
2897: @*/
2898: PetscErrorCode VecHIPReplaceArray(Vec vin,const PetscScalar a[])
2899: {
2900: #if defined(PETSC_HAVE_HIP)
2901: hipError_t err;
2902: #endif
2907: #if defined(PETSC_HAVE_HIP)
2908: err = hipFree(((Vec_HIP*)vin->spptr)->GPUarray);CHKERRHIP(err);
2909: ((Vec_HIP*)vin->spptr)->GPUarray = (PetscScalar*)a;
2910: vin->offloadmask = PETSC_OFFLOAD_GPU;
2911: #endif
2912: PetscObjectStateIncrease((PetscObject)vin);
2913: return(0);
2914: }
2916: /*@C
2917: VecHIPResetArray - Resets a vector to use its default memory. Call this
2918: after the use of VecHIPPlaceArray().
2920: Not Collective
2922: Input Parameters:
2923: . vec - the vector
2925: Level: developer
2927: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecPlaceArray(), VecResetArray(), VecHIPPlaceArray(), VecHIPReplaceArray()
2929: @*/
2930: PetscErrorCode VecHIPResetArray(Vec vin)
2931: {
2936: #if defined(PETSC_HAVE_HIP)
2937: VecHIPCopyToGPU(vin);
2938: ((Vec_HIP*)vin->spptr)->GPUarray = (PetscScalar *) ((Vec_Seq*)vin->data)->unplacedarray;
2939: ((Vec_Seq*)vin->data)->unplacedarray = 0;
2940: vin->offloadmask = PETSC_OFFLOAD_GPU;
2941: #endif
2942: PetscObjectStateIncrease((PetscObject)vin);
2943: return(0);
2944: }
2946: /*MC
2947: VecDuplicateVecsF90 - Creates several vectors of the same type as an existing vector
2948: and makes them accessible via a Fortran90 pointer.
2950: Synopsis:
2951: VecDuplicateVecsF90(Vec x,PetscInt n,{Vec, pointer :: y(:)},integer ierr)
2953: Collective on Vec
2955: Input Parameters:
2956: + x - a vector to mimic
2957: - n - the number of vectors to obtain
2959: Output Parameters:
2960: + y - Fortran90 pointer to the array of vectors
2961: - ierr - error code
2963: Example of Usage:
2964: .vb
2965: #include <petsc/finclude/petscvec.h>
2966: use petscvec
2968: Vec x
2969: Vec, pointer :: y(:)
2970: ....
2971: call VecDuplicateVecsF90(x,2,y,ierr)
2972: call VecSet(y(2),alpha,ierr)
2973: call VecSet(y(2),alpha,ierr)
2974: ....
2975: call VecDestroyVecsF90(2,y,ierr)
2976: .ve
2978: Notes:
2979: Not yet supported for all F90 compilers
2981: Use VecDestroyVecsF90() to free the space.
2983: Level: beginner
2985: .seealso: VecDestroyVecsF90(), VecDuplicateVecs()
2987: M*/
2989: /*MC
2990: VecRestoreArrayF90 - Restores a vector to a usable state after a call to
2991: VecGetArrayF90().
2993: Synopsis:
2994: VecRestoreArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2996: Logically Collective on Vec
2998: Input Parameters:
2999: + x - vector
3000: - xx_v - the Fortran90 pointer to the array
3002: Output Parameter:
3003: . ierr - error code
3005: Example of Usage:
3006: .vb
3007: #include <petsc/finclude/petscvec.h>
3008: use petscvec
3010: PetscScalar, pointer :: xx_v(:)
3011: ....
3012: call VecGetArrayF90(x,xx_v,ierr)
3013: xx_v(3) = a
3014: call VecRestoreArrayF90(x,xx_v,ierr)
3015: .ve
3017: Level: beginner
3019: .seealso: VecGetArrayF90(), VecGetArray(), VecRestoreArray(), VecRestoreArrayReadF90()
3021: M*/
3023: /*MC
3024: VecDestroyVecsF90 - Frees a block of vectors obtained with VecDuplicateVecsF90().
3026: Synopsis:
3027: VecDestroyVecsF90(PetscInt n,{Vec, pointer :: x(:)},PetscErrorCode ierr)
3029: Collective on Vec
3031: Input Parameters:
3032: + n - the number of vectors previously obtained
3033: - x - pointer to array of vector pointers
3035: Output Parameter:
3036: . ierr - error code
3038: Notes:
3039: Not yet supported for all F90 compilers
3041: Level: beginner
3043: .seealso: VecDestroyVecs(), VecDuplicateVecsF90()
3045: M*/
3047: /*MC
3048: VecGetArrayF90 - Accesses a vector array from Fortran90. For default PETSc
3049: vectors, VecGetArrayF90() returns a pointer to the local data array. Otherwise,
3050: this routine is implementation dependent. You MUST call VecRestoreArrayF90()
3051: when you no longer need access to the array.
3053: Synopsis:
3054: VecGetArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
3056: Logically Collective on Vec
3058: Input Parameter:
3059: . x - vector
3061: Output Parameters:
3062: + xx_v - the Fortran90 pointer to the array
3063: - ierr - error code
3065: Example of Usage:
3066: .vb
3067: #include <petsc/finclude/petscvec.h>
3068: use petscvec
3070: PetscScalar, pointer :: xx_v(:)
3071: ....
3072: call VecGetArrayF90(x,xx_v,ierr)
3073: xx_v(3) = a
3074: call VecRestoreArrayF90(x,xx_v,ierr)
3075: .ve
3077: If you ONLY intend to read entries from the array and not change any entries you should use VecGetArrayReadF90().
3079: Level: beginner
3081: .seealso: VecRestoreArrayF90(), VecGetArray(), VecRestoreArray(), VecGetArrayReadF90()
3083: M*/
3085: /*MC
3086: VecGetArrayReadF90 - Accesses a read only array from Fortran90. For default PETSc
3087: vectors, VecGetArrayF90() returns a pointer to the local data array. Otherwise,
3088: this routine is implementation dependent. You MUST call VecRestoreArrayReadF90()
3089: when you no longer need access to the array.
3091: Synopsis:
3092: VecGetArrayReadF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
3094: Logically Collective on Vec
3096: Input Parameter:
3097: . x - vector
3099: Output Parameters:
3100: + xx_v - the Fortran90 pointer to the array
3101: - ierr - error code
3103: Example of Usage:
3104: .vb
3105: #include <petsc/finclude/petscvec.h>
3106: use petscvec
3108: PetscScalar, pointer :: xx_v(:)
3109: ....
3110: call VecGetArrayReadF90(x,xx_v,ierr)
3111: a = xx_v(3)
3112: call VecRestoreArrayReadF90(x,xx_v,ierr)
3113: .ve
3115: If you intend to write entries into the array you must use VecGetArrayF90().
3117: Level: beginner
3119: .seealso: VecRestoreArrayReadF90(), VecGetArray(), VecRestoreArray(), VecGetArrayRead(), VecRestoreArrayRead(), VecGetArrayF90()
3121: M*/
3123: /*MC
3124: VecRestoreArrayReadF90 - Restores a readonly vector to a usable state after a call to
3125: VecGetArrayReadF90().
3127: Synopsis:
3128: VecRestoreArrayReadF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
3130: Logically Collective on Vec
3132: Input Parameters:
3133: + x - vector
3134: - xx_v - the Fortran90 pointer to the array
3136: Output Parameter:
3137: . ierr - error code
3139: Example of Usage:
3140: .vb
3141: #include <petsc/finclude/petscvec.h>
3142: use petscvec
3144: PetscScalar, pointer :: xx_v(:)
3145: ....
3146: call VecGetArrayReadF90(x,xx_v,ierr)
3147: a = xx_v(3)
3148: call VecRestoreArrayReadF90(x,xx_v,ierr)
3149: .ve
3151: Level: beginner
3153: .seealso: VecGetArrayReadF90(), VecGetArray(), VecRestoreArray(), VecGetArrayRead(), VecRestoreArrayRead(), VecRestoreArrayF90()
3155: M*/
3157: /*@C
3158: VecGetArray2d - Returns a pointer to a 2d contiguous array that contains this
3159: processor's portion of the vector data. You MUST call VecRestoreArray2d()
3160: when you no longer need access to the array.
3162: Logically Collective
3164: Input Parameter:
3165: + x - the vector
3166: . m - first dimension of two dimensional array
3167: . n - second dimension of two dimensional array
3168: . mstart - first index you will use in first coordinate direction (often 0)
3169: - nstart - first index in the second coordinate direction (often 0)
3171: Output Parameter:
3172: . a - location to put pointer to the array
3174: Level: developer
3176: Notes:
3177: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
3178: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3179: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3180: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
3182: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3184: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3185: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3186: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3187: @*/
3188: PetscErrorCode VecGetArray2d(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3189: {
3191: PetscInt i,N;
3192: PetscScalar *aa;
3198: VecGetLocalSize(x,&N);
3199: 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);
3200: VecGetArray(x,&aa);
3202: PetscMalloc1(m,a);
3203: for (i=0; i<m; i++) (*a)[i] = aa + i*n - nstart;
3204: *a -= mstart;
3205: return(0);
3206: }
3208: /*@C
3209: VecGetArray2dWrite - Returns a pointer to a 2d contiguous array that will contain this
3210: processor's portion of the vector data. You MUST call VecRestoreArray2dWrite()
3211: when you no longer need access to the array.
3213: Logically Collective
3215: Input Parameter:
3216: + x - the vector
3217: . m - first dimension of two dimensional array
3218: . n - second dimension of two dimensional array
3219: . mstart - first index you will use in first coordinate direction (often 0)
3220: - nstart - first index in the second coordinate direction (often 0)
3222: Output Parameter:
3223: . a - location to put pointer to the array
3225: Level: developer
3227: Notes:
3228: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
3229: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3230: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3231: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
3233: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3235: Concepts: vector^accessing local values as 2d array
3237: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3238: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3239: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3240: @*/
3241: PetscErrorCode VecGetArray2dWrite(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3242: {
3244: PetscInt i,N;
3245: PetscScalar *aa;
3251: VecGetLocalSize(x,&N);
3252: 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);
3253: VecGetArrayWrite(x,&aa);
3255: PetscMalloc1(m,a);
3256: for (i=0; i<m; i++) (*a)[i] = aa + i*n - nstart;
3257: *a -= mstart;
3258: return(0);
3259: }
3261: /*@C
3262: VecRestoreArray2d - Restores a vector after VecGetArray2d() has been called.
3264: Logically Collective
3266: Input Parameters:
3267: + x - the vector
3268: . m - first dimension of two dimensional array
3269: . n - second dimension of the two dimensional array
3270: . mstart - first index you will use in first coordinate direction (often 0)
3271: . nstart - first index in the second coordinate direction (often 0)
3272: - a - location of pointer to array obtained from VecGetArray2d()
3274: Level: developer
3276: Notes:
3277: For regular PETSc vectors this routine does not involve any copies. For
3278: any special vectors that do not store local vector data in a contiguous
3279: array, this routine will copy the data back into the underlying
3280: vector data structure from the array obtained with VecGetArray().
3282: This routine actually zeros out the a pointer.
3284: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3285: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3286: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3287: @*/
3288: PetscErrorCode VecRestoreArray2d(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3289: {
3291: void *dummy;
3297: dummy = (void*)(*a + mstart);
3298: PetscFree(dummy);
3299: VecRestoreArray(x,NULL);
3300: return(0);
3301: }
3303: /*@C
3304: VecRestoreArray2dWrite - Restores a vector after VecGetArray2dWrite() has been called.
3306: Logically Collective
3308: Input Parameters:
3309: + x - the vector
3310: . m - first dimension of two dimensional array
3311: . n - second dimension of the two dimensional array
3312: . mstart - first index you will use in first coordinate direction (often 0)
3313: . nstart - first index in the second coordinate direction (often 0)
3314: - a - location of pointer to array obtained from VecGetArray2d()
3316: Level: developer
3318: Notes:
3319: For regular PETSc vectors this routine does not involve any copies. For
3320: any special vectors that do not store local vector data in a contiguous
3321: array, this routine will copy the data back into the underlying
3322: vector data structure from the array obtained with VecGetArray().
3324: This routine actually zeros out the a pointer.
3326: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3327: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3328: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3329: @*/
3330: PetscErrorCode VecRestoreArray2dWrite(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3331: {
3333: void *dummy;
3339: dummy = (void*)(*a + mstart);
3340: PetscFree(dummy);
3341: VecRestoreArrayWrite(x,NULL);
3342: return(0);
3343: }
3345: /*@C
3346: VecGetArray1d - Returns a pointer to a 1d contiguous array that contains this
3347: processor's portion of the vector data. You MUST call VecRestoreArray1d()
3348: when you no longer need access to the array.
3350: Logically Collective
3352: Input Parameter:
3353: + x - the vector
3354: . m - first dimension of two dimensional array
3355: - mstart - first index you will use in first coordinate direction (often 0)
3357: Output Parameter:
3358: . a - location to put pointer to the array
3360: Level: developer
3362: Notes:
3363: For a vector obtained from DMCreateLocalVector() mstart are likely
3364: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3365: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
3367: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3369: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3370: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3371: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3372: @*/
3373: PetscErrorCode VecGetArray1d(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
3374: {
3376: PetscInt N;
3382: VecGetLocalSize(x,&N);
3383: if (m != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Local array size %D does not match 1d array dimensions %D",N,m);
3384: VecGetArray(x,a);
3385: *a -= mstart;
3386: return(0);
3387: }
3389: /*@C
3390: VecGetArray1dWrite - Returns a pointer to a 1d contiguous array that will contain this
3391: processor's portion of the vector data. You MUST call VecRestoreArray1dWrite()
3392: when you no longer need access to the array.
3394: Logically Collective
3396: Input Parameter:
3397: + x - the vector
3398: . m - first dimension of two dimensional array
3399: - mstart - first index you will use in first coordinate direction (often 0)
3401: Output Parameter:
3402: . a - location to put pointer to the array
3404: Level: developer
3406: Notes:
3407: For a vector obtained from DMCreateLocalVector() mstart are likely
3408: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3409: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
3411: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3413: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3414: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3415: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3416: @*/
3417: PetscErrorCode VecGetArray1dWrite(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
3418: {
3420: PetscInt N;
3426: VecGetLocalSize(x,&N);
3427: if (m != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Local array size %D does not match 1d array dimensions %D",N,m);
3428: VecGetArrayWrite(x,a);
3429: *a -= mstart;
3430: return(0);
3431: }
3433: /*@C
3434: VecRestoreArray1d - Restores a vector after VecGetArray1d() has been called.
3436: Logically Collective
3438: Input Parameters:
3439: + x - the vector
3440: . m - first dimension of two dimensional array
3441: . mstart - first index you will use in first coordinate direction (often 0)
3442: - a - location of pointer to array obtained from VecGetArray21()
3444: Level: developer
3446: Notes:
3447: For regular PETSc vectors this routine does not involve any copies. For
3448: any special vectors that do not store local vector data in a contiguous
3449: array, this routine will copy the data back into the underlying
3450: vector data structure from the array obtained with VecGetArray1d().
3452: This routine actually zeros out the a pointer.
3454: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3455: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3456: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
3457: @*/
3458: PetscErrorCode VecRestoreArray1d(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
3459: {
3465: VecRestoreArray(x,NULL);
3466: return(0);
3467: }
3469: /*@C
3470: VecRestoreArray1dWrite - Restores a vector after VecGetArray1dWrite() has been called.
3472: Logically Collective
3474: Input Parameters:
3475: + x - the vector
3476: . m - first dimension of two dimensional array
3477: . mstart - first index you will use in first coordinate direction (often 0)
3478: - a - location of pointer to array obtained from VecGetArray21()
3480: Level: developer
3482: Notes:
3483: For regular PETSc vectors this routine does not involve any copies. For
3484: any special vectors that do not store local vector data in a contiguous
3485: array, this routine will copy the data back into the underlying
3486: vector data structure from the array obtained with VecGetArray1d().
3488: This routine actually zeros out the a pointer.
3490: Concepts: vector^accessing local values as 1d array
3492: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3493: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3494: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
3495: @*/
3496: PetscErrorCode VecRestoreArray1dWrite(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
3497: {
3503: VecRestoreArrayWrite(x,NULL);
3504: return(0);
3505: }
3507: /*@C
3508: VecGetArray3d - Returns a pointer to a 3d contiguous array that contains this
3509: processor's portion of the vector data. You MUST call VecRestoreArray3d()
3510: when you no longer need access to the array.
3512: Logically Collective
3514: Input Parameter:
3515: + x - the vector
3516: . m - first dimension of three dimensional array
3517: . n - second dimension of three dimensional array
3518: . p - third dimension of three dimensional array
3519: . mstart - first index you will use in first coordinate direction (often 0)
3520: . nstart - first index in the second coordinate direction (often 0)
3521: - pstart - first index in the third coordinate direction (often 0)
3523: Output Parameter:
3524: . a - location to put pointer to the array
3526: Level: developer
3528: Notes:
3529: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3530: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3531: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3532: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3534: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3536: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3537: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3538: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3539: @*/
3540: PetscErrorCode VecGetArray3d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3541: {
3543: PetscInt i,N,j;
3544: PetscScalar *aa,**b;
3550: VecGetLocalSize(x,&N);
3551: 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);
3552: VecGetArray(x,&aa);
3554: PetscMalloc1(m*sizeof(PetscScalar**)+m*n,a);
3555: b = (PetscScalar**)((*a) + m);
3556: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3557: for (i=0; i<m; i++)
3558: for (j=0; j<n; j++)
3559: b[i*n+j] = aa + i*n*p + j*p - pstart;
3561: *a -= mstart;
3562: return(0);
3563: }
3565: /*@C
3566: VecGetArray3dWrite - Returns a pointer to a 3d contiguous array that will contain this
3567: processor's portion of the vector data. You MUST call VecRestoreArray3dWrite()
3568: when you no longer need access to the array.
3570: Logically Collective
3572: Input Parameter:
3573: + x - the vector
3574: . m - first dimension of three dimensional array
3575: . n - second dimension of three dimensional array
3576: . p - third dimension of three dimensional array
3577: . mstart - first index you will use in first coordinate direction (often 0)
3578: . nstart - first index in the second coordinate direction (often 0)
3579: - pstart - first index in the third coordinate direction (often 0)
3581: Output Parameter:
3582: . a - location to put pointer to the array
3584: Level: developer
3586: Notes:
3587: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3588: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3589: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3590: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3592: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3594: Concepts: vector^accessing local values as 3d array
3596: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3597: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3598: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3599: @*/
3600: PetscErrorCode VecGetArray3dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3601: {
3603: PetscInt i,N,j;
3604: PetscScalar *aa,**b;
3610: VecGetLocalSize(x,&N);
3611: 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);
3612: VecGetArrayWrite(x,&aa);
3614: PetscMalloc1(m*sizeof(PetscScalar**)+m*n,a);
3615: b = (PetscScalar**)((*a) + m);
3616: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3617: for (i=0; i<m; i++)
3618: for (j=0; j<n; j++)
3619: b[i*n+j] = aa + i*n*p + j*p - pstart;
3621: *a -= mstart;
3622: return(0);
3623: }
3625: /*@C
3626: VecRestoreArray3d - Restores a vector after VecGetArray3d() has been called.
3628: Logically Collective
3630: Input Parameters:
3631: + x - the vector
3632: . m - first dimension of three dimensional array
3633: . n - second dimension of the three dimensional array
3634: . p - third dimension of the three dimensional array
3635: . mstart - first index you will use in first coordinate direction (often 0)
3636: . nstart - first index in the second coordinate direction (often 0)
3637: . pstart - first index in the third coordinate direction (often 0)
3638: - a - location of pointer to array obtained from VecGetArray3d()
3640: Level: developer
3642: Notes:
3643: For regular PETSc vectors this routine does not involve any copies. For
3644: any special vectors that do not store local vector data in a contiguous
3645: array, this routine will copy the data back into the underlying
3646: vector data structure from the array obtained with VecGetArray().
3648: This routine actually zeros out the a pointer.
3650: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3651: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3652: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3653: @*/
3654: PetscErrorCode VecRestoreArray3d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3655: {
3657: void *dummy;
3663: dummy = (void*)(*a + mstart);
3664: PetscFree(dummy);
3665: VecRestoreArray(x,NULL);
3666: return(0);
3667: }
3669: /*@C
3670: VecRestoreArray3dWrite - Restores a vector after VecGetArray3dWrite() has been called.
3672: Logically Collective
3674: Input Parameters:
3675: + x - the vector
3676: . m - first dimension of three dimensional array
3677: . n - second dimension of the three dimensional array
3678: . p - third dimension of the three dimensional array
3679: . mstart - first index you will use in first coordinate direction (often 0)
3680: . nstart - first index in the second coordinate direction (often 0)
3681: . pstart - first index in the third coordinate direction (often 0)
3682: - a - location of pointer to array obtained from VecGetArray3d()
3684: Level: developer
3686: Notes:
3687: For regular PETSc vectors this routine does not involve any copies. For
3688: any special vectors that do not store local vector data in a contiguous
3689: array, this routine will copy the data back into the underlying
3690: vector data structure from the array obtained with VecGetArray().
3692: This routine actually zeros out the a pointer.
3694: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3695: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3696: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3697: @*/
3698: PetscErrorCode VecRestoreArray3dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3699: {
3701: void *dummy;
3707: dummy = (void*)(*a + mstart);
3708: PetscFree(dummy);
3709: VecRestoreArrayWrite(x,NULL);
3710: return(0);
3711: }
3713: /*@C
3714: VecGetArray4d - Returns a pointer to a 4d contiguous array that contains this
3715: processor's portion of the vector data. You MUST call VecRestoreArray4d()
3716: when you no longer need access to the array.
3718: Logically Collective
3720: Input Parameter:
3721: + x - the vector
3722: . m - first dimension of four dimensional array
3723: . n - second dimension of four dimensional array
3724: . p - third dimension of four dimensional array
3725: . q - fourth dimension of four dimensional array
3726: . mstart - first index you will use in first coordinate direction (often 0)
3727: . nstart - first index in the second coordinate direction (often 0)
3728: . pstart - first index in the third coordinate direction (often 0)
3729: - qstart - first index in the fourth coordinate direction (often 0)
3731: Output Parameter:
3732: . a - location to put pointer to the array
3734: Level: beginner
3736: Notes:
3737: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3738: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3739: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3740: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3742: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3744: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3745: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3746: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3747: @*/
3748: PetscErrorCode VecGetArray4d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3749: {
3751: PetscInt i,N,j,k;
3752: PetscScalar *aa,***b,**c;
3758: VecGetLocalSize(x,&N);
3759: 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);
3760: VecGetArray(x,&aa);
3762: PetscMalloc1(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p,a);
3763: b = (PetscScalar***)((*a) + m);
3764: c = (PetscScalar**)(b + m*n);
3765: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3766: for (i=0; i<m; i++)
3767: for (j=0; j<n; j++)
3768: b[i*n+j] = c + i*n*p + j*p - pstart;
3769: for (i=0; i<m; i++)
3770: for (j=0; j<n; j++)
3771: for (k=0; k<p; k++)
3772: c[i*n*p+j*p+k] = aa + i*n*p*q + j*p*q + k*q - qstart;
3773: *a -= mstart;
3774: return(0);
3775: }
3777: /*@C
3778: VecGetArray4dWrite - Returns a pointer to a 4d contiguous array that will contain this
3779: processor's portion of the vector data. You MUST call VecRestoreArray4dWrite()
3780: when you no longer need access to the array.
3782: Logically Collective
3784: Input Parameter:
3785: + x - the vector
3786: . m - first dimension of four dimensional array
3787: . n - second dimension of four dimensional array
3788: . p - third dimension of four dimensional array
3789: . q - fourth dimension of four dimensional array
3790: . mstart - first index you will use in first coordinate direction (often 0)
3791: . nstart - first index in the second coordinate direction (often 0)
3792: . pstart - first index in the third coordinate direction (often 0)
3793: - qstart - first index in the fourth coordinate direction (often 0)
3795: Output Parameter:
3796: . a - location to put pointer to the array
3798: Level: beginner
3800: Notes:
3801: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3802: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3803: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3804: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3806: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3808: Concepts: vector^accessing local values as 3d array
3810: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3811: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3812: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3813: @*/
3814: PetscErrorCode VecGetArray4dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3815: {
3817: PetscInt i,N,j,k;
3818: PetscScalar *aa,***b,**c;
3824: VecGetLocalSize(x,&N);
3825: 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);
3826: VecGetArrayWrite(x,&aa);
3828: PetscMalloc1(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p,a);
3829: b = (PetscScalar***)((*a) + m);
3830: c = (PetscScalar**)(b + m*n);
3831: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3832: for (i=0; i<m; i++)
3833: for (j=0; j<n; j++)
3834: b[i*n+j] = c + i*n*p + j*p - pstart;
3835: for (i=0; i<m; i++)
3836: for (j=0; j<n; j++)
3837: for (k=0; k<p; k++)
3838: c[i*n*p+j*p+k] = aa + i*n*p*q + j*p*q + k*q - qstart;
3839: *a -= mstart;
3840: return(0);
3841: }
3843: /*@C
3844: VecRestoreArray4d - Restores a vector after VecGetArray3d() has been called.
3846: Logically Collective
3848: Input Parameters:
3849: + x - the vector
3850: . m - first dimension of four dimensional array
3851: . n - second dimension of the four dimensional array
3852: . p - third dimension of the four dimensional array
3853: . q - fourth dimension of the four dimensional array
3854: . mstart - first index you will use in first coordinate direction (often 0)
3855: . nstart - first index in the second coordinate direction (often 0)
3856: . pstart - first index in the third coordinate direction (often 0)
3857: . qstart - first index in the fourth coordinate direction (often 0)
3858: - a - location of pointer to array obtained from VecGetArray4d()
3860: Level: beginner
3862: Notes:
3863: For regular PETSc vectors this routine does not involve any copies. For
3864: any special vectors that do not store local vector data in a contiguous
3865: array, this routine will copy the data back into the underlying
3866: vector data structure from the array obtained with VecGetArray().
3868: This routine actually zeros out the a pointer.
3870: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3871: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3872: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3873: @*/
3874: PetscErrorCode VecRestoreArray4d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3875: {
3877: void *dummy;
3883: dummy = (void*)(*a + mstart);
3884: PetscFree(dummy);
3885: VecRestoreArray(x,NULL);
3886: return(0);
3887: }
3889: /*@C
3890: VecRestoreArray4dWrite - Restores a vector after VecGetArray3dWrite() has been called.
3892: Logically Collective
3894: Input Parameters:
3895: + x - the vector
3896: . m - first dimension of four dimensional array
3897: . n - second dimension of the four dimensional array
3898: . p - third dimension of the four dimensional array
3899: . q - fourth dimension of the four dimensional array
3900: . mstart - first index you will use in first coordinate direction (often 0)
3901: . nstart - first index in the second coordinate direction (often 0)
3902: . pstart - first index in the third coordinate direction (often 0)
3903: . qstart - first index in the fourth coordinate direction (often 0)
3904: - a - location of pointer to array obtained from VecGetArray4d()
3906: Level: beginner
3908: Notes:
3909: For regular PETSc vectors this routine does not involve any copies. For
3910: any special vectors that do not store local vector data in a contiguous
3911: array, this routine will copy the data back into the underlying
3912: vector data structure from the array obtained with VecGetArray().
3914: This routine actually zeros out the a pointer.
3916: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3917: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3918: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3919: @*/
3920: PetscErrorCode VecRestoreArray4dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3921: {
3923: void *dummy;
3929: dummy = (void*)(*a + mstart);
3930: PetscFree(dummy);
3931: VecRestoreArrayWrite(x,NULL);
3932: return(0);
3933: }
3935: /*@C
3936: VecGetArray2dRead - Returns a pointer to a 2d contiguous array that contains this
3937: processor's portion of the vector data. You MUST call VecRestoreArray2dRead()
3938: when you no longer need access to the array.
3940: Logically Collective
3942: Input Parameter:
3943: + x - the vector
3944: . m - first dimension of two dimensional array
3945: . n - second dimension of two dimensional array
3946: . mstart - first index you will use in first coordinate direction (often 0)
3947: - nstart - first index in the second coordinate direction (often 0)
3949: Output Parameter:
3950: . a - location to put pointer to the array
3952: Level: developer
3954: Notes:
3955: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
3956: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3957: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3958: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
3960: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3962: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3963: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3964: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3965: @*/
3966: PetscErrorCode VecGetArray2dRead(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3967: {
3968: PetscErrorCode ierr;
3969: PetscInt i,N;
3970: const PetscScalar *aa;
3976: VecGetLocalSize(x,&N);
3977: 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);
3978: VecGetArrayRead(x,&aa);
3980: PetscMalloc1(m,a);
3981: for (i=0; i<m; i++) (*a)[i] = (PetscScalar*) aa + i*n - nstart;
3982: *a -= mstart;
3983: return(0);
3984: }
3986: /*@C
3987: VecRestoreArray2dRead - Restores a vector after VecGetArray2dRead() has been called.
3989: Logically Collective
3991: Input Parameters:
3992: + x - the vector
3993: . m - first dimension of two dimensional array
3994: . n - second dimension of the two dimensional array
3995: . mstart - first index you will use in first coordinate direction (often 0)
3996: . nstart - first index in the second coordinate direction (often 0)
3997: - a - location of pointer to array obtained from VecGetArray2d()
3999: Level: developer
4001: Notes:
4002: For regular PETSc vectors this routine does not involve any copies. For
4003: any special vectors that do not store local vector data in a contiguous
4004: array, this routine will copy the data back into the underlying
4005: vector data structure from the array obtained with VecGetArray().
4007: This routine actually zeros out the a pointer.
4009: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
4010: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
4011: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
4012: @*/
4013: PetscErrorCode VecRestoreArray2dRead(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
4014: {
4016: void *dummy;
4022: dummy = (void*)(*a + mstart);
4023: PetscFree(dummy);
4024: VecRestoreArrayRead(x,NULL);
4025: return(0);
4026: }
4028: /*@C
4029: VecGetArray1dRead - Returns a pointer to a 1d contiguous array that contains this
4030: processor's portion of the vector data. You MUST call VecRestoreArray1dRead()
4031: when you no longer need access to the array.
4033: Logically Collective
4035: Input Parameter:
4036: + x - the vector
4037: . m - first dimension of two dimensional array
4038: - mstart - first index you will use in first coordinate direction (often 0)
4040: Output Parameter:
4041: . a - location to put pointer to the array
4043: Level: developer
4045: Notes:
4046: For a vector obtained from DMCreateLocalVector() mstart are likely
4047: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
4048: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
4050: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
4052: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
4053: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
4054: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
4055: @*/
4056: PetscErrorCode VecGetArray1dRead(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
4057: {
4059: PetscInt N;
4065: VecGetLocalSize(x,&N);
4066: if (m != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Local array size %D does not match 1d array dimensions %D",N,m);
4067: VecGetArrayRead(x,(const PetscScalar**)a);
4068: *a -= mstart;
4069: return(0);
4070: }
4072: /*@C
4073: VecRestoreArray1dRead - Restores a vector after VecGetArray1dRead() has been called.
4075: Logically Collective
4077: Input Parameters:
4078: + x - the vector
4079: . m - first dimension of two dimensional array
4080: . mstart - first index you will use in first coordinate direction (often 0)
4081: - a - location of pointer to array obtained from VecGetArray21()
4083: Level: developer
4085: Notes:
4086: For regular PETSc vectors this routine does not involve any copies. For
4087: any special vectors that do not store local vector data in a contiguous
4088: array, this routine will copy the data back into the underlying
4089: vector data structure from the array obtained with VecGetArray1dRead().
4091: This routine actually zeros out the a pointer.
4093: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
4094: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
4095: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
4096: @*/
4097: PetscErrorCode VecRestoreArray1dRead(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
4098: {
4104: VecRestoreArrayRead(x,NULL);
4105: return(0);
4106: }
4109: /*@C
4110: VecGetArray3dRead - Returns a pointer to a 3d contiguous array that contains this
4111: processor's portion of the vector data. You MUST call VecRestoreArray3dRead()
4112: when you no longer need access to the array.
4114: Logically Collective
4116: Input Parameter:
4117: + x - the vector
4118: . m - first dimension of three dimensional array
4119: . n - second dimension of three dimensional array
4120: . p - third dimension of three dimensional array
4121: . mstart - first index you will use in first coordinate direction (often 0)
4122: . nstart - first index in the second coordinate direction (often 0)
4123: - pstart - first index in the third coordinate direction (often 0)
4125: Output Parameter:
4126: . a - location to put pointer to the array
4128: Level: developer
4130: Notes:
4131: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
4132: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
4133: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
4134: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3dRead().
4136: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
4138: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
4139: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
4140: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
4141: @*/
4142: PetscErrorCode VecGetArray3dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
4143: {
4144: PetscErrorCode ierr;
4145: PetscInt i,N,j;
4146: const PetscScalar *aa;
4147: PetscScalar **b;
4153: VecGetLocalSize(x,&N);
4154: 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);
4155: VecGetArrayRead(x,&aa);
4157: PetscMalloc1(m*sizeof(PetscScalar**)+m*n,a);
4158: b = (PetscScalar**)((*a) + m);
4159: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
4160: for (i=0; i<m; i++)
4161: for (j=0; j<n; j++)
4162: b[i*n+j] = (PetscScalar *)aa + i*n*p + j*p - pstart;
4164: *a -= mstart;
4165: return(0);
4166: }
4168: /*@C
4169: VecRestoreArray3dRead - Restores a vector after VecGetArray3dRead() has been called.
4171: Logically Collective
4173: Input Parameters:
4174: + x - the vector
4175: . m - first dimension of three dimensional array
4176: . n - second dimension of the three dimensional array
4177: . p - third dimension of the three dimensional array
4178: . mstart - first index you will use in first coordinate direction (often 0)
4179: . nstart - first index in the second coordinate direction (often 0)
4180: . pstart - first index in the third coordinate direction (often 0)
4181: - a - location of pointer to array obtained from VecGetArray3dRead()
4183: Level: developer
4185: Notes:
4186: For regular PETSc vectors this routine does not involve any copies. For
4187: any special vectors that do not store local vector data in a contiguous
4188: array, this routine will copy the data back into the underlying
4189: vector data structure from the array obtained with VecGetArray().
4191: This routine actually zeros out the a pointer.
4193: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
4194: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
4195: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
4196: @*/
4197: PetscErrorCode VecRestoreArray3dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
4198: {
4200: void *dummy;
4206: dummy = (void*)(*a + mstart);
4207: PetscFree(dummy);
4208: VecRestoreArrayRead(x,NULL);
4209: return(0);
4210: }
4212: /*@C
4213: VecGetArray4dRead - Returns a pointer to a 4d contiguous array that contains this
4214: processor's portion of the vector data. You MUST call VecRestoreArray4dRead()
4215: when you no longer need access to the array.
4217: Logically Collective
4219: Input Parameter:
4220: + x - the vector
4221: . m - first dimension of four dimensional array
4222: . n - second dimension of four dimensional array
4223: . p - third dimension of four dimensional array
4224: . q - fourth dimension of four dimensional array
4225: . mstart - first index you will use in first coordinate direction (often 0)
4226: . nstart - first index in the second coordinate direction (often 0)
4227: . pstart - first index in the third coordinate direction (often 0)
4228: - qstart - first index in the fourth coordinate direction (often 0)
4230: Output Parameter:
4231: . a - location to put pointer to the array
4233: Level: beginner
4235: Notes:
4236: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
4237: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
4238: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
4239: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
4241: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
4243: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
4244: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
4245: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
4246: @*/
4247: PetscErrorCode VecGetArray4dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
4248: {
4249: PetscErrorCode ierr;
4250: PetscInt i,N,j,k;
4251: const PetscScalar *aa;
4252: PetscScalar ***b,**c;
4258: VecGetLocalSize(x,&N);
4259: 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);
4260: VecGetArrayRead(x,&aa);
4262: PetscMalloc1(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p,a);
4263: b = (PetscScalar***)((*a) + m);
4264: c = (PetscScalar**)(b + m*n);
4265: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
4266: for (i=0; i<m; i++)
4267: for (j=0; j<n; j++)
4268: b[i*n+j] = c + i*n*p + j*p - pstart;
4269: for (i=0; i<m; i++)
4270: for (j=0; j<n; j++)
4271: for (k=0; k<p; k++)
4272: c[i*n*p+j*p+k] = (PetscScalar*) aa + i*n*p*q + j*p*q + k*q - qstart;
4273: *a -= mstart;
4274: return(0);
4275: }
4277: /*@C
4278: VecRestoreArray4dRead - Restores a vector after VecGetArray3d() has been called.
4280: Logically Collective
4282: Input Parameters:
4283: + x - the vector
4284: . m - first dimension of four dimensional array
4285: . n - second dimension of the four dimensional array
4286: . p - third dimension of the four dimensional array
4287: . q - fourth dimension of the four dimensional array
4288: . mstart - first index you will use in first coordinate direction (often 0)
4289: . nstart - first index in the second coordinate direction (often 0)
4290: . pstart - first index in the third coordinate direction (often 0)
4291: . qstart - first index in the fourth coordinate direction (often 0)
4292: - a - location of pointer to array obtained from VecGetArray4dRead()
4294: Level: beginner
4296: Notes:
4297: For regular PETSc vectors this routine does not involve any copies. For
4298: any special vectors that do not store local vector data in a contiguous
4299: array, this routine will copy the data back into the underlying
4300: vector data structure from the array obtained with VecGetArray().
4302: This routine actually zeros out the a pointer.
4304: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
4305: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
4306: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
4307: @*/
4308: PetscErrorCode VecRestoreArray4dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
4309: {
4311: void *dummy;
4317: dummy = (void*)(*a + mstart);
4318: PetscFree(dummy);
4319: VecRestoreArrayRead(x,NULL);
4320: return(0);
4321: }
4323: #if defined(PETSC_USE_DEBUG)
4325: /*@
4326: VecLockGet - Gets the current lock status of a vector
4328: Logically Collective on Vec
4330: Input Parameter:
4331: . x - the vector
4333: Output Parameter:
4334: . state - greater than zero indicates the vector is locked for read; less then zero indicates the vector is
4335: locked for write; equal to zero means the vector is unlocked, that is, it is free to read or write.
4337: Level: beginner
4339: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPush(), VecLockReadPop()
4340: @*/
4341: PetscErrorCode VecLockGet(Vec x,PetscInt *state)
4342: {
4345: *state = x->lock;
4346: return(0);
4347: }
4349: /*@
4350: VecLockReadPush - Pushes a read-only lock on a vector to prevent it from writing
4352: Logically Collective on Vec
4354: Input Parameter:
4355: . x - the vector
4357: Notes:
4358: If this is set then calls to VecGetArray() or VecSetValues() or any other routines that change the vectors values will fail.
4360: The call can be nested, i.e., called multiple times on the same vector, but each VecLockReadPush(x) has to have one matching
4361: VecLockReadPop(x), which removes the latest read-only lock.
4363: Level: beginner
4365: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPop(), VecLockGet()
4366: @*/
4367: PetscErrorCode VecLockReadPush(Vec x)
4368: {
4371: if (x->lock < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Vector is already locked for exclusive write access but you want to read it");
4372: x->lock++;
4373: return(0);
4374: }
4376: /*@
4377: VecLockReadPop - Pops a read-only lock from a vector
4379: Logically Collective on Vec
4381: Input Parameter:
4382: . x - the vector
4384: Level: beginner
4386: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPush(), VecLockGet()
4387: @*/
4388: PetscErrorCode VecLockReadPop(Vec x)
4389: {
4392: x->lock--;
4393: if (x->lock < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Vector has been unlocked from read-only access too many times");
4394: return(0);
4395: }
4397: /*@C
4398: VecLockWriteSet_Private - Lock or unlock a vector for exclusive read/write access
4400: Logically Collective on Vec
4402: Input Parameter:
4403: + x - the vector
4404: - flg - PETSC_TRUE to lock the vector for writing; PETSC_FALSE to unlock it.
4406: Notes:
4407: The function is usefull in split-phase computations, which usually have a begin phase and an end phase.
4408: One can call VecLockWriteSet_Private(x,PETSC_TRUE) in the begin phase to lock a vector for exclusive
4409: access, and call VecLockWriteSet_Private(x,PETSC_FALSE) in the end phase to unlock the vector from exclusive
4410: access. In this way, one is ensured no other operations can access the vector in between. The code may like
4413: VecGetArray(x,&xdata); // begin phase
4414: VecLockWriteSet_Private(v,PETSC_TRUE);
4416: Other operations, which can not acceess x anymore (they can access xdata, of course)
4418: VecRestoreArray(x,&vdata); // end phase
4419: VecLockWriteSet_Private(v,PETSC_FALSE);
4421: The call can not be nested on the same vector, in other words, one can not call VecLockWriteSet_Private(x,PETSC_TRUE)
4422: again before calling VecLockWriteSet_Private(v,PETSC_FALSE).
4424: Level: beginner
4426: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPush(), VecLockReadPop(), VecLockGet()
4427: @*/
4428: PetscErrorCode VecLockWriteSet_Private(Vec x,PetscBool flg)
4429: {
4432: if (flg) {
4433: if (x->lock > 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Vector is already locked for read-only access but you want to write it");
4434: else if (x->lock < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Vector is already locked for exclusive write access but you want to write it");
4435: else x->lock = -1;
4436: } else {
4437: if (x->lock != -1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Vector is not locked for exclusive write access but you want to unlock it from that");
4438: x->lock = 0;
4439: }
4440: return(0);
4441: }
4443: /*@
4444: VecLockPush - Pushes a read-only lock on a vector to prevent it from writing
4446: Level: deprecated
4448: .seealso: VecLockReadPush()
4449: @*/
4450: PetscErrorCode VecLockPush(Vec x)
4451: {
4454: VecLockReadPush(x);
4455: return(0);
4456: }
4458: /*@
4459: VecLockPop - Pops a read-only lock from a vector
4461: Level: deprecated
4463: .seealso: VecLockReadPop()
4464: @*/
4465: PetscErrorCode VecLockPop(Vec x)
4466: {
4469: VecLockReadPop(x);
4470: return(0);
4471: }
4473: #endif