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 imbalance 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
112: .seealso: VecMDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecDotRealPart()
113: @*/
114: PetscErrorCode VecDot(Vec x,Vec y,PetscScalar *val)
115: {
125: VecCheckSameSize(x,1,y,2);
127: PetscLogEventBegin(VEC_Dot,x,y,0,0);
128: (*x->ops->dot)(x,y,val);
129: PetscLogEventEnd(VEC_Dot,x,y,0,0);
130: return(0);
131: }
133: /*@
134: VecDotRealPart - Computes the real part of the vector dot product.
136: Collective on Vec
138: Input Parameters:
139: . x, y - the vectors
141: Output Parameter:
142: . val - the real part of the dot product;
144: Performance Issues:
145: $ per-processor memory bandwidth
146: $ interprocessor latency
147: $ work load imbalance that causes certain processes to arrive much earlier than others
149: Notes for Users of Complex Numbers:
150: See VecDot() for more details on the definition of the dot product for complex numbers
152: For real numbers this returns the same value as VecDot()
154: 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
155: the space R^{2n} (that is a vector of 2n components with the real or imaginary part of the complex numbers for components)
157: Developer Note: This is not currently optimized to compute only the real part of the dot product.
159: Level: intermediate
161: .seealso: VecMDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecDot(), VecDotNorm2()
162: @*/
163: PetscErrorCode VecDotRealPart(Vec x,Vec y,PetscReal *val)
164: {
166: PetscScalar fdot;
169: VecDot(x,y,&fdot);
170: *val = PetscRealPart(fdot);
171: return(0);
172: }
174: /*@
175: VecNorm - Computes the vector norm.
177: Collective on Vec
179: Input Parameters:
180: + x - the vector
181: - type - one of NORM_1, NORM_2, NORM_INFINITY. Also available
182: NORM_1_AND_2, which computes both norms and stores them
183: in a two element array.
185: Output Parameter:
186: . val - the norm
188: Notes:
189: $ NORM_1 denotes sum_i |x_i|
190: $ NORM_2 denotes sqrt(sum_i |x_i|^2)
191: $ NORM_INFINITY denotes max_i |x_i|
193: For complex numbers NORM_1 will return the traditional 1 norm of the 2 norm of the complex numbers; that is the 1
194: norm of the absolute values of the complex entries. In PETSc 3.6 and earlier releases it returned the 1 norm of
195: the 1 norm of the complex entries (what is returned by the BLAS routine asum()). Both are valid norms but most
196: people expect the former.
198: Level: intermediate
200: Performance Issues:
201: $ per-processor memory bandwidth
202: $ interprocessor latency
203: $ work load imbalance that causes certain processes to arrive much earlier than others
205: .seealso: VecDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecNormAvailable(),
206: VecNormBegin(), VecNormEnd()
208: @*/
210: PetscErrorCode VecNorm(Vec x,NormType type,PetscReal *val)
211: {
212: PetscBool flg;
220: /*
221: * Cached data?
222: */
223: if (type!=NORM_1_AND_2) {
224: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[type],*val,flg);
225: if (flg) return(0);
226: }
227: PetscLogEventBegin(VEC_Norm,x,0,0,0);
228: (*x->ops->norm)(x,type,val);
229: PetscLogEventEnd(VEC_Norm,x,0,0,0);
230: if (type!=NORM_1_AND_2) {
231: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[type],*val);
232: }
233: return(0);
234: }
236: /*@
237: VecNormAvailable - Returns the vector norm if it is already known.
239: Not Collective
241: Input Parameters:
242: + x - the vector
243: - type - one of NORM_1, NORM_2, NORM_INFINITY. Also available
244: NORM_1_AND_2, which computes both norms and stores them
245: in a two element array.
247: Output Parameters:
248: + available - PETSC_TRUE if the val returned is valid
249: - val - the norm
251: Notes:
252: $ NORM_1 denotes sum_i |x_i|
253: $ NORM_2 denotes sqrt(sum_i (x_i)^2)
254: $ NORM_INFINITY denotes max_i |x_i|
256: Level: intermediate
258: Performance Issues:
259: $ per-processor memory bandwidth
260: $ interprocessor latency
261: $ work load imbalance that causes certain processes to arrive much earlier than others
263: Compile Option:
264: PETSC_HAVE_SLOW_BLAS_NORM2 will cause a C (loop unrolled) version of the norm to be used, rather
265: than the BLAS. This should probably only be used when one is using the FORTRAN BLAS routines
266: (as opposed to vendor provided) because the FORTRAN BLAS NRM2() routine is very slow.
268: .seealso: VecDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecNorm()
269: VecNormBegin(), VecNormEnd()
271: @*/
272: PetscErrorCode VecNormAvailable(Vec x,NormType type,PetscBool *available,PetscReal *val)
273: {
281: *available = PETSC_FALSE;
282: if (type!=NORM_1_AND_2) {
283: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[type],*val,*available);
284: }
285: return(0);
286: }
288: /*@
289: VecNormalize - Normalizes a vector by 2-norm.
291: Collective on Vec
293: Input Parameter:
294: . x - the vector
296: Output Parameter:
297: . val - the vector norm before normalization
299: Level: intermediate
301: @*/
302: PetscErrorCode VecNormalize(Vec x,PetscReal *val)
303: {
305: PetscReal norm;
310: PetscLogEventBegin(VEC_Normalize,x,0,0,0);
311: VecNorm(x,NORM_2,&norm);
312: if (norm == 0.0) {
313: PetscInfo(x,"Vector of zero norm can not be normalized; Returning only the zero norm\n");
314: } else if (norm != 1.0) {
315: PetscScalar tmp = 1.0/norm;
316: VecScale(x,tmp);
317: }
318: if (val) *val = norm;
319: PetscLogEventEnd(VEC_Normalize,x,0,0,0);
320: return(0);
321: }
323: /*@C
324: VecMax - Determines the vector component with maximum real part and its location.
326: Collective on Vec
328: Input Parameter:
329: . x - the vector
331: Output Parameters:
332: + p - the location of val (pass NULL if you don't want this)
333: - val - the maximum component
335: Notes:
336: Returns the value PETSC_MIN_REAL and negative p if the vector is of length 0.
338: Returns the smallest index with the maximum value
339: Level: intermediate
341: .seealso: VecNorm(), VecMin()
342: @*/
343: PetscErrorCode VecMax(Vec x,PetscInt *p,PetscReal *val)
344: {
351: PetscLogEventBegin(VEC_Max,x,0,0,0);
352: (*x->ops->max)(x,p,val);
353: PetscLogEventEnd(VEC_Max,x,0,0,0);
354: return(0);
355: }
357: /*@C
358: VecMin - Determines the vector component with minimum real part and its location.
360: Collective on Vec
362: Input Parameter:
363: . x - the vector
365: Output Parameters:
366: + p - the location of val (pass NULL if you don't want this location)
367: - val - the minimum component
369: Level: intermediate
371: Notes:
372: Returns the value PETSC_MAX_REAL and negative p if the vector is of length 0.
374: This returns the smallest index with the minumum value
376: .seealso: VecMax()
377: @*/
378: PetscErrorCode VecMin(Vec x,PetscInt *p,PetscReal *val)
379: {
386: PetscLogEventBegin(VEC_Min,x,0,0,0);
387: (*x->ops->min)(x,p,val);
388: PetscLogEventEnd(VEC_Min,x,0,0,0);
389: return(0);
390: }
392: /*@
393: VecTDot - Computes an indefinite vector dot product. That is, this
394: routine does NOT use the complex conjugate.
396: Collective on Vec
398: Input Parameters:
399: . x, y - the vectors
401: Output Parameter:
402: . val - the dot product
404: Notes for Users of Complex Numbers:
405: For complex vectors, VecTDot() computes the indefinite form
406: $ val = (x,y) = y^T x,
407: where y^T denotes the transpose of y.
409: Use VecDot() for the inner product
410: $ val = (x,y) = y^H x,
411: where y^H denotes the conjugate transpose of y.
413: Level: intermediate
415: .seealso: VecDot(), VecMTDot()
416: @*/
417: PetscErrorCode VecTDot(Vec x,Vec y,PetscScalar *val)
418: {
428: VecCheckSameSize(x,1,y,2);
430: PetscLogEventBegin(VEC_TDot,x,y,0,0);
431: (*x->ops->tdot)(x,y,val);
432: PetscLogEventEnd(VEC_TDot,x,y,0,0);
433: return(0);
434: }
436: /*@
437: VecScale - Scales a vector.
439: Not collective on Vec
441: Input Parameters:
442: + x - the vector
443: - alpha - the scalar
445: Note:
446: For a vector with n components, VecScale() computes
447: $ x[i] = alpha * x[i], for i=1,...,n.
449: Level: intermediate
451: @*/
452: PetscErrorCode VecScale(Vec x, PetscScalar alpha)
453: {
454: PetscReal norms[4] = {0.0,0.0,0.0, 0.0};
455: PetscBool flgs[4];
457: PetscInt i;
462: if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
463: PetscLogEventBegin(VEC_Scale,x,0,0,0);
464: if (alpha != (PetscScalar)1.0) {
465: VecSetErrorIfLocked(x,1);
466: /* get current stashed norms */
467: for (i=0; i<4; i++) {
468: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[i],norms[i],flgs[i]);
469: }
470: (*x->ops->scale)(x,alpha);
471: PetscObjectStateIncrease((PetscObject)x);
472: /* put the scaled stashed norms back into the Vec */
473: for (i=0; i<4; i++) {
474: if (flgs[i]) {
475: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[i],PetscAbsScalar(alpha)*norms[i]);
476: }
477: }
478: }
479: PetscLogEventEnd(VEC_Scale,x,0,0,0);
480: return(0);
481: }
483: /*@
484: VecSet - Sets all components of a vector to a single scalar value.
486: Logically Collective on Vec
488: Input Parameters:
489: + x - the vector
490: - alpha - the scalar
492: Output Parameter:
493: . x - the vector
495: Note:
496: For a vector of dimension n, VecSet() computes
497: $ x[i] = alpha, for i=1,...,n,
498: so that all vector entries then equal the identical
499: scalar value, alpha. Use the more general routine
500: VecSetValues() to set different vector entries.
502: You CANNOT call this after you have called VecSetValues() but before you call
503: VecAssemblyBegin/End().
505: Level: beginner
507: .seealso VecSetValues(), VecSetValuesBlocked(), VecSetRandom()
509: @*/
510: PetscErrorCode VecSet(Vec x,PetscScalar alpha)
511: {
512: PetscReal val;
518: 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()");
520: VecSetErrorIfLocked(x,1);
522: PetscLogEventBegin(VEC_Set,x,0,0,0);
523: (*x->ops->set)(x,alpha);
524: PetscLogEventEnd(VEC_Set,x,0,0,0);
525: PetscObjectStateIncrease((PetscObject)x);
527: /* norms can be simply set (if |alpha|*N not too large) */
528: val = PetscAbsScalar(alpha);
529: if (x->map->N == 0) {
530: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_1],0.0l);
531: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],0.0);
532: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_2],0.0);
533: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_FROBENIUS],0.0);
534: } else if (val > PETSC_MAX_REAL/x->map->N) {
535: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],val);
536: } else {
537: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_1],x->map->N * val);
538: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],val);
539: val = PetscSqrtReal((PetscReal)x->map->N) * val;
540: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_2],val);
541: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_FROBENIUS],val);
542: }
543: return(0);
544: }
546: /*@
547: VecAXPY - Computes y = alpha x + y.
549: Logically Collective on Vec
551: Input Parameters:
552: + alpha - the scalar
553: - x, y - the vectors
555: Output Parameter:
556: . y - output vector
558: Level: intermediate
560: Notes:
561: x and y MUST be different vectors
562: This routine is optimized for alpha of 0.0, otherwise it calls the BLAS routine
564: $ VecAXPY(y,alpha,x) y = alpha x + y
565: $ VecAYPX(y,beta,x) y = x + beta y
566: $ VecAXPBY(y,alpha,beta,x) y = alpha x + beta y
567: $ VecWAXPY(w,alpha,x,y) w = alpha x + y
568: $ VecAXPBYPCZ(w,alpha,beta,gamma,x,y) z = alpha x + beta y + gamma z
569: $ VecMAXPY(y,nv,alpha[],x[]) y = sum alpha[i] x[i] + y
571: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPBYPCZ(), VecAXPBY()
572: @*/
573: PetscErrorCode VecAXPY(Vec y,PetscScalar alpha,Vec x)
574: {
583: VecCheckSameSize(x,3,y,1);
584: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y cannot be the same vector");
586: if (alpha == (PetscScalar)0.0) return(0);
587: VecSetErrorIfLocked(y,1);
589: VecLockReadPush(x);
590: PetscLogEventBegin(VEC_AXPY,x,y,0,0);
591: (*y->ops->axpy)(y,alpha,x);
592: PetscLogEventEnd(VEC_AXPY,x,y,0,0);
593: VecLockReadPop(x);
594: PetscObjectStateIncrease((PetscObject)y);
595: return(0);
596: }
598: /*@
599: VecAXPBY - Computes y = alpha x + beta y.
601: Logically Collective on Vec
603: Input Parameters:
604: + alpha,beta - the scalars
605: - x, y - the vectors
607: Output Parameter:
608: . y - output vector
610: Level: intermediate
612: Notes:
613: x and y MUST be different vectors
614: The implementation is optimized for alpha and/or beta values of 0.0 and 1.0
616: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPY(), VecAXPBYPCZ()
617: @*/
618: PetscErrorCode VecAXPBY(Vec y,PetscScalar alpha,PetscScalar beta,Vec x)
619: {
628: VecCheckSameSize(y,1,x,4);
629: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y cannot be the same vector");
632: if (alpha == (PetscScalar)0.0 && beta == (PetscScalar)1.0) return(0);
633: VecSetErrorIfLocked(y,1);
634: PetscLogEventBegin(VEC_AXPY,x,y,0,0);
635: (*y->ops->axpby)(y,alpha,beta,x);
636: PetscLogEventEnd(VEC_AXPY,x,y,0,0);
637: PetscObjectStateIncrease((PetscObject)y);
638: return(0);
639: }
641: /*@
642: VecAXPBYPCZ - Computes z = alpha x + beta y + gamma z
644: Logically Collective on Vec
646: Input Parameters:
647: + alpha,beta, gamma - the scalars
648: - x, y, z - the vectors
650: Output Parameter:
651: . z - output vector
653: Level: intermediate
655: Notes:
656: x, y and z must be different vectors
657: The implementation is optimized for alpha of 1.0 and gamma of 1.0 or 0.0
659: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPY(), VecAXPBY()
660: @*/
661: PetscErrorCode VecAXPBYPCZ(Vec z,PetscScalar alpha,PetscScalar beta,PetscScalar gamma,Vec x,Vec y)
662: {
674: VecCheckSameSize(x,1,y,5);
675: VecCheckSameSize(x,1,z,6);
676: if (x == y || x == z) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x, y, and z must be different vectors");
677: if (y == z) SETERRQ(PetscObjectComm((PetscObject)y),PETSC_ERR_ARG_IDN,"x, y, and z must be different vectors");
681: if (alpha == (PetscScalar)0.0 && beta == (PetscScalar)0.0 && gamma == (PetscScalar)1.0) return(0);
682: VecSetErrorIfLocked(z,1);
684: PetscLogEventBegin(VEC_AXPBYPCZ,x,y,z,0);
685: (*y->ops->axpbypcz)(z,alpha,beta,gamma,x,y);
686: PetscLogEventEnd(VEC_AXPBYPCZ,x,y,z,0);
687: PetscObjectStateIncrease((PetscObject)z);
688: return(0);
689: }
691: /*@
692: VecAYPX - Computes y = x + beta y.
694: Logically Collective on Vec
696: Input Parameters:
697: + beta - the scalar
698: - x, y - the vectors
700: Output Parameter:
701: . y - output vector
703: Level: intermediate
705: Notes:
706: x and y MUST be different vectors
707: The implementation is optimized for beta of -1.0, 0.0, and 1.0
709: .seealso: VecMAXPY(), VecWAXPY(), VecAXPY(), VecAXPBYPCZ(), VecAXPBY()
710: @*/
711: PetscErrorCode VecAYPX(Vec y,PetscScalar beta,Vec x)
712: {
721: VecCheckSameSize(x,1,y,3);
722: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y must be different vectors");
724: VecSetErrorIfLocked(y,1);
726: PetscLogEventBegin(VEC_AYPX,x,y,0,0);
727: (*y->ops->aypx)(y,beta,x);
728: PetscLogEventEnd(VEC_AYPX,x,y,0,0);
729: PetscObjectStateIncrease((PetscObject)y);
730: return(0);
731: }
733: /*@
734: VecWAXPY - Computes w = alpha x + y.
736: Logically Collective on Vec
738: Input Parameters:
739: + alpha - the scalar
740: - x, y - the vectors
742: Output Parameter:
743: . w - the result
745: Level: intermediate
747: Notes:
748: w cannot be either x or y, but x and y can be the same
749: The implementation is optimzed for alpha of -1.0, 0.0, and 1.0
751: .seealso: VecAXPY(), VecAYPX(), VecAXPBY(), VecMAXPY(), VecAXPBYPCZ()
752: @*/
753: PetscErrorCode VecWAXPY(Vec w,PetscScalar alpha,Vec x,Vec y)
754: {
766: VecCheckSameSize(x,3,y,4);
767: VecCheckSameSize(x,3,w,1);
768: if (w == y) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Result vector w cannot be same as input vector y, suggest VecAXPY()");
769: if (w == x) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Result vector w cannot be same as input vector x, suggest VecAYPX()");
771: VecSetErrorIfLocked(w,1);
773: PetscLogEventBegin(VEC_WAXPY,x,y,w,0);
774: (*w->ops->waxpy)(w,alpha,x,y);
775: PetscLogEventEnd(VEC_WAXPY,x,y,w,0);
776: PetscObjectStateIncrease((PetscObject)w);
777: return(0);
778: }
780: /*@C
781: VecSetValues - Inserts or adds values into certain locations of a vector.
783: Not Collective
785: Input Parameters:
786: + x - vector to insert in
787: . ni - number of elements to add
788: . ix - indices where to add
789: . y - array of values
790: - iora - either INSERT_VALUES or ADD_VALUES, where
791: ADD_VALUES adds values to any existing entries, and
792: INSERT_VALUES replaces existing entries with new values
794: Notes:
795: VecSetValues() sets x[ix[i]] = y[i], for i=0,...,ni-1.
797: Calls to VecSetValues() with the INSERT_VALUES and ADD_VALUES
798: options cannot be mixed without intervening calls to the assembly
799: routines.
801: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
802: MUST be called after all calls to VecSetValues() have been completed.
804: VecSetValues() uses 0-based indices in Fortran as well as in C.
806: If you call VecSetOption(x, VEC_IGNORE_NEGATIVE_INDICES,PETSC_TRUE),
807: negative indices may be passed in ix. These rows are
808: simply ignored. This allows easily inserting element load matrices
809: with homogeneous Dirchlet boundary conditions that you don't want represented
810: in the vector.
812: Level: beginner
814: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValuesLocal(),
815: VecSetValue(), VecSetValuesBlocked(), InsertMode, INSERT_VALUES, ADD_VALUES, VecGetValues()
816: @*/
817: PetscErrorCode VecSetValues(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
818: {
823: if (!ni) return(0);
828: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
829: (*x->ops->setvalues)(x,ni,ix,y,iora);
830: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
831: PetscObjectStateIncrease((PetscObject)x);
832: return(0);
833: }
835: /*@C
836: VecGetValues - Gets values from certain locations of a vector. Currently
837: can only get values on the same processor
839: Not Collective
841: Input Parameters:
842: + x - vector to get values from
843: . ni - number of elements to get
844: - ix - indices where to get them from (in global 1d numbering)
846: Output Parameter:
847: . y - array of values
849: Notes:
850: The user provides the allocated array y; it is NOT allocated in this routine
852: VecGetValues() gets y[i] = x[ix[i]], for i=0,...,ni-1.
854: VecAssemblyBegin() and VecAssemblyEnd() MUST be called before calling this
856: VecGetValues() uses 0-based indices in Fortran as well as in C.
858: If you call VecSetOption(x, VEC_IGNORE_NEGATIVE_INDICES,PETSC_TRUE),
859: negative indices may be passed in ix. These rows are
860: simply ignored.
862: Level: beginner
864: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues()
865: @*/
866: PetscErrorCode VecGetValues(Vec x,PetscInt ni,const PetscInt ix[],PetscScalar y[])
867: {
872: if (!ni) return(0);
876: (*x->ops->getvalues)(x,ni,ix,y);
877: return(0);
878: }
880: /*@C
881: VecSetValuesBlocked - Inserts or adds blocks of values into certain locations of a vector.
883: Not Collective
885: Input Parameters:
886: + x - vector to insert in
887: . ni - number of blocks to add
888: . ix - indices where to add in block count, rather than element count
889: . y - array of values
890: - iora - either INSERT_VALUES or ADD_VALUES, where
891: ADD_VALUES adds values to any existing entries, and
892: INSERT_VALUES replaces existing entries with new values
894: Notes:
895: VecSetValuesBlocked() sets x[bs*ix[i]+j] = y[bs*i+j],
896: for j=0,...,bs-1, for i=0,...,ni-1. where bs was set with VecSetBlockSize().
898: Calls to VecSetValuesBlocked() with the INSERT_VALUES and ADD_VALUES
899: options cannot be mixed without intervening calls to the assembly
900: routines.
902: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
903: MUST be called after all calls to VecSetValuesBlocked() have been completed.
905: VecSetValuesBlocked() uses 0-based indices in Fortran as well as in C.
907: Negative indices may be passed in ix, these rows are
908: simply ignored. This allows easily inserting element load matrices
909: with homogeneous Dirchlet boundary conditions that you don't want represented
910: in the vector.
912: Level: intermediate
914: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValuesBlockedLocal(),
915: VecSetValues()
916: @*/
917: PetscErrorCode VecSetValuesBlocked(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
918: {
923: if (!ni) return(0);
928: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
929: (*x->ops->setvaluesblocked)(x,ni,ix,y,iora);
930: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
931: PetscObjectStateIncrease((PetscObject)x);
932: return(0);
933: }
935: /*@C
936: VecSetValuesLocal - Inserts or adds values into certain locations of a vector,
937: using a local ordering of the nodes.
939: Not Collective
941: Input Parameters:
942: + x - vector to insert in
943: . ni - number of elements to add
944: . ix - indices where to add
945: . y - array of values
946: - iora - either INSERT_VALUES or ADD_VALUES, where
947: ADD_VALUES adds values to any existing entries, and
948: INSERT_VALUES replaces existing entries with new values
950: Level: intermediate
952: Notes:
953: VecSetValuesLocal() sets x[ix[i]] = y[i], for i=0,...,ni-1.
955: Calls to VecSetValues() with the INSERT_VALUES and ADD_VALUES
956: options cannot be mixed without intervening calls to the assembly
957: routines.
959: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
960: MUST be called after all calls to VecSetValuesLocal() have been completed.
962: VecSetValuesLocal() uses 0-based indices in Fortran as well as in C.
964: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetLocalToGlobalMapping(),
965: VecSetValuesBlockedLocal()
966: @*/
967: PetscErrorCode VecSetValuesLocal(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
968: {
970: PetscInt lixp[128],*lix = lixp;
974: if (!ni) return(0);
979: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
980: if (!x->ops->setvalueslocal) {
981: if (!x->map->mapping) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Local to global never set with VecSetLocalToGlobalMapping()");
982: if (ni > 128) {
983: PetscMalloc1(ni,&lix);
984: }
985: ISLocalToGlobalMappingApply(x->map->mapping,ni,(PetscInt*)ix,lix);
986: (*x->ops->setvalues)(x,ni,lix,y,iora);
987: if (ni > 128) {
988: PetscFree(lix);
989: }
990: } else {
991: (*x->ops->setvalueslocal)(x,ni,ix,y,iora);
992: }
993: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
994: PetscObjectStateIncrease((PetscObject)x);
995: return(0);
996: }
998: /*@
999: VecSetValuesBlockedLocal - Inserts or adds values into certain locations of a vector,
1000: using a local ordering of the nodes.
1002: Not Collective
1004: Input Parameters:
1005: + x - vector to insert in
1006: . ni - number of blocks to add
1007: . ix - indices where to add in block count, not element count
1008: . y - array of values
1009: - iora - either INSERT_VALUES or ADD_VALUES, where
1010: ADD_VALUES adds values to any existing entries, and
1011: INSERT_VALUES replaces existing entries with new values
1013: Level: intermediate
1015: Notes:
1016: VecSetValuesBlockedLocal() sets x[bs*ix[i]+j] = y[bs*i+j],
1017: for j=0,..bs-1, for i=0,...,ni-1, where bs has been set with VecSetBlockSize().
1019: Calls to VecSetValuesBlockedLocal() with the INSERT_VALUES and ADD_VALUES
1020: options cannot be mixed without intervening calls to the assembly
1021: routines.
1023: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
1024: MUST be called after all calls to VecSetValuesBlockedLocal() have been completed.
1026: VecSetValuesBlockedLocal() uses 0-based indices in Fortran as well as in C.
1028: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetValuesBlocked(),
1029: VecSetLocalToGlobalMapping()
1030: @*/
1031: PetscErrorCode VecSetValuesBlockedLocal(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
1032: {
1034: PetscInt lixp[128],*lix = lixp;
1038: if (!ni) return(0);
1042: if (ni > 128) {
1043: PetscMalloc1(ni,&lix);
1044: }
1046: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
1047: ISLocalToGlobalMappingApplyBlock(x->map->mapping,ni,(PetscInt*)ix,lix);
1048: (*x->ops->setvaluesblocked)(x,ni,lix,y,iora);
1049: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
1050: if (ni > 128) {
1051: PetscFree(lix);
1052: }
1053: PetscObjectStateIncrease((PetscObject)x);
1054: return(0);
1055: }
1057: /*@
1058: VecMTDot - Computes indefinite vector multiple dot products.
1059: That is, it does NOT use the complex conjugate.
1061: Collective on Vec
1063: Input Parameters:
1064: + x - one vector
1065: . nv - number of vectors
1066: - y - array of vectors. Note that vectors are pointers
1068: Output Parameter:
1069: . val - array of the dot products
1071: Notes for Users of Complex Numbers:
1072: For complex vectors, VecMTDot() computes the indefinite form
1073: $ val = (x,y) = y^T x,
1074: where y^T denotes the transpose of y.
1076: Use VecMDot() for the inner product
1077: $ val = (x,y) = y^H x,
1078: where y^H denotes the conjugate transpose of y.
1080: Level: intermediate
1082: .seealso: VecMDot(), VecTDot()
1083: @*/
1084: PetscErrorCode VecMTDot(Vec x,PetscInt nv,const Vec y[],PetscScalar val[])
1085: {
1091: if (!nv) return(0);
1098: VecCheckSameSize(x,1,*y,3);
1100: PetscLogEventBegin(VEC_MTDot,x,*y,0,0);
1101: (*x->ops->mtdot)(x,nv,y,val);
1102: PetscLogEventEnd(VEC_MTDot,x,*y,0,0);
1103: return(0);
1104: }
1106: /*@
1107: VecMDot - Computes vector multiple dot products.
1109: Collective on Vec
1111: Input Parameters:
1112: + x - one vector
1113: . nv - number of vectors
1114: - y - array of vectors.
1116: Output Parameter:
1117: . val - array of the dot products (does not allocate the array)
1119: Notes for Users of Complex Numbers:
1120: For complex vectors, VecMDot() computes
1121: $ val = (x,y) = y^H x,
1122: where y^H denotes the conjugate transpose of y.
1124: Use VecMTDot() for the indefinite form
1125: $ val = (x,y) = y^T x,
1126: where y^T denotes the transpose of y.
1128: Level: intermediate
1130: .seealso: VecMTDot(), VecDot()
1131: @*/
1132: PetscErrorCode VecMDot(Vec x,PetscInt nv,const Vec y[],PetscScalar val[])
1133: {
1139: if (!nv) return(0);
1140: if (nv < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Number of vectors (given %D) cannot be negative",nv);
1147: VecCheckSameSize(x,1,*y,3);
1149: PetscLogEventBegin(VEC_MDot,x,*y,0,0);
1150: (*x->ops->mdot)(x,nv,y,val);
1151: PetscLogEventEnd(VEC_MDot,x,*y,0,0);
1152: return(0);
1153: }
1155: /*@
1156: VecMAXPY - Computes y = y + sum alpha[i] x[i]
1158: Logically Collective on Vec
1160: Input Parameters:
1161: + nv - number of scalars and x-vectors
1162: . alpha - array of scalars
1163: . y - one vector
1164: - x - array of vectors
1166: Level: intermediate
1168: Notes:
1169: y cannot be any of the x vectors
1171: .seealso: VecAYPX(), VecWAXPY(), VecAXPY(), VecAXPBYPCZ(), VecAXPBY()
1172: @*/
1173: PetscErrorCode VecMAXPY(Vec y,PetscInt nv,const PetscScalar alpha[],Vec x[])
1174: {
1176: PetscInt i;
1177: PetscBool nonzero;
1182: if (!nv) return(0);
1183: if (nv < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Number of vectors (given %D) cannot be negative",nv);
1190: VecCheckSameSize(y,1,*x,4);
1192: for (i=0, nonzero = PETSC_FALSE; i<nv && !nonzero; i++) nonzero = (PetscBool)(nonzero || alpha[i] != (PetscScalar)0.0);
1193: if (!nonzero) return(0);
1194: VecSetErrorIfLocked(y,1);
1195: PetscLogEventBegin(VEC_MAXPY,*x,y,0,0);
1196: (*y->ops->maxpy)(y,nv,alpha,x);
1197: PetscLogEventEnd(VEC_MAXPY,*x,y,0,0);
1198: PetscObjectStateIncrease((PetscObject)y);
1199: return(0);
1200: }
1202: /*@
1203: VecConcatenate - Creates a new vector that is a vertical concatenation of all the given array of vectors
1204: in the order they appear in the array. The concatenated vector resides on the same
1205: communicator and is the same type as the source vectors.
1207: Collective on X
1209: Input Parameters:
1210: + nx - number of vectors to be concatenated
1211: - X - array containing the vectors to be concatenated in the order of concatenation
1213: Output Parameters:
1214: + Y - concatenated vector
1215: - x_is - array of index sets corresponding to the concatenated components of Y (NULL if not needed)
1217: Notes:
1218: Concatenation is similar to the functionality of a VecNest object; they both represent combination of
1219: different vector spaces. However, concatenated vectors do not store any information about their
1220: sub-vectors and own their own data. Consequently, this function provides index sets to enable the
1221: manipulation of data in the concatenated vector that corresponds to the original components at creation.
1223: This is a useful tool for outer loop algorithms, particularly constrained optimizers, where the solver
1224: has to operate on combined vector spaces and cannot utilize VecNest objects due to incompatibility with
1225: bound projections.
1227: Level: advanced
1229: .seealso: VECNEST, VECSCATTER, VecScatterCreate()
1230: @*/
1231: PetscErrorCode VecConcatenate(PetscInt nx, const Vec X[], Vec *Y, IS *x_is[])
1232: {
1233: MPI_Comm comm;
1234: VecType vec_type;
1235: Vec Ytmp, Xtmp;
1236: IS *is_tmp;
1237: PetscInt i, shift=0, Xnl, Xng, Xbegin;
1246: if ((*X)->ops->concatenate) {
1247: /* use the dedicated concatenation function if available */
1248: (*(*X)->ops->concatenate)(nx,X,Y,x_is);
1249: } else {
1250: /* loop over vectors and start creating IS */
1251: comm = PetscObjectComm((PetscObject)(*X));
1252: VecGetType(*X, &vec_type);
1253: PetscMalloc1(nx, &is_tmp);
1254: for (i=0; i<nx; i++) {
1255: VecGetSize(X[i], &Xng);
1256: VecGetLocalSize(X[i], &Xnl);
1257: VecGetOwnershipRange(X[i], &Xbegin, NULL);
1258: ISCreateStride(comm, Xnl, shift + Xbegin, 1, &is_tmp[i]);
1259: shift += Xng;
1260: }
1261: /* create the concatenated vector */
1262: VecCreate(comm, &Ytmp);
1263: VecSetType(Ytmp, vec_type);
1264: VecSetSizes(Ytmp, PETSC_DECIDE, shift);
1265: VecSetUp(Ytmp);
1266: /* copy data from X array to Y and return */
1267: for (i=0; i<nx; i++) {
1268: VecGetSubVector(Ytmp, is_tmp[i], &Xtmp);
1269: VecCopy(X[i], Xtmp);
1270: VecRestoreSubVector(Ytmp, is_tmp[i], &Xtmp);
1271: }
1272: *Y = Ytmp;
1273: if (x_is) {
1274: *x_is = is_tmp;
1275: } else {
1276: for (i=0; i<nx; i++) {
1277: ISDestroy(&is_tmp[i]);
1278: }
1279: PetscFree(is_tmp);
1280: }
1281: }
1282: return(0);
1283: }
1285: /*@
1286: VecGetSubVector - Gets a vector representing part of another vector
1288: Collective on X and IS
1290: Input Parameters:
1291: + X - vector from which to extract a subvector
1292: - is - index set representing portion of X to extract
1294: Output Parameter:
1295: . Y - subvector corresponding to is
1297: Level: advanced
1299: Notes:
1300: The subvector Y should be returned with VecRestoreSubVector().
1301: X and is must be defined on the same communicator
1303: This function may return a subvector without making a copy, therefore it is not safe to use the original vector while
1304: modifying the subvector. Other non-overlapping subvectors can still be obtained from X using this function.
1305: 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.
1307: .seealso: MatCreateSubMatrix()
1308: @*/
1309: PetscErrorCode VecGetSubVector(Vec X,IS is,Vec *Y)
1310: {
1311: PetscErrorCode ierr;
1312: Vec Z;
1319: if (X->ops->getsubvector) {
1320: (*X->ops->getsubvector)(X,is,&Z);
1321: } else { /* Default implementation currently does no caching */
1322: PetscInt gstart,gend,start;
1323: PetscBool red[2] = { PETSC_TRUE, PETSC_TRUE };
1324: PetscInt n,N,ibs,vbs,bs = -1;
1326: ISGetLocalSize(is,&n);
1327: ISGetSize(is,&N);
1328: ISGetBlockSize(is,&ibs);
1329: VecGetBlockSize(X,&vbs);
1330: VecGetOwnershipRange(X,&gstart,&gend);
1331: ISContiguousLocal(is,gstart,gend,&start,&red[0]);
1332: /* block size is given by IS if ibs > 1; otherwise, check the vector */
1333: if (ibs > 1) {
1334: MPIU_Allreduce(MPI_IN_PLACE,red,1,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)is));
1335: bs = ibs;
1336: } else {
1337: if (n%vbs || vbs == 1) red[1] = PETSC_FALSE; /* this process invalidate the collectiveness of block size */
1338: MPIU_Allreduce(MPI_IN_PLACE,red,2,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)is));
1339: if (red[0] && red[1]) bs = vbs; /* all processes have a valid block size and the access will be contiguous */
1340: }
1341: if (red[0]) { /* We can do a no-copy implementation */
1342: const PetscScalar *x;
1343: PetscInt state = 0;
1344: PetscBool isstd,iscuda,iship;
1346: PetscObjectTypeCompareAny((PetscObject)X,&isstd,VECSEQ,VECMPI,VECSTANDARD,"");
1347: PetscObjectTypeCompareAny((PetscObject)X,&iscuda,VECSEQCUDA,VECMPICUDA,"");
1348: PetscObjectTypeCompareAny((PetscObject)X,&iship,VECSEQHIP,VECMPIHIP,"");
1349: if (iscuda) {
1350: #if defined(PETSC_HAVE_CUDA)
1351: const PetscScalar *x_d;
1352: PetscMPIInt size;
1353: PetscOffloadMask flg;
1355: VecCUDAGetArrays_Private(X,&x,&x_d,&flg);
1356: if (flg == PETSC_OFFLOAD_UNALLOCATED) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Not for PETSC_OFFLOAD_UNALLOCATED");
1357: if (n && !x && !x_d) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Missing vector data");
1358: if (x) x += start;
1359: if (x_d) x_d += start;
1360: MPI_Comm_size(PetscObjectComm((PetscObject)X),&size);
1361: if (size == 1) {
1362: VecCreateSeqCUDAWithArrays(PetscObjectComm((PetscObject)X),bs,n,x,x_d,&Z);
1363: } else {
1364: VecCreateMPICUDAWithArrays(PetscObjectComm((PetscObject)X),bs,n,N,x,x_d,&Z);
1365: }
1366: Z->offloadmask = flg;
1367: #endif
1368: } else if (iship) {
1369: #if defined(PETSC_HAVE_HIP)
1370: const PetscScalar *x_d;
1371: PetscMPIInt size;
1372: PetscOffloadMask flg;
1374: VecHIPGetArrays_Private(X,&x,&x_d,&flg);
1375: if (flg == PETSC_OFFLOAD_UNALLOCATED) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Not for PETSC_OFFLOAD_UNALLOCATED");
1376: if (n && !x && !x_d) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Missing vector data");
1377: if (x) x += start;
1378: if (x_d) x_d += start;
1379: MPI_Comm_size(PetscObjectComm((PetscObject)X),&size);
1380: if (size == 1) {
1381: VecCreateSeqHIPWithArrays(PetscObjectComm((PetscObject)X),bs,n,x,x_d,&Z);
1382: } else {
1383: VecCreateMPIHIPWithArrays(PetscObjectComm((PetscObject)X),bs,n,N,x,x_d,&Z);
1384: }
1385: Z->offloadmask = flg;
1386: #endif
1387: } else if (isstd) {
1388: PetscMPIInt size;
1390: MPI_Comm_size(PetscObjectComm((PetscObject)X),&size);
1391: VecGetArrayRead(X,&x);
1392: if (x) x += start;
1393: if (size == 1) {
1394: VecCreateSeqWithArray(PetscObjectComm((PetscObject)X),bs,n,x,&Z);
1395: } else {
1396: VecCreateMPIWithArray(PetscObjectComm((PetscObject)X),bs,n,N,x,&Z);
1397: }
1398: VecRestoreArrayRead(X,&x);
1399: } else { /* default implementation: use place array */
1400: VecGetArrayRead(X,&x);
1401: VecCreate(PetscObjectComm((PetscObject)X),&Z);
1402: VecSetType(Z,((PetscObject)X)->type_name);
1403: VecSetSizes(Z,n,N);
1404: VecSetBlockSize(Z,bs);
1405: VecPlaceArray(Z,x ? x+start : NULL);
1406: VecRestoreArrayRead(X,&x);
1407: }
1409: /* this is relevant only in debug mode */
1410: VecLockGet(X,&state);
1411: if (state) {
1412: VecLockReadPush(Z);
1413: }
1414: Z->ops->placearray = NULL;
1415: Z->ops->replacearray = NULL;
1416: } else { /* Have to create a scatter and do a copy */
1417: VecScatter scatter;
1419: VecCreate(PetscObjectComm((PetscObject)is),&Z);
1420: VecSetSizes(Z,n,N);
1421: VecSetBlockSize(Z,bs);
1422: VecSetType(Z,((PetscObject)X)->type_name);
1423: VecScatterCreate(X,is,Z,NULL,&scatter);
1424: VecScatterBegin(scatter,X,Z,INSERT_VALUES,SCATTER_FORWARD);
1425: VecScatterEnd(scatter,X,Z,INSERT_VALUES,SCATTER_FORWARD);
1426: PetscObjectCompose((PetscObject)Z,"VecGetSubVector_Scatter",(PetscObject)scatter);
1427: VecScatterDestroy(&scatter);
1428: }
1429: }
1430: /* Record the state when the subvector was gotten so we know whether its values need to be put back */
1431: if (VecGetSubVectorSavedStateId < 0) {PetscObjectComposedDataRegister(&VecGetSubVectorSavedStateId);}
1432: PetscObjectComposedDataSetInt((PetscObject)Z,VecGetSubVectorSavedStateId,1);
1433: *Y = Z;
1434: return(0);
1435: }
1437: /*@
1438: VecRestoreSubVector - Restores a subvector extracted using VecGetSubVector()
1440: Collective on IS
1442: Input Parameters:
1443: + X - vector from which subvector was obtained
1444: . is - index set representing the subset of X
1445: - Y - subvector being restored
1447: Level: advanced
1449: .seealso: VecGetSubVector()
1450: @*/
1451: PetscErrorCode VecRestoreSubVector(Vec X,IS is,Vec *Y)
1452: {
1461: if (X->ops->restoresubvector) {
1462: (*X->ops->restoresubvector)(X,is,Y);
1463: } else {
1464: PETSC_UNUSED PetscObjectState dummystate = 0;
1465: PetscBool valid;
1467: PetscObjectComposedDataGetInt((PetscObject)*Y,VecGetSubVectorSavedStateId,dummystate,valid);
1468: if (!valid) {
1469: VecScatter scatter;
1470: PetscInt state;
1472: VecLockGet(X,&state);
1473: if (state != 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Vec X is locked for read-only or read/write access");
1475: PetscObjectQuery((PetscObject)*Y,"VecGetSubVector_Scatter",(PetscObject*)&scatter);
1476: if (scatter) {
1477: VecScatterBegin(scatter,*Y,X,INSERT_VALUES,SCATTER_REVERSE);
1478: VecScatterEnd(scatter,*Y,X,INSERT_VALUES,SCATTER_REVERSE);
1479: } else {
1480: PetscBool iscuda,iship;
1481: PetscObjectTypeCompareAny((PetscObject)X,&iscuda,VECSEQCUDA,VECMPICUDA,"");
1482: PetscObjectTypeCompareAny((PetscObject)X,&iship,VECSEQHIP,VECMPIHIP,"");
1484: if (iscuda) {
1485: #if defined(PETSC_HAVE_CUDA)
1486: PetscOffloadMask ymask = (*Y)->offloadmask;
1488: /* The offloadmask of X dictates where to move memory
1489: If X GPU data is valid, then move Y data on GPU if needed
1490: Otherwise, move back to the CPU */
1491: switch (X->offloadmask) {
1492: case PETSC_OFFLOAD_BOTH:
1493: if (ymask == PETSC_OFFLOAD_CPU) {
1494: VecCUDAResetArray(*Y);
1495: } else if (ymask == PETSC_OFFLOAD_GPU) {
1496: X->offloadmask = PETSC_OFFLOAD_GPU;
1497: }
1498: break;
1499: case PETSC_OFFLOAD_GPU:
1500: if (ymask == PETSC_OFFLOAD_CPU) {
1501: VecCUDAResetArray(*Y);
1502: }
1503: break;
1504: case PETSC_OFFLOAD_CPU:
1505: if (ymask == PETSC_OFFLOAD_GPU) {
1506: VecResetArray(*Y);
1507: }
1508: break;
1509: case PETSC_OFFLOAD_UNALLOCATED:
1510: case PETSC_OFFLOAD_VECKOKKOS:
1511: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"This should not happen");
1512: }
1513: #endif
1514: } else if (iship) {
1515: #if defined(PETSC_HAVE_HIP)
1516: PetscOffloadMask ymask = (*Y)->offloadmask;
1518: /* The offloadmask of X dictates where to move memory
1519: If X GPU data is valid, then move Y data on GPU if needed
1520: Otherwise, move back to the CPU */
1521: switch (X->offloadmask) {
1522: case PETSC_OFFLOAD_BOTH:
1523: if (ymask == PETSC_OFFLOAD_CPU) {
1524: VecHIPResetArray(*Y);
1525: } else if (ymask == PETSC_OFFLOAD_GPU) {
1526: X->offloadmask = PETSC_OFFLOAD_GPU;
1527: }
1528: break;
1529: case PETSC_OFFLOAD_GPU:
1530: if (ymask == PETSC_OFFLOAD_CPU) {
1531: VecHIPResetArray(*Y);
1532: }
1533: break;
1534: case PETSC_OFFLOAD_CPU:
1535: if (ymask == PETSC_OFFLOAD_GPU) {
1536: VecResetArray(*Y);
1537: }
1538: break;
1539: case PETSC_OFFLOAD_UNALLOCATED:
1540: case PETSC_OFFLOAD_VECKOKKOS:
1541: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"This should not happen");
1542: }
1543: #endif
1544: } else {
1545: /* If OpenCL vecs updated the device memory, this triggers a copy on the CPU */
1546: VecResetArray(*Y);
1547: }
1548: PetscObjectStateIncrease((PetscObject)X);
1549: }
1550: }
1551: VecDestroy(Y);
1552: }
1553: return(0);
1554: }
1556: /*@
1557: VecGetLocalVectorRead - Maps the local portion of a vector into a
1558: vector. You must call VecRestoreLocalVectorRead() when the local
1559: vector is no longer needed.
1561: Not collective.
1563: Input parameter:
1564: . v - The vector for which the local vector is desired.
1566: Output parameter:
1567: . w - Upon exit this contains the local vector.
1569: Level: beginner
1571: Notes:
1572: This function is similar to VecGetArrayRead() which maps the local
1573: portion into a raw pointer. VecGetLocalVectorRead() is usually
1574: almost as efficient as VecGetArrayRead() but in certain circumstances
1575: VecGetLocalVectorRead() can be much more efficient than
1576: VecGetArrayRead(). This is because the construction of a contiguous
1577: array representing the vector data required by VecGetArrayRead() can
1578: be an expensive operation for certain vector types. For example, for
1579: GPU vectors VecGetArrayRead() requires that the data between device
1580: and host is synchronized.
1582: Unlike VecGetLocalVector(), this routine is not collective and
1583: preserves cached information.
1585: .seealso: VecRestoreLocalVectorRead(), VecGetLocalVector(), VecGetArrayRead(), VecGetArray()
1586: @*/
1587: PetscErrorCode VecGetLocalVectorRead(Vec v,Vec w)
1588: {
1590: PetscScalar *a;
1595: VecCheckSameLocalSize(v,1,w,2);
1596: if (v->ops->getlocalvectorread) {
1597: (*v->ops->getlocalvectorread)(v,w);
1598: } else {
1599: VecGetArrayRead(v,(const PetscScalar**)&a);
1600: VecPlaceArray(w,a);
1601: }
1602: PetscObjectStateIncrease((PetscObject)w);
1603: VecLockReadPush(v);
1604: VecLockReadPush(w);
1605: return(0);
1606: }
1608: /*@
1609: VecRestoreLocalVectorRead - Unmaps the local portion of a vector
1610: previously mapped into a vector using VecGetLocalVectorRead().
1612: Not collective.
1614: Input parameter:
1615: + v - The local portion of this vector was previously mapped into w using VecGetLocalVectorRead().
1616: - w - The vector into which the local portion of v was mapped.
1618: Level: beginner
1620: .seealso: VecGetLocalVectorRead(), VecGetLocalVector(), VecGetArrayRead(), VecGetArray()
1621: @*/
1622: PetscErrorCode VecRestoreLocalVectorRead(Vec v,Vec w)
1623: {
1625: PetscScalar *a;
1630: if (v->ops->restorelocalvectorread) {
1631: (*v->ops->restorelocalvectorread)(v,w);
1632: } else {
1633: VecGetArrayRead(w,(const PetscScalar**)&a);
1634: VecRestoreArrayRead(v,(const PetscScalar**)&a);
1635: VecResetArray(w);
1636: }
1637: VecLockReadPop(v);
1638: VecLockReadPop(w);
1639: PetscObjectStateIncrease((PetscObject)w);
1640: return(0);
1641: }
1643: /*@
1644: VecGetLocalVector - Maps the local portion of a vector into a
1645: vector.
1647: Collective on v, not collective on w.
1649: Input parameter:
1650: . v - The vector for which the local vector is desired.
1652: Output parameter:
1653: . w - Upon exit this contains the local vector.
1655: Level: beginner
1657: Notes:
1658: This function is similar to VecGetArray() which maps the local
1659: portion into a raw pointer. VecGetLocalVector() is usually about as
1660: efficient as VecGetArray() but in certain circumstances
1661: VecGetLocalVector() can be much more efficient than VecGetArray().
1662: This is because the construction of a contiguous array representing
1663: the vector data required by VecGetArray() can be an expensive
1664: operation for certain vector types. For example, for GPU vectors
1665: VecGetArray() requires that the data between device and host is
1666: synchronized.
1668: .seealso: VecRestoreLocalVector(), VecGetLocalVectorRead(), VecGetArrayRead(), VecGetArray()
1669: @*/
1670: PetscErrorCode VecGetLocalVector(Vec v,Vec w)
1671: {
1673: PetscScalar *a;
1678: VecCheckSameLocalSize(v,1,w,2);
1679: if (v->ops->getlocalvector) {
1680: (*v->ops->getlocalvector)(v,w);
1681: } else {
1682: VecGetArray(v,&a);
1683: VecPlaceArray(w,a);
1684: }
1685: PetscObjectStateIncrease((PetscObject)w);
1686: return(0);
1687: }
1689: /*@
1690: VecRestoreLocalVector - Unmaps the local portion of a vector
1691: previously mapped into a vector using VecGetLocalVector().
1693: Logically collective.
1695: Input parameter:
1696: + v - The local portion of this vector was previously mapped into w using VecGetLocalVector().
1697: - w - The vector into which the local portion of v was mapped.
1699: Level: beginner
1701: .seealso: VecGetLocalVector(), VecGetLocalVectorRead(), VecRestoreLocalVectorRead(), LocalVectorRead(), VecGetArrayRead(), VecGetArray()
1702: @*/
1703: PetscErrorCode VecRestoreLocalVector(Vec v,Vec w)
1704: {
1706: PetscScalar *a;
1711: if (v->ops->restorelocalvector) {
1712: (*v->ops->restorelocalvector)(v,w);
1713: } else {
1714: VecGetArray(w,&a);
1715: VecRestoreArray(v,&a);
1716: VecResetArray(w);
1717: }
1718: PetscObjectStateIncrease((PetscObject)w);
1719: PetscObjectStateIncrease((PetscObject)v);
1720: return(0);
1721: }
1723: /*@C
1724: VecGetArray - Returns a pointer to a contiguous array that contains this
1725: processor's portion of the vector data. For the standard PETSc
1726: vectors, VecGetArray() returns a pointer to the local data array and
1727: does not use any copies. If the underlying vector data is not stored
1728: in a contiguous array this routine will copy the data to a contiguous
1729: array and return a pointer to that. You MUST call VecRestoreArray()
1730: when you no longer need access to the array.
1732: Logically Collective on Vec
1734: Input Parameter:
1735: . x - the vector
1737: Output Parameter:
1738: . a - location to put pointer to the array
1740: Fortran Note:
1741: This routine is used differently from Fortran 77
1742: $ Vec x
1743: $ PetscScalar x_array(1)
1744: $ PetscOffset i_x
1745: $ PetscErrorCode ierr
1746: $ call VecGetArray(x,x_array,i_x,ierr)
1747: $
1748: $ Access first local entry in vector with
1749: $ value = x_array(i_x + 1)
1750: $
1751: $ ...... other code
1752: $ call VecRestoreArray(x,x_array,i_x,ierr)
1753: For Fortran 90 see VecGetArrayF90()
1755: See the Fortran chapter of the users manual and
1756: petsc/src/snes/tutorials/ex5f.F for details.
1758: Level: beginner
1760: .seealso: VecRestoreArray(), VecGetArrayRead(), VecGetArrays(), VecGetArrayF90(), VecGetArrayReadF90(), VecPlaceArray(), VecGetArray2d(),
1761: VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayWrite(), VecRestoreArrayWrite()
1762: @*/
1763: PetscErrorCode VecGetArray(Vec x,PetscScalar **a)
1764: {
1769: VecSetErrorIfLocked(x,1);
1770: if (x->ops->getarray) { /* The if-else order matters! VECNEST, VECCUDA etc should have ops->getarray while VECCUDA etc are petscnative */
1771: (*x->ops->getarray)(x,a);
1772: } else if (x->petscnative) { /* VECSTANDARD */
1773: *a = *((PetscScalar**)x->data);
1774: } else SETERRQ1(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Cannot get array for vector type \"%s\"",((PetscObject)x)->type_name);
1775: return(0);
1776: }
1778: /*@C
1779: VecRestoreArray - Restores a vector after VecGetArray() has been called.
1781: Logically Collective on Vec
1783: Input Parameters:
1784: + x - the vector
1785: - a - location of pointer to array obtained from VecGetArray()
1787: Level: beginner
1789: .seealso: VecGetArray(), VecRestoreArrayRead(), VecRestoreArrays(), VecRestoreArrayF90(), VecRestoreArrayReadF90(), VecPlaceArray(), VecRestoreArray2d(),
1790: VecGetArrayPair(), VecRestoreArrayPair()
1791: @*/
1792: PetscErrorCode VecRestoreArray(Vec x,PetscScalar **a)
1793: {
1798: if (x->ops->restorearray) { /* VECNEST, VECCUDA etc */
1799: (*x->ops->restorearray)(x,a);
1800: } else if (x->petscnative) { /* VECSTANDARD */
1801: /* nothing */
1802: } else SETERRQ1(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Cannot restore array for vector type \"%s\"",((PetscObject)x)->type_name);
1803: if (a) *a = NULL;
1804: PetscObjectStateIncrease((PetscObject)x);
1805: return(0);
1806: }
1807: /*@C
1808: VecGetArrayRead - Get read-only pointer to contiguous array containing this processor's portion of the vector data.
1810: Not Collective
1812: Input Parameter:
1813: . x - the vector
1815: Output Parameter:
1816: . a - the array
1818: Level: beginner
1820: Notes:
1821: The array must be returned using a matching call to VecRestoreArrayRead().
1823: Unlike VecGetArray(), this routine is not collective and preserves cached information like vector norms.
1825: Standard PETSc vectors use contiguous storage so that this routine does not perform a copy. Other vector
1826: implementations may require a copy, but must such implementations should cache the contiguous representation so that
1827: only one copy is performed when this routine is called multiple times in sequence.
1829: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair()
1830: @*/
1831: PetscErrorCode VecGetArrayRead(Vec x,const PetscScalar **a)
1832: {
1837: if (x->ops->getarray) { /* VECNEST, VECCUDA, VECKOKKOS etc */
1838: (*x->ops->getarray)(x,(PetscScalar**)a);
1839: } else if (x->petscnative) { /* VECSTANDARD */
1840: *a = *((PetscScalar**)x->data);
1841: } else SETERRQ1(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Cannot get array read for vector type \"%s\"",((PetscObject)x)->type_name);
1842: return(0);
1843: }
1845: /*@C
1846: VecRestoreArrayRead - Restore array obtained with VecGetArrayRead()
1848: Not Collective
1850: Input Parameters:
1851: + vec - the vector
1852: - array - the array
1854: Level: beginner
1856: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair()
1857: @*/
1858: PetscErrorCode VecRestoreArrayRead(Vec x,const PetscScalar **a)
1859: {
1864: if (x->petscnative) { /* VECSTANDARD, VECCUDA, VECKOKKOS etc */
1865: /* nothing */
1866: } else if (x->ops->restorearrayread) { /* VECNEST */
1867: (*x->ops->restorearrayread)(x,a);
1868: } else { /* No one? */
1869: (*x->ops->restorearray)(x,(PetscScalar**)a);
1870: }
1871: if (a) *a = NULL;
1872: return(0);
1873: }
1875: /*@C
1876: VecGetArrayWrite - Returns a pointer to a contiguous array that WILL contains this
1877: processor's portion of the vector data. The values in this array are NOT valid, the routine calling this
1878: routine is responsible for putting values into the array; any values it does not set will be invalid
1880: Logically Collective on Vec
1882: Input Parameter:
1883: . x - the vector
1885: Output Parameter:
1886: . a - location to put pointer to the array
1888: Level: intermediate
1890: This is for vectors associate with GPUs, the vector is not copied up before giving access. If you need correct
1891: values in the array use VecGetArray()
1893: Concepts: vector^accessing local values
1895: .seealso: VecRestoreArray(), VecGetArrayRead(), VecGetArrays(), VecGetArrayF90(), VecGetArrayReadF90(), VecPlaceArray(), VecGetArray2d(),
1896: VecGetArrayPair(), VecRestoreArrayPair(), VecGetArray(), VecRestoreArrayWrite()
1897: @*/
1898: PetscErrorCode VecGetArrayWrite(Vec x,PetscScalar **a)
1899: {
1904: VecSetErrorIfLocked(x,1);
1905: if (x->ops->getarraywrite) {
1906: (*x->ops->getarraywrite)(x,a);
1907: } else {
1908: VecGetArray(x,a);
1909: }
1910: return(0);
1911: }
1913: /*@C
1914: VecRestoreArrayWrite - Restores a vector after VecGetArrayWrite() has been called.
1916: Logically Collective on Vec
1918: Input Parameters:
1919: + x - the vector
1920: - a - location of pointer to array obtained from VecGetArray()
1922: Level: beginner
1924: .seealso: VecGetArray(), VecRestoreArrayRead(), VecRestoreArrays(), VecRestoreArrayF90(), VecRestoreArrayReadF90(), VecPlaceArray(), VecRestoreArray2d(),
1925: VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayWrite()
1926: @*/
1927: PetscErrorCode VecRestoreArrayWrite(Vec x,PetscScalar **a)
1928: {
1933: if (x->ops->restorearraywrite) {
1934: (*x->ops->restorearraywrite)(x,a);
1935: } else if (x->ops->restorearray) {
1936: (*x->ops->restorearray)(x,a);
1937: }
1938: if (a) *a = NULL;
1939: PetscObjectStateIncrease((PetscObject)x);
1940: return(0);
1941: }
1943: /*@C
1944: VecGetArrays - Returns a pointer to the arrays in a set of vectors
1945: that were created by a call to VecDuplicateVecs(). You MUST call
1946: VecRestoreArrays() when you no longer need access to the array.
1948: Logically Collective on Vec
1950: Input Parameters:
1951: + x - the vectors
1952: - n - the number of vectors
1954: Output Parameter:
1955: . a - location to put pointer to the array
1957: Fortran Note:
1958: This routine is not supported in Fortran.
1960: Level: intermediate
1962: .seealso: VecGetArray(), VecRestoreArrays()
1963: @*/
1964: PetscErrorCode VecGetArrays(const Vec x[],PetscInt n,PetscScalar **a[])
1965: {
1967: PetscInt i;
1968: PetscScalar **q;
1974: if (n <= 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must get at least one array n = %D",n);
1975: PetscMalloc1(n,&q);
1976: for (i=0; i<n; ++i) {
1977: VecGetArray(x[i],&q[i]);
1978: }
1979: *a = q;
1980: return(0);
1981: }
1983: /*@C
1984: VecRestoreArrays - Restores a group of vectors after VecGetArrays()
1985: has been called.
1987: Logically Collective on Vec
1989: Input Parameters:
1990: + x - the vector
1991: . n - the number of vectors
1992: - a - location of pointer to arrays obtained from VecGetArrays()
1994: Notes:
1995: For regular PETSc vectors this routine does not involve any copies. For
1996: any special vectors that do not store local vector data in a contiguous
1997: array, this routine will copy the data back into the underlying
1998: vector data structure from the arrays obtained with VecGetArrays().
2000: Fortran Note:
2001: This routine is not supported in Fortran.
2003: Level: intermediate
2005: .seealso: VecGetArrays(), VecRestoreArray()
2006: @*/
2007: PetscErrorCode VecRestoreArrays(const Vec x[],PetscInt n,PetscScalar **a[])
2008: {
2010: PetscInt i;
2011: PetscScalar **q = *a;
2018: for (i=0; i<n; ++i) {
2019: VecRestoreArray(x[i],&q[i]);
2020: }
2021: PetscFree(q);
2022: return(0);
2023: }
2025: /*@C
2026: VecGetArrayAndMemType - Like VecGetArray(), but if this is a device vector (e.g., VECCUDA) and the device has up-to-date data,
2027: the returned pointer will be a device pointer to the device memory that contains this processor's portion of the vector data.
2028: Otherwise, when this is a host vector (e.g., VECMPI), or a device vector with the host having newer data than device,
2029: it functions as VecGetArray() and returns a host pointer.
2031: Logically Collective on Vec
2033: Input Parameter:
2034: . x - the vector
2036: Output Parameters:
2037: + a - location to put pointer to the array
2038: - mtype - memory type of the array
2040: Level: beginner
2042: .seealso: VecRestoreArrayAndMemType(), VecRestoreArrayAndMemType(), VecRestoreArray(), VecGetArrayRead(), VecGetArrays(), VecGetArrayF90(), VecGetArrayReadF90(),
2043: VecPlaceArray(), VecGetArray2d(), VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayWrite(), VecRestoreArrayWrite()
2044: @*/
2045: PetscErrorCode VecGetArrayAndMemType(Vec x,PetscScalar **a,PetscMemType *mtype)
2046: {
2052: VecSetErrorIfLocked(x,1);
2053: if (x->ops->getarrayandmemtype) { /* VECCUDA, VECKOKKOS etc */
2054: (*x->ops->getarrayandmemtype)(x,a,mtype);
2055: } else { /* VECSTANDARD, VECNEST, VECVIENNACL */
2056: VecGetArray(x,a);
2057: if (mtype) *mtype = PETSC_MEMTYPE_HOST;
2058: }
2059: return(0);
2060: }
2062: /*@C
2063: VecRestoreArrayAndMemType - Restores a vector after VecGetArrayAndMemType() has been called.
2065: Logically Collective on Vec
2067: Input Parameters:
2068: + x - the vector
2069: - a - location of pointer to array obtained from VecGetArrayAndMemType()
2071: Level: beginner
2073: .seealso: VecGetArrayAndMemType(), VecGetArray(), VecRestoreArrayRead(), VecRestoreArrays(), VecRestoreArrayF90(), VecRestoreArrayReadF90(),
2074: VecPlaceArray(), VecRestoreArray2d(), VecGetArrayPair(), VecRestoreArrayPair()
2075: @*/
2076: PetscErrorCode VecRestoreArrayAndMemType(Vec x,PetscScalar **a)
2077: {
2083: if (x->ops->restorearrayandmemtype) { /* VECCUDA, VECKOKKOS etc */
2084: (*x->ops->restorearrayandmemtype)(x,a);
2085: } else if (x->ops->restorearray) { /* VECNEST, VECVIENNACL */
2086: (*x->ops->restorearray)(x,a);
2087: } /* VECSTANDARD does nothing */
2088: if (a) *a = NULL;
2089: PetscObjectStateIncrease((PetscObject)x);
2090: return(0);
2091: }
2093: /*@C
2094: VecGetArrayReadAndMemType - Like VecGetArrayRead(), but if this is a CUDA vector and it is currently offloaded to GPU,
2095: the returned pointer will be a GPU pointer to the GPU memory that contains this processor's portion of the
2096: vector data. Otherwise, it functions as VecGetArrayRead().
2098: Not Collective
2100: Input Parameter:
2101: . x - the vector
2103: Output Parameters:
2104: + a - the array
2105: - mtype - memory type of the array
2107: Level: beginner
2109: Notes:
2110: The array must be returned using a matching call to VecRestoreArrayReadAndMemType().
2112: .seealso: VecRestoreArrayReadAndMemType(), VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayAndMemType()
2113: @*/
2114: PetscErrorCode VecGetArrayReadAndMemType(Vec x,const PetscScalar **a,PetscMemType *mtype)
2115: {
2121: #if defined(PETSC_USE_DEBUG)
2122: if (x->ops->getarrayreadandmemtype) SETERRQ1(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Not expected vector type \"%s\" has ops->getarrayreadandmemtype",((PetscObject)x)->type_name);
2123: #endif
2125: if (x->ops->getarrayandmemtype) { /* VECCUDA, VECKOKKOS etc, though they are also petscnative */
2126: (*x->ops->getarrayandmemtype)(x,(PetscScalar**)a,mtype);
2127: } else if (x->ops->getarray) { /* VECNEST, VECVIENNACL */
2128: (*x->ops->getarray)(x,(PetscScalar**)a);
2129: if (mtype) *mtype = PETSC_MEMTYPE_HOST;
2130: } else if (x->petscnative) { /* VECSTANDARD */
2131: *a = *((PetscScalar**)x->data);
2132: if (mtype) *mtype = PETSC_MEMTYPE_HOST;
2133: } else SETERRQ1(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Cannot get array read in place for vector type \"%s\"",((PetscObject)x)->type_name);
2134: return(0);
2135: }
2137: /*@C
2138: VecRestoreArrayReadAndMemType - Restore array obtained with VecGetArrayReadAndMemType()
2140: Not Collective
2142: Input Parameters:
2143: + vec - the vector
2144: - array - the array
2146: Level: beginner
2148: .seealso: VecGetArrayReadAndMemType(), VecRestoreArrayAndMemType(), VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair()
2149: @*/
2150: PetscErrorCode VecRestoreArrayReadAndMemType(Vec x,const PetscScalar **a)
2151: {
2157: if (x->petscnative) { /* VECSTANDARD, VECCUDA, VECKOKKOS, VECVIENNACL etc */
2158: /* nothing */
2159: } else if (x->ops->restorearrayread) { /* VECNEST */
2160: (*x->ops->restorearrayread)(x,a);
2161: } else SETERRQ1(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Cannot restore array read in place for vector type \"%s\"",((PetscObject)x)->type_name);
2162: if (a) *a = NULL;
2163: return(0);
2164: }
2166: /*@
2167: VecPlaceArray - Allows one to replace the array in a vector with an
2168: array provided by the user. This is useful to avoid copying an array
2169: into a vector.
2171: Not Collective
2173: Input Parameters:
2174: + vec - the vector
2175: - array - the array
2177: Notes:
2178: You can return to the original array with a call to VecResetArray()
2180: Level: developer
2182: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecResetArray()
2184: @*/
2185: PetscErrorCode VecPlaceArray(Vec vec,const PetscScalar array[])
2186: {
2193: if (vec->ops->placearray) {
2194: (*vec->ops->placearray)(vec,array);
2195: } else SETERRQ(PetscObjectComm((PetscObject)vec),PETSC_ERR_SUP,"Cannot place array in this type of vector");
2196: PetscObjectStateIncrease((PetscObject)vec);
2197: return(0);
2198: }
2200: /*@C
2201: VecReplaceArray - Allows one to replace the array in a vector with an
2202: array provided by the user. This is useful to avoid copying an array
2203: into a vector.
2205: Not Collective
2207: Input Parameters:
2208: + vec - the vector
2209: - array - the array
2211: Notes:
2212: This permanently replaces the array and frees the memory associated
2213: with the old array.
2215: The memory passed in MUST be obtained with PetscMalloc() and CANNOT be
2216: freed by the user. It will be freed when the vector is destroyed.
2218: Not supported from Fortran
2220: Level: developer
2222: .seealso: VecGetArray(), VecRestoreArray(), VecPlaceArray(), VecResetArray()
2224: @*/
2225: PetscErrorCode VecReplaceArray(Vec vec,const PetscScalar array[])
2226: {
2232: if (vec->ops->replacearray) {
2233: (*vec->ops->replacearray)(vec,array);
2234: } else SETERRQ(PetscObjectComm((PetscObject)vec),PETSC_ERR_SUP,"Cannot replace array in this type of vector");
2235: PetscObjectStateIncrease((PetscObject)vec);
2236: return(0);
2237: }
2239: /*@C
2240: VecCUDAGetArray - Provides access to the CUDA buffer inside a vector.
2242: This function has semantics similar to VecGetArray(): the pointer
2243: returned by this function points to a consistent view of the vector
2244: data. This may involve a copy operation of data from the host to the
2245: device if the data on the device is out of date. If the device
2246: memory hasn't been allocated previously it will be allocated as part
2247: of this function call. VecCUDAGetArray() assumes that
2248: the user will modify the vector data. This is similar to
2249: intent(inout) in fortran.
2251: The CUDA device pointer has to be released by calling
2252: VecCUDARestoreArray(). Upon restoring the vector data
2253: the data on the host will be marked as out of date. A subsequent
2254: access of the host data will thus incur a data transfer from the
2255: device to the host.
2257: Input Parameter:
2258: . v - the vector
2260: Output Parameter:
2261: . a - the CUDA device pointer
2263: Fortran note:
2264: This function is not currently available from Fortran.
2266: Level: intermediate
2268: .seealso: VecCUDARestoreArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2269: @*/
2270: PETSC_EXTERN PetscErrorCode VecCUDAGetArray(Vec v, PetscScalar **a)
2271: {
2274: #if defined(PETSC_HAVE_CUDA)
2275: {
2277: VecCUDACopyToGPU(v);
2278: *a = ((Vec_CUDA*)v->spptr)->GPUarray;
2279: }
2280: #endif
2281: return(0);
2282: }
2284: /*@C
2285: VecCUDARestoreArray - Restore a CUDA device pointer previously acquired with VecCUDAGetArray().
2287: This marks the host data as out of date. Subsequent access to the
2288: vector data on the host side with for instance VecGetArray() incurs a
2289: data transfer.
2291: Input Parameters:
2292: + v - the vector
2293: - a - the CUDA device pointer. This pointer is invalid after
2294: VecCUDARestoreArray() returns.
2296: Fortran note:
2297: This function is not currently available from Fortran.
2299: Level: intermediate
2301: .seealso: VecCUDAGetArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2302: @*/
2303: PETSC_EXTERN PetscErrorCode VecCUDARestoreArray(Vec v, PetscScalar **a)
2304: {
2309: #if defined(PETSC_HAVE_CUDA)
2310: v->offloadmask = PETSC_OFFLOAD_GPU;
2311: #endif
2312: PetscObjectStateIncrease((PetscObject)v);
2313: return(0);
2314: }
2316: /*@C
2317: VecCUDAGetArrayRead - Provides read access to the CUDA buffer inside a vector.
2319: This function is analogous to VecGetArrayRead(): The pointer
2320: returned by this function points to a consistent view of the vector
2321: data. This may involve a copy operation of data from the host to the
2322: device if the data on the device is out of date. If the device
2323: memory hasn't been allocated previously it will be allocated as part
2324: of this function call. VecCUDAGetArrayRead() assumes that the
2325: user will not modify the vector data. This is analgogous to
2326: intent(in) in Fortran.
2328: The CUDA device pointer has to be released by calling
2329: VecCUDARestoreArrayRead(). If the data on the host side was
2330: previously up to date it will remain so, i.e. data on both the device
2331: and the host is up to date. Accessing data on the host side does not
2332: incur a device to host data transfer.
2334: Input Parameter:
2335: . v - the vector
2337: Output Parameter:
2338: . a - the CUDA pointer.
2340: Fortran note:
2341: This function is not currently available from Fortran.
2343: Level: intermediate
2345: .seealso: VecCUDARestoreArrayRead(), VecCUDAGetArray(), VecCUDAGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2346: @*/
2347: PETSC_EXTERN PetscErrorCode VecCUDAGetArrayRead(Vec v,const PetscScalar** a)
2348: {
2351: VecCUDAGetArray(v,(PetscScalar**)a);
2352: return(0);
2353: }
2355: /*@C
2356: VecCUDARestoreArrayRead - Restore a CUDA device pointer previously acquired with VecCUDAGetArrayRead().
2358: If the data on the host side was previously up to date it will remain
2359: so, i.e. data on both the device and the host is up to date.
2360: Accessing data on the host side e.g. with VecGetArray() does not
2361: incur a device to host data transfer.
2363: Input Parameters:
2364: + v - the vector
2365: - a - the CUDA device pointer. This pointer is invalid after
2366: VecCUDARestoreArrayRead() returns.
2368: Fortran note:
2369: This function is not currently available from Fortran.
2371: Level: intermediate
2373: .seealso: VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecCUDAGetArray(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2374: @*/
2375: PETSC_EXTERN PetscErrorCode VecCUDARestoreArrayRead(Vec v, const PetscScalar **a)
2376: {
2379: *a = NULL;
2380: return(0);
2381: }
2383: /*@C
2384: VecCUDAGetArrayWrite - Provides write access to the CUDA buffer inside a vector.
2386: The data pointed to by the device pointer is uninitialized. The user
2387: may not read from this data. Furthermore, the entire array needs to
2388: be filled by the user to obtain well-defined behaviour. The device
2389: memory will be allocated by this function if it hasn't been allocated
2390: previously. This is analogous to intent(out) in Fortran.
2392: The device pointer needs to be released with
2393: VecCUDARestoreArrayWrite(). When the pointer is released the
2394: host data of the vector is marked as out of data. Subsequent access
2395: of the host data with e.g. VecGetArray() incurs a device to host data
2396: transfer.
2398: Input Parameter:
2399: . v - the vector
2401: Output Parameter:
2402: . a - the CUDA pointer
2404: Fortran note:
2405: This function is not currently available from Fortran.
2407: Level: advanced
2409: .seealso: VecCUDARestoreArrayWrite(), VecCUDAGetArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2410: @*/
2411: PETSC_EXTERN PetscErrorCode VecCUDAGetArrayWrite(Vec v, PetscScalar **a)
2412: {
2415: #if defined(PETSC_HAVE_CUDA)
2416: {
2418: VecCUDAAllocateCheck(v);
2419: *a = ((Vec_CUDA*)v->spptr)->GPUarray;
2420: }
2421: #endif
2422: return(0);
2423: }
2425: /*@C
2426: VecCUDARestoreArrayWrite - Restore a CUDA device pointer previously acquired with VecCUDAGetArrayWrite().
2428: Data on the host will be marked as out of date. Subsequent access of
2429: the data on the host side e.g. with VecGetArray() will incur a device
2430: to host data transfer.
2432: Input Parameters:
2433: + v - the vector
2434: - a - the CUDA device pointer. This pointer is invalid after
2435: VecCUDARestoreArrayWrite() returns.
2437: Fortran note:
2438: This function is not currently available from Fortran.
2440: Level: intermediate
2442: .seealso: VecCUDAGetArrayWrite(), VecCUDAGetArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2443: @*/
2444: PETSC_EXTERN PetscErrorCode VecCUDARestoreArrayWrite(Vec v, PetscScalar **a)
2445: {
2450: #if defined(PETSC_HAVE_CUDA)
2451: v->offloadmask = PETSC_OFFLOAD_GPU;
2452: if (a) *a = NULL;
2453: #endif
2454: PetscObjectStateIncrease((PetscObject)v);
2455: return(0);
2456: }
2458: /*@C
2459: VecCUDAPlaceArray - Allows one to replace the GPU array in a vector with a
2460: GPU array provided by the user. This is useful to avoid copying an
2461: array into a vector.
2463: Not Collective
2465: Input Parameters:
2466: + vec - the vector
2467: - array - the GPU array
2469: Notes:
2470: You can return to the original GPU array with a call to VecCUDAResetArray()
2471: It is not possible to use VecCUDAPlaceArray() and VecPlaceArray() at the
2472: same time on the same vector.
2474: Level: developer
2476: .seealso: VecPlaceArray(), VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecResetArray(), VecCUDAResetArray(), VecCUDAReplaceArray()
2478: @*/
2479: PetscErrorCode VecCUDAPlaceArray(Vec vin,const PetscScalar a[])
2480: {
2485: #if defined(PETSC_HAVE_CUDA)
2486: VecCUDACopyToGPU(vin);
2487: 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()");
2488: ((Vec_Seq*)vin->data)->unplacedarray = (PetscScalar *) ((Vec_CUDA*)vin->spptr)->GPUarray; /* save previous GPU array so reset can bring it back */
2489: ((Vec_CUDA*)vin->spptr)->GPUarray = (PetscScalar*)a;
2490: vin->offloadmask = PETSC_OFFLOAD_GPU;
2491: #endif
2492: PetscObjectStateIncrease((PetscObject)vin);
2493: return(0);
2494: }
2496: /*@C
2497: VecCUDAReplaceArray - Allows one to replace the GPU array in a vector
2498: with a GPU array provided by the user. This is useful to avoid copying
2499: a GPU array into a vector.
2501: Not Collective
2503: Input Parameters:
2504: + vec - the vector
2505: - array - the GPU array
2507: Notes:
2508: This permanently replaces the GPU array and frees the memory associated
2509: with the old GPU array.
2511: The memory passed in CANNOT be freed by the user. It will be freed
2512: when the vector is destroyed.
2514: Not supported from Fortran
2516: Level: developer
2518: .seealso: VecGetArray(), VecRestoreArray(), VecPlaceArray(), VecResetArray(), VecCUDAResetArray(), VecCUDAPlaceArray(), VecReplaceArray()
2520: @*/
2521: PetscErrorCode VecCUDAReplaceArray(Vec vin,const PetscScalar a[])
2522: {
2523: #if defined(PETSC_HAVE_CUDA)
2524: cudaError_t err;
2525: #endif
2530: #if defined(PETSC_HAVE_CUDA)
2531: if (((Vec_CUDA*)vin->spptr)->GPUarray_allocated) {
2532: err = cudaFree(((Vec_CUDA*)vin->spptr)->GPUarray_allocated);CHKERRCUDA(err);
2533: }
2534: ((Vec_CUDA*)vin->spptr)->GPUarray_allocated = ((Vec_CUDA*)vin->spptr)->GPUarray = (PetscScalar*)a;
2535: vin->offloadmask = PETSC_OFFLOAD_GPU;
2536: #endif
2537: PetscObjectStateIncrease((PetscObject)vin);
2538: return(0);
2539: }
2541: /*@C
2542: VecCUDAResetArray - Resets a vector to use its default memory. Call this
2543: after the use of VecCUDAPlaceArray().
2545: Not Collective
2547: Input Parameters:
2548: . vec - the vector
2550: Level: developer
2552: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecPlaceArray(), VecResetArray(), VecCUDAPlaceArray(), VecCUDAReplaceArray()
2554: @*/
2555: PetscErrorCode VecCUDAResetArray(Vec vin)
2556: {
2561: #if defined(PETSC_HAVE_CUDA)
2562: VecCUDACopyToGPU(vin);
2563: ((Vec_CUDA*)vin->spptr)->GPUarray = (PetscScalar *) ((Vec_Seq*)vin->data)->unplacedarray;
2564: ((Vec_Seq*)vin->data)->unplacedarray = 0;
2565: vin->offloadmask = PETSC_OFFLOAD_GPU;
2566: #endif
2567: PetscObjectStateIncrease((PetscObject)vin);
2568: return(0);
2569: }
2571: /*@C
2572: VecHIPGetArray - Provides access to the HIP buffer inside a vector.
2574: This function has semantics similar to VecGetArray(): the pointer
2575: returned by this function points to a consistent view of the vector
2576: data. This may involve a copy operation of data from the host to the
2577: device if the data on the device is out of date. If the device
2578: memory hasn't been allocated previously it will be allocated as part
2579: of this function call. VecHIPGetArray() assumes that
2580: the user will modify the vector data. This is similar to
2581: intent(inout) in fortran.
2583: The HIP device pointer has to be released by calling
2584: VecHIPRestoreArray(). Upon restoring the vector data
2585: the data on the host will be marked as out of date. A subsequent
2586: access of the host data will thus incur a data transfer from the
2587: device to the host.
2589: Input Parameter:
2590: . v - the vector
2592: Output Parameter:
2593: . a - the HIP device pointer
2595: Fortran note:
2596: This function is not currently available from Fortran.
2598: Level: intermediate
2600: .seealso: VecHIPRestoreArray(), VecHIPGetArrayRead(), VecHIPGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2601: @*/
2602: PETSC_EXTERN PetscErrorCode VecHIPGetArray(Vec v, PetscScalar **a)
2603: {
2604: #if defined(PETSC_HAVE_HIP)
2606: #endif
2610: #if defined(PETSC_HAVE_HIP)
2611: *a = 0;
2612: VecHIPCopyToGPU(v);
2613: *a = ((Vec_HIP*)v->spptr)->GPUarray;
2614: #endif
2615: return(0);
2616: }
2618: /*@C
2619: VecHIPRestoreArray - Restore a HIP device pointer previously acquired with VecHIPGetArray().
2621: This marks the host data as out of date. Subsequent access to the
2622: vector data on the host side with for instance VecGetArray() incurs a
2623: data transfer.
2625: Input Parameters:
2626: + v - the vector
2627: - a - the HIP device pointer. This pointer is invalid after
2628: VecHIPRestoreArray() returns.
2630: Fortran note:
2631: This function is not currently available from Fortran.
2633: Level: intermediate
2635: .seealso: VecHIPGetArray(), VecHIPGetArrayRead(), VecHIPGetArrayWrite(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2636: @*/
2637: PETSC_EXTERN PetscErrorCode VecHIPRestoreArray(Vec v, PetscScalar **a)
2638: {
2643: #if defined(PETSC_HAVE_HIP)
2644: v->offloadmask = PETSC_OFFLOAD_GPU;
2645: #endif
2647: PetscObjectStateIncrease((PetscObject)v);
2648: return(0);
2649: }
2651: /*@C
2652: VecHIPGetArrayRead - Provides read access to the HIP buffer inside a vector.
2654: This function is analogous to VecGetArrayRead(): The pointer
2655: returned by this function points to a consistent view of the vector
2656: data. This may involve a copy operation of data from the host to the
2657: device if the data on the device is out of date. If the device
2658: memory hasn't been allocated previously it will be allocated as part
2659: of this function call. VecHIPGetArrayRead() assumes that the
2660: user will not modify the vector data. This is analgogous to
2661: intent(in) in Fortran.
2663: The HIP device pointer has to be released by calling
2664: VecHIPRestoreArrayRead(). If the data on the host side was
2665: previously up to date it will remain so, i.e. data on both the device
2666: and the host is up to date. Accessing data on the host side does not
2667: incur a device to host data transfer.
2669: Input Parameter:
2670: . v - the vector
2672: Output Parameter:
2673: . a - the HIP pointer.
2675: Fortran note:
2676: This function is not currently available from Fortran.
2678: Level: intermediate
2680: .seealso: VecHIPRestoreArrayRead(), VecHIPGetArray(), VecHIPGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2681: @*/
2682: PETSC_EXTERN PetscErrorCode VecHIPGetArrayRead(Vec v, const PetscScalar **a)
2683: {
2684: #if defined(PETSC_HAVE_HIP)
2686: #endif
2690: #if defined(PETSC_HAVE_HIP)
2691: *a = 0;
2692: VecHIPCopyToGPU(v);
2693: *a = ((Vec_HIP*)v->spptr)->GPUarray;
2694: #endif
2695: return(0);
2696: }
2698: /*@C
2699: VecHIPRestoreArrayRead - Restore a HIP device pointer previously acquired with VecHIPGetArrayRead().
2701: If the data on the host side was previously up to date it will remain
2702: so, i.e. data on both the device and the host is up to date.
2703: Accessing data on the host side e.g. with VecGetArray() does not
2704: incur a device to host data transfer.
2706: Input Parameters:
2707: + v - the vector
2708: - a - the HIP device pointer. This pointer is invalid after
2709: VecHIPRestoreArrayRead() returns.
2711: Fortran note:
2712: This function is not currently available from Fortran.
2714: Level: intermediate
2716: .seealso: VecHIPGetArrayRead(), VecHIPGetArrayWrite(), VecHIPGetArray(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2717: @*/
2718: PETSC_EXTERN PetscErrorCode VecHIPRestoreArrayRead(Vec v, const PetscScalar **a)
2719: {
2722: *a = NULL;
2723: return(0);
2724: }
2726: /*@C
2727: VecHIPGetArrayWrite - Provides write access to the HIP buffer inside a vector.
2729: The data pointed to by the device pointer is uninitialized. The user
2730: may not read from this data. Furthermore, the entire array needs to
2731: be filled by the user to obtain well-defined behaviour. The device
2732: memory will be allocated by this function if it hasn't been allocated
2733: previously. This is analogous to intent(out) in Fortran.
2735: The device pointer needs to be released with
2736: VecHIPRestoreArrayWrite(). When the pointer is released the
2737: host data of the vector is marked as out of data. Subsequent access
2738: of the host data with e.g. VecGetArray() incurs a device to host data
2739: transfer.
2741: Input Parameter:
2742: . v - the vector
2744: Output Parameter:
2745: . a - the HIP pointer
2747: Fortran note:
2748: This function is not currently available from Fortran.
2750: Level: advanced
2752: .seealso: VecHIPRestoreArrayWrite(), VecHIPGetArray(), VecHIPGetArrayRead(), VecHIPGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2753: @*/
2754: PETSC_EXTERN PetscErrorCode VecHIPGetArrayWrite(Vec v, PetscScalar **a)
2755: {
2756: #if defined(PETSC_HAVE_HIP)
2758: #endif
2762: #if defined(PETSC_HAVE_HIP)
2763: *a = 0;
2764: VecHIPAllocateCheck(v);
2765: *a = ((Vec_HIP*)v->spptr)->GPUarray;
2766: #endif
2767: return(0);
2768: }
2770: /*@C
2771: VecHIPRestoreArrayWrite - Restore a HIP device pointer previously acquired with VecHIPGetArrayWrite().
2773: Data on the host will be marked as out of date. Subsequent access of
2774: the data on the host side e.g. with VecGetArray() will incur a device
2775: to host data transfer.
2777: Input Parameters:
2778: + v - the vector
2779: - a - the HIP device pointer. This pointer is invalid after
2780: VecHIPRestoreArrayWrite() returns.
2782: Fortran note:
2783: This function is not currently available from Fortran.
2785: Level: intermediate
2787: .seealso: VecHIPGetArrayWrite(), VecHIPGetArray(), VecHIPGetArrayRead(), VecHIPGetArrayWrite(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2788: @*/
2789: PETSC_EXTERN PetscErrorCode VecHIPRestoreArrayWrite(Vec v, PetscScalar **a)
2790: {
2795: #if defined(PETSC_HAVE_HIP)
2796: v->offloadmask = PETSC_OFFLOAD_GPU;
2797: #endif
2799: PetscObjectStateIncrease((PetscObject)v);
2800: return(0);
2801: }
2803: /*@C
2804: VecHIPPlaceArray - Allows one to replace the GPU array in a vector with a
2805: GPU array provided by the user. This is useful to avoid copying an
2806: array into a vector.
2808: Not Collective
2810: Input Parameters:
2811: + vec - the vector
2812: - array - the GPU array
2814: Notes:
2815: You can return to the original GPU array with a call to VecHIPResetArray()
2816: It is not possible to use VecHIPPlaceArray() and VecPlaceArray() at the
2817: same time on the same vector.
2819: Level: developer
2821: .seealso: VecPlaceArray(), VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecResetArray(), VecHIPResetArray(), VecHIPReplaceArray()
2823: @*/
2824: PetscErrorCode VecHIPPlaceArray(Vec vin,const PetscScalar a[])
2825: {
2830: #if defined(PETSC_HAVE_HIP)
2831: VecHIPCopyToGPU(vin);
2832: 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()");
2833: ((Vec_Seq*)vin->data)->unplacedarray = (PetscScalar *) ((Vec_HIP*)vin->spptr)->GPUarray; /* save previous GPU array so reset can bring it back */
2834: ((Vec_HIP*)vin->spptr)->GPUarray = (PetscScalar*)a;
2835: vin->offloadmask = PETSC_OFFLOAD_GPU;
2836: #endif
2837: PetscObjectStateIncrease((PetscObject)vin);
2838: return(0);
2839: }
2841: /*@C
2842: VecHIPReplaceArray - Allows one to replace the GPU array in a vector
2843: with a GPU array provided by the user. This is useful to avoid copying
2844: a GPU array into a vector.
2846: Not Collective
2848: Input Parameters:
2849: + vec - the vector
2850: - array - the GPU array
2852: Notes:
2853: This permanently replaces the GPU array and frees the memory associated
2854: with the old GPU array.
2856: The memory passed in CANNOT be freed by the user. It will be freed
2857: when the vector is destroyed.
2859: Not supported from Fortran
2861: Level: developer
2863: .seealso: VecGetArray(), VecRestoreArray(), VecPlaceArray(), VecResetArray(), VecHIPResetArray(), VecHIPPlaceArray(), VecReplaceArray()
2865: @*/
2866: PetscErrorCode VecHIPReplaceArray(Vec vin,const PetscScalar a[])
2867: {
2868: #if defined(PETSC_HAVE_HIP)
2869: hipError_t err;
2870: #endif
2875: #if defined(PETSC_HAVE_HIP)
2876: err = hipFree(((Vec_HIP*)vin->spptr)->GPUarray);CHKERRHIP(err);
2877: ((Vec_HIP*)vin->spptr)->GPUarray = (PetscScalar*)a;
2878: vin->offloadmask = PETSC_OFFLOAD_GPU;
2879: #endif
2880: PetscObjectStateIncrease((PetscObject)vin);
2881: return(0);
2882: }
2884: /*@C
2885: VecHIPResetArray - Resets a vector to use its default memory. Call this
2886: after the use of VecHIPPlaceArray().
2888: Not Collective
2890: Input Parameters:
2891: . vec - the vector
2893: Level: developer
2895: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecPlaceArray(), VecResetArray(), VecHIPPlaceArray(), VecHIPReplaceArray()
2897: @*/
2898: PetscErrorCode VecHIPResetArray(Vec vin)
2899: {
2904: #if defined(PETSC_HAVE_HIP)
2905: VecHIPCopyToGPU(vin);
2906: ((Vec_HIP*)vin->spptr)->GPUarray = (PetscScalar *) ((Vec_Seq*)vin->data)->unplacedarray;
2907: ((Vec_Seq*)vin->data)->unplacedarray = 0;
2908: vin->offloadmask = PETSC_OFFLOAD_GPU;
2909: #endif
2910: PetscObjectStateIncrease((PetscObject)vin);
2911: return(0);
2912: }
2914: /*MC
2915: VecDuplicateVecsF90 - Creates several vectors of the same type as an existing vector
2916: and makes them accessible via a Fortran90 pointer.
2918: Synopsis:
2919: VecDuplicateVecsF90(Vec x,PetscInt n,{Vec, pointer :: y(:)},integer ierr)
2921: Collective on Vec
2923: Input Parameters:
2924: + x - a vector to mimic
2925: - n - the number of vectors to obtain
2927: Output Parameters:
2928: + y - Fortran90 pointer to the array of vectors
2929: - ierr - error code
2931: Example of Usage:
2932: .vb
2933: #include <petsc/finclude/petscvec.h>
2934: use petscvec
2936: Vec x
2937: Vec, pointer :: y(:)
2938: ....
2939: call VecDuplicateVecsF90(x,2,y,ierr)
2940: call VecSet(y(2),alpha,ierr)
2941: call VecSet(y(2),alpha,ierr)
2942: ....
2943: call VecDestroyVecsF90(2,y,ierr)
2944: .ve
2946: Notes:
2947: Not yet supported for all F90 compilers
2949: Use VecDestroyVecsF90() to free the space.
2951: Level: beginner
2953: .seealso: VecDestroyVecsF90(), VecDuplicateVecs()
2955: M*/
2957: /*MC
2958: VecRestoreArrayF90 - Restores a vector to a usable state after a call to
2959: VecGetArrayF90().
2961: Synopsis:
2962: VecRestoreArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2964: Logically Collective on Vec
2966: Input Parameters:
2967: + x - vector
2968: - xx_v - the Fortran90 pointer to the array
2970: Output Parameter:
2971: . ierr - error code
2973: Example of Usage:
2974: .vb
2975: #include <petsc/finclude/petscvec.h>
2976: use petscvec
2978: PetscScalar, pointer :: xx_v(:)
2979: ....
2980: call VecGetArrayF90(x,xx_v,ierr)
2981: xx_v(3) = a
2982: call VecRestoreArrayF90(x,xx_v,ierr)
2983: .ve
2985: Level: beginner
2987: .seealso: VecGetArrayF90(), VecGetArray(), VecRestoreArray(), VecRestoreArrayReadF90()
2989: M*/
2991: /*MC
2992: VecDestroyVecsF90 - Frees a block of vectors obtained with VecDuplicateVecsF90().
2994: Synopsis:
2995: VecDestroyVecsF90(PetscInt n,{Vec, pointer :: x(:)},PetscErrorCode ierr)
2997: Collective on Vec
2999: Input Parameters:
3000: + n - the number of vectors previously obtained
3001: - x - pointer to array of vector pointers
3003: Output Parameter:
3004: . ierr - error code
3006: Notes:
3007: Not yet supported for all F90 compilers
3009: Level: beginner
3011: .seealso: VecDestroyVecs(), VecDuplicateVecsF90()
3013: M*/
3015: /*MC
3016: VecGetArrayF90 - Accesses a vector array from Fortran90. For default PETSc
3017: vectors, VecGetArrayF90() returns a pointer to the local data array. Otherwise,
3018: this routine is implementation dependent. You MUST call VecRestoreArrayF90()
3019: when you no longer need access to the array.
3021: Synopsis:
3022: VecGetArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
3024: Logically Collective on Vec
3026: Input Parameter:
3027: . x - vector
3029: Output Parameters:
3030: + xx_v - the Fortran90 pointer to the array
3031: - ierr - error code
3033: Example of Usage:
3034: .vb
3035: #include <petsc/finclude/petscvec.h>
3036: use petscvec
3038: PetscScalar, pointer :: xx_v(:)
3039: ....
3040: call VecGetArrayF90(x,xx_v,ierr)
3041: xx_v(3) = a
3042: call VecRestoreArrayF90(x,xx_v,ierr)
3043: .ve
3045: If you ONLY intend to read entries from the array and not change any entries you should use VecGetArrayReadF90().
3047: Level: beginner
3049: .seealso: VecRestoreArrayF90(), VecGetArray(), VecRestoreArray(), VecGetArrayReadF90()
3051: M*/
3053: /*MC
3054: VecGetArrayReadF90 - Accesses a read only array from Fortran90. For default PETSc
3055: vectors, VecGetArrayF90() returns a pointer to the local data array. Otherwise,
3056: this routine is implementation dependent. You MUST call VecRestoreArrayReadF90()
3057: when you no longer need access to the array.
3059: Synopsis:
3060: VecGetArrayReadF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
3062: Logically Collective on Vec
3064: Input Parameter:
3065: . x - vector
3067: Output Parameters:
3068: + xx_v - the Fortran90 pointer to the array
3069: - ierr - error code
3071: Example of Usage:
3072: .vb
3073: #include <petsc/finclude/petscvec.h>
3074: use petscvec
3076: PetscScalar, pointer :: xx_v(:)
3077: ....
3078: call VecGetArrayReadF90(x,xx_v,ierr)
3079: a = xx_v(3)
3080: call VecRestoreArrayReadF90(x,xx_v,ierr)
3081: .ve
3083: If you intend to write entries into the array you must use VecGetArrayF90().
3085: Level: beginner
3087: .seealso: VecRestoreArrayReadF90(), VecGetArray(), VecRestoreArray(), VecGetArrayRead(), VecRestoreArrayRead(), VecGetArrayF90()
3089: M*/
3091: /*MC
3092: VecRestoreArrayReadF90 - Restores a readonly vector to a usable state after a call to
3093: VecGetArrayReadF90().
3095: Synopsis:
3096: VecRestoreArrayReadF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
3098: Logically Collective on Vec
3100: Input Parameters:
3101: + x - vector
3102: - xx_v - the Fortran90 pointer to the array
3104: Output Parameter:
3105: . ierr - error code
3107: Example of Usage:
3108: .vb
3109: #include <petsc/finclude/petscvec.h>
3110: use petscvec
3112: PetscScalar, pointer :: xx_v(:)
3113: ....
3114: call VecGetArrayReadF90(x,xx_v,ierr)
3115: a = xx_v(3)
3116: call VecRestoreArrayReadF90(x,xx_v,ierr)
3117: .ve
3119: Level: beginner
3121: .seealso: VecGetArrayReadF90(), VecGetArray(), VecRestoreArray(), VecGetArrayRead(), VecRestoreArrayRead(), VecRestoreArrayF90()
3123: M*/
3125: /*@C
3126: VecGetArray2d - Returns a pointer to a 2d contiguous array that contains this
3127: processor's portion of the vector data. You MUST call VecRestoreArray2d()
3128: when you no longer need access to the array.
3130: Logically Collective
3132: Input Parameters:
3133: + x - the vector
3134: . m - first dimension of two dimensional array
3135: . n - second dimension of two dimensional array
3136: . mstart - first index you will use in first coordinate direction (often 0)
3137: - nstart - first index in the second coordinate direction (often 0)
3139: Output Parameter:
3140: . a - location to put pointer to the array
3142: Level: developer
3144: Notes:
3145: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
3146: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3147: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3148: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
3150: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3152: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3153: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3154: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3155: @*/
3156: PetscErrorCode VecGetArray2d(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3157: {
3159: PetscInt i,N;
3160: PetscScalar *aa;
3166: VecGetLocalSize(x,&N);
3167: 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);
3168: VecGetArray(x,&aa);
3170: PetscMalloc1(m,a);
3171: for (i=0; i<m; i++) (*a)[i] = aa + i*n - nstart;
3172: *a -= mstart;
3173: return(0);
3174: }
3176: /*@C
3177: VecGetArray2dWrite - Returns a pointer to a 2d contiguous array that will contain this
3178: processor's portion of the vector data. You MUST call VecRestoreArray2dWrite()
3179: when you no longer need access to the array.
3181: Logically Collective
3183: Input Parameters:
3184: + x - the vector
3185: . m - first dimension of two dimensional array
3186: . n - second dimension of two dimensional array
3187: . mstart - first index you will use in first coordinate direction (often 0)
3188: - nstart - first index in the second coordinate direction (often 0)
3190: Output Parameter:
3191: . a - location to put pointer to the array
3193: Level: developer
3195: Notes:
3196: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
3197: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3198: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3199: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
3201: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3203: Concepts: vector^accessing local values as 2d array
3205: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3206: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3207: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3208: @*/
3209: PetscErrorCode VecGetArray2dWrite(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3210: {
3212: PetscInt i,N;
3213: PetscScalar *aa;
3219: VecGetLocalSize(x,&N);
3220: 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);
3221: VecGetArrayWrite(x,&aa);
3223: PetscMalloc1(m,a);
3224: for (i=0; i<m; i++) (*a)[i] = aa + i*n - nstart;
3225: *a -= mstart;
3226: return(0);
3227: }
3229: /*@C
3230: VecRestoreArray2d - Restores a vector after VecGetArray2d() has been called.
3232: Logically Collective
3234: Input Parameters:
3235: + x - the vector
3236: . m - first dimension of two dimensional array
3237: . n - second dimension of the two dimensional array
3238: . mstart - first index you will use in first coordinate direction (often 0)
3239: . nstart - first index in the second coordinate direction (often 0)
3240: - a - location of pointer to array obtained from VecGetArray2d()
3242: Level: developer
3244: Notes:
3245: For regular PETSc vectors this routine does not involve any copies. For
3246: any special vectors that do not store local vector data in a contiguous
3247: array, this routine will copy the data back into the underlying
3248: vector data structure from the array obtained with VecGetArray().
3250: This routine actually zeros out the a pointer.
3252: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3253: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3254: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3255: @*/
3256: PetscErrorCode VecRestoreArray2d(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3257: {
3259: void *dummy;
3265: dummy = (void*)(*a + mstart);
3266: PetscFree(dummy);
3267: VecRestoreArray(x,NULL);
3268: return(0);
3269: }
3271: /*@C
3272: VecRestoreArray2dWrite - Restores a vector after VecGetArray2dWrite() has been called.
3274: Logically Collective
3276: Input Parameters:
3277: + x - the vector
3278: . m - first dimension of two dimensional array
3279: . n - second dimension of the two dimensional array
3280: . mstart - first index you will use in first coordinate direction (often 0)
3281: . nstart - first index in the second coordinate direction (often 0)
3282: - a - location of pointer to array obtained from VecGetArray2d()
3284: Level: developer
3286: Notes:
3287: For regular PETSc vectors this routine does not involve any copies. For
3288: any special vectors that do not store local vector data in a contiguous
3289: array, this routine will copy the data back into the underlying
3290: vector data structure from the array obtained with VecGetArray().
3292: This routine actually zeros out the a pointer.
3294: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3295: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3296: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3297: @*/
3298: PetscErrorCode VecRestoreArray2dWrite(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3299: {
3301: void *dummy;
3307: dummy = (void*)(*a + mstart);
3308: PetscFree(dummy);
3309: VecRestoreArrayWrite(x,NULL);
3310: return(0);
3311: }
3313: /*@C
3314: VecGetArray1d - Returns a pointer to a 1d contiguous array that contains this
3315: processor's portion of the vector data. You MUST call VecRestoreArray1d()
3316: when you no longer need access to the array.
3318: Logically Collective
3320: Input Parameters:
3321: + x - the vector
3322: . m - first dimension of two dimensional array
3323: - mstart - first index you will use in first coordinate direction (often 0)
3325: Output Parameter:
3326: . a - location to put pointer to the array
3328: Level: developer
3330: Notes:
3331: For a vector obtained from DMCreateLocalVector() mstart are likely
3332: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3333: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
3335: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3337: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3338: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3339: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3340: @*/
3341: PetscErrorCode VecGetArray1d(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
3342: {
3344: PetscInt N;
3350: VecGetLocalSize(x,&N);
3351: if (m != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Local array size %D does not match 1d array dimensions %D",N,m);
3352: VecGetArray(x,a);
3353: *a -= mstart;
3354: return(0);
3355: }
3357: /*@C
3358: VecGetArray1dWrite - Returns a pointer to a 1d contiguous array that will contain this
3359: processor's portion of the vector data. You MUST call VecRestoreArray1dWrite()
3360: when you no longer need access to the array.
3362: Logically Collective
3364: Input Parameters:
3365: + x - the vector
3366: . m - first dimension of two dimensional array
3367: - mstart - first index you will use in first coordinate direction (often 0)
3369: Output Parameter:
3370: . a - location to put pointer to the array
3372: Level: developer
3374: Notes:
3375: For a vector obtained from DMCreateLocalVector() mstart are likely
3376: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3377: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
3379: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3381: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3382: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3383: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3384: @*/
3385: PetscErrorCode VecGetArray1dWrite(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
3386: {
3388: PetscInt N;
3394: VecGetLocalSize(x,&N);
3395: if (m != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Local array size %D does not match 1d array dimensions %D",N,m);
3396: VecGetArrayWrite(x,a);
3397: *a -= mstart;
3398: return(0);
3399: }
3401: /*@C
3402: VecRestoreArray1d - Restores a vector after VecGetArray1d() has been called.
3404: Logically Collective
3406: Input Parameters:
3407: + x - the vector
3408: . m - first dimension of two dimensional array
3409: . mstart - first index you will use in first coordinate direction (often 0)
3410: - a - location of pointer to array obtained from VecGetArray21()
3412: Level: developer
3414: Notes:
3415: For regular PETSc vectors this routine does not involve any copies. For
3416: any special vectors that do not store local vector data in a contiguous
3417: array, this routine will copy the data back into the underlying
3418: vector data structure from the array obtained with VecGetArray1d().
3420: This routine actually zeros out the a pointer.
3422: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3423: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3424: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
3425: @*/
3426: PetscErrorCode VecRestoreArray1d(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
3427: {
3433: VecRestoreArray(x,NULL);
3434: return(0);
3435: }
3437: /*@C
3438: VecRestoreArray1dWrite - Restores a vector after VecGetArray1dWrite() has been called.
3440: Logically Collective
3442: Input Parameters:
3443: + x - the vector
3444: . m - first dimension of two dimensional array
3445: . mstart - first index you will use in first coordinate direction (often 0)
3446: - a - location of pointer to array obtained from VecGetArray21()
3448: Level: developer
3450: Notes:
3451: For regular PETSc vectors this routine does not involve any copies. For
3452: any special vectors that do not store local vector data in a contiguous
3453: array, this routine will copy the data back into the underlying
3454: vector data structure from the array obtained with VecGetArray1d().
3456: This routine actually zeros out the a pointer.
3458: Concepts: vector^accessing local values as 1d array
3460: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3461: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3462: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
3463: @*/
3464: PetscErrorCode VecRestoreArray1dWrite(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
3465: {
3471: VecRestoreArrayWrite(x,NULL);
3472: return(0);
3473: }
3475: /*@C
3476: VecGetArray3d - Returns a pointer to a 3d contiguous array that contains this
3477: processor's portion of the vector data. You MUST call VecRestoreArray3d()
3478: when you no longer need access to the array.
3480: Logically Collective
3482: Input Parameters:
3483: + x - the vector
3484: . m - first dimension of three dimensional array
3485: . n - second dimension of three dimensional array
3486: . p - third dimension of three dimensional array
3487: . mstart - first index you will use in first coordinate direction (often 0)
3488: . nstart - first index in the second coordinate direction (often 0)
3489: - pstart - first index in the third coordinate direction (often 0)
3491: Output Parameter:
3492: . a - location to put pointer to the array
3494: Level: developer
3496: Notes:
3497: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3498: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3499: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3500: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3502: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3504: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3505: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3506: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3507: @*/
3508: PetscErrorCode VecGetArray3d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3509: {
3511: PetscInt i,N,j;
3512: PetscScalar *aa,**b;
3518: VecGetLocalSize(x,&N);
3519: 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);
3520: VecGetArray(x,&aa);
3522: PetscMalloc1(m*sizeof(PetscScalar**)+m*n,a);
3523: b = (PetscScalar**)((*a) + m);
3524: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3525: for (i=0; i<m; i++)
3526: for (j=0; j<n; j++)
3527: b[i*n+j] = aa + i*n*p + j*p - pstart;
3529: *a -= mstart;
3530: return(0);
3531: }
3533: /*@C
3534: VecGetArray3dWrite - Returns a pointer to a 3d contiguous array that will contain this
3535: processor's portion of the vector data. You MUST call VecRestoreArray3dWrite()
3536: when you no longer need access to the array.
3538: Logically Collective
3540: Input Parameters:
3541: + x - the vector
3542: . m - first dimension of three dimensional array
3543: . n - second dimension of three dimensional array
3544: . p - third dimension of three dimensional array
3545: . mstart - first index you will use in first coordinate direction (often 0)
3546: . nstart - first index in the second coordinate direction (often 0)
3547: - pstart - first index in the third coordinate direction (often 0)
3549: Output Parameter:
3550: . a - location to put pointer to the array
3552: Level: developer
3554: Notes:
3555: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3556: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3557: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3558: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3560: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3562: Concepts: vector^accessing local values as 3d array
3564: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3565: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3566: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3567: @*/
3568: PetscErrorCode VecGetArray3dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3569: {
3571: PetscInt i,N,j;
3572: PetscScalar *aa,**b;
3578: VecGetLocalSize(x,&N);
3579: 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);
3580: VecGetArrayWrite(x,&aa);
3582: PetscMalloc1(m*sizeof(PetscScalar**)+m*n,a);
3583: b = (PetscScalar**)((*a) + m);
3584: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3585: for (i=0; i<m; i++)
3586: for (j=0; j<n; j++)
3587: b[i*n+j] = aa + i*n*p + j*p - pstart;
3589: *a -= mstart;
3590: return(0);
3591: }
3593: /*@C
3594: VecRestoreArray3d - Restores a vector after VecGetArray3d() has been called.
3596: Logically Collective
3598: Input Parameters:
3599: + x - the vector
3600: . m - first dimension of three dimensional array
3601: . n - second dimension of the three dimensional array
3602: . p - third dimension of the three dimensional array
3603: . mstart - first index you will use in first coordinate direction (often 0)
3604: . nstart - first index in the second coordinate direction (often 0)
3605: . pstart - first index in the third coordinate direction (often 0)
3606: - a - location of pointer to array obtained from VecGetArray3d()
3608: Level: developer
3610: Notes:
3611: For regular PETSc vectors this routine does not involve any copies. For
3612: any special vectors that do not store local vector data in a contiguous
3613: array, this routine will copy the data back into the underlying
3614: vector data structure from the array obtained with VecGetArray().
3616: This routine actually zeros out the a pointer.
3618: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3619: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3620: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3621: @*/
3622: PetscErrorCode VecRestoreArray3d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3623: {
3625: void *dummy;
3631: dummy = (void*)(*a + mstart);
3632: PetscFree(dummy);
3633: VecRestoreArray(x,NULL);
3634: return(0);
3635: }
3637: /*@C
3638: VecRestoreArray3dWrite - Restores a vector after VecGetArray3dWrite() has been called.
3640: Logically Collective
3642: Input Parameters:
3643: + x - the vector
3644: . m - first dimension of three dimensional array
3645: . n - second dimension of the three dimensional array
3646: . p - third dimension of the three dimensional array
3647: . mstart - first index you will use in first coordinate direction (often 0)
3648: . nstart - first index in the second coordinate direction (often 0)
3649: . pstart - first index in the third coordinate direction (often 0)
3650: - a - location of pointer to array obtained from VecGetArray3d()
3652: Level: developer
3654: Notes:
3655: For regular PETSc vectors this routine does not involve any copies. For
3656: any special vectors that do not store local vector data in a contiguous
3657: array, this routine will copy the data back into the underlying
3658: vector data structure from the array obtained with VecGetArray().
3660: This routine actually zeros out the a pointer.
3662: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3663: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3664: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3665: @*/
3666: PetscErrorCode VecRestoreArray3dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3667: {
3669: void *dummy;
3675: dummy = (void*)(*a + mstart);
3676: PetscFree(dummy);
3677: VecRestoreArrayWrite(x,NULL);
3678: return(0);
3679: }
3681: /*@C
3682: VecGetArray4d - Returns a pointer to a 4d contiguous array that contains this
3683: processor's portion of the vector data. You MUST call VecRestoreArray4d()
3684: when you no longer need access to the array.
3686: Logically Collective
3688: Input Parameters:
3689: + x - the vector
3690: . m - first dimension of four dimensional array
3691: . n - second dimension of four dimensional array
3692: . p - third dimension of four dimensional array
3693: . q - fourth dimension of four dimensional array
3694: . mstart - first index you will use in first coordinate direction (often 0)
3695: . nstart - first index in the second coordinate direction (often 0)
3696: . pstart - first index in the third coordinate direction (often 0)
3697: - qstart - first index in the fourth coordinate direction (often 0)
3699: Output Parameter:
3700: . a - location to put pointer to the array
3702: Level: beginner
3704: Notes:
3705: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3706: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3707: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3708: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3710: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3712: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3713: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3714: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3715: @*/
3716: PetscErrorCode VecGetArray4d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3717: {
3719: PetscInt i,N,j,k;
3720: PetscScalar *aa,***b,**c;
3726: VecGetLocalSize(x,&N);
3727: 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);
3728: VecGetArray(x,&aa);
3730: PetscMalloc1(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p,a);
3731: b = (PetscScalar***)((*a) + m);
3732: c = (PetscScalar**)(b + m*n);
3733: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3734: for (i=0; i<m; i++)
3735: for (j=0; j<n; j++)
3736: b[i*n+j] = c + i*n*p + j*p - pstart;
3737: for (i=0; i<m; i++)
3738: for (j=0; j<n; j++)
3739: for (k=0; k<p; k++)
3740: c[i*n*p+j*p+k] = aa + i*n*p*q + j*p*q + k*q - qstart;
3741: *a -= mstart;
3742: return(0);
3743: }
3745: /*@C
3746: VecGetArray4dWrite - Returns a pointer to a 4d contiguous array that will contain this
3747: processor's portion of the vector data. You MUST call VecRestoreArray4dWrite()
3748: when you no longer need access to the array.
3750: Logically Collective
3752: Input Parameters:
3753: + x - the vector
3754: . m - first dimension of four dimensional array
3755: . n - second dimension of four dimensional array
3756: . p - third dimension of four dimensional array
3757: . q - fourth dimension of four dimensional array
3758: . mstart - first index you will use in first coordinate direction (often 0)
3759: . nstart - first index in the second coordinate direction (often 0)
3760: . pstart - first index in the third coordinate direction (often 0)
3761: - qstart - first index in the fourth coordinate direction (often 0)
3763: Output Parameter:
3764: . a - location to put pointer to the array
3766: Level: beginner
3768: Notes:
3769: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3770: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3771: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3772: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3774: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3776: Concepts: vector^accessing local values as 3d array
3778: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3779: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3780: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3781: @*/
3782: PetscErrorCode VecGetArray4dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3783: {
3785: PetscInt i,N,j,k;
3786: PetscScalar *aa,***b,**c;
3792: VecGetLocalSize(x,&N);
3793: 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);
3794: VecGetArrayWrite(x,&aa);
3796: PetscMalloc1(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p,a);
3797: b = (PetscScalar***)((*a) + m);
3798: c = (PetscScalar**)(b + m*n);
3799: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3800: for (i=0; i<m; i++)
3801: for (j=0; j<n; j++)
3802: b[i*n+j] = c + i*n*p + j*p - pstart;
3803: for (i=0; i<m; i++)
3804: for (j=0; j<n; j++)
3805: for (k=0; k<p; k++)
3806: c[i*n*p+j*p+k] = aa + i*n*p*q + j*p*q + k*q - qstart;
3807: *a -= mstart;
3808: return(0);
3809: }
3811: /*@C
3812: VecRestoreArray4d - Restores a vector after VecGetArray3d() has been called.
3814: Logically Collective
3816: Input Parameters:
3817: + x - the vector
3818: . m - first dimension of four dimensional array
3819: . n - second dimension of the four dimensional array
3820: . p - third dimension of the four dimensional array
3821: . q - fourth dimension of the four dimensional array
3822: . mstart - first index you will use in first coordinate direction (often 0)
3823: . nstart - first index in the second coordinate direction (often 0)
3824: . pstart - first index in the third coordinate direction (often 0)
3825: . qstart - first index in the fourth coordinate direction (often 0)
3826: - a - location of pointer to array obtained from VecGetArray4d()
3828: Level: beginner
3830: Notes:
3831: For regular PETSc vectors this routine does not involve any copies. For
3832: any special vectors that do not store local vector data in a contiguous
3833: array, this routine will copy the data back into the underlying
3834: vector data structure from the array obtained with VecGetArray().
3836: This routine actually zeros out the a pointer.
3838: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3839: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3840: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3841: @*/
3842: PetscErrorCode VecRestoreArray4d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3843: {
3845: void *dummy;
3851: dummy = (void*)(*a + mstart);
3852: PetscFree(dummy);
3853: VecRestoreArray(x,NULL);
3854: return(0);
3855: }
3857: /*@C
3858: VecRestoreArray4dWrite - Restores a vector after VecGetArray3dWrite() has been called.
3860: Logically Collective
3862: Input Parameters:
3863: + x - the vector
3864: . m - first dimension of four dimensional array
3865: . n - second dimension of the four dimensional array
3866: . p - third dimension of the four dimensional array
3867: . q - fourth dimension of the four dimensional array
3868: . mstart - first index you will use in first coordinate direction (often 0)
3869: . nstart - first index in the second coordinate direction (often 0)
3870: . pstart - first index in the third coordinate direction (often 0)
3871: . qstart - first index in the fourth coordinate direction (often 0)
3872: - a - location of pointer to array obtained from VecGetArray4d()
3874: Level: beginner
3876: Notes:
3877: For regular PETSc vectors this routine does not involve any copies. For
3878: any special vectors that do not store local vector data in a contiguous
3879: array, this routine will copy the data back into the underlying
3880: vector data structure from the array obtained with VecGetArray().
3882: This routine actually zeros out the a pointer.
3884: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3885: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3886: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3887: @*/
3888: PetscErrorCode VecRestoreArray4dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3889: {
3891: void *dummy;
3897: dummy = (void*)(*a + mstart);
3898: PetscFree(dummy);
3899: VecRestoreArrayWrite(x,NULL);
3900: return(0);
3901: }
3903: /*@C
3904: VecGetArray2dRead - Returns a pointer to a 2d contiguous array that contains this
3905: processor's portion of the vector data. You MUST call VecRestoreArray2dRead()
3906: when you no longer need access to the array.
3908: Logically Collective
3910: Input Parameters:
3911: + x - the vector
3912: . m - first dimension of two dimensional array
3913: . n - second dimension of two dimensional array
3914: . mstart - first index you will use in first coordinate direction (often 0)
3915: - nstart - first index in the second coordinate direction (often 0)
3917: Output Parameter:
3918: . a - location to put pointer to the array
3920: Level: developer
3922: Notes:
3923: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
3924: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3925: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3926: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
3928: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3930: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3931: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3932: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3933: @*/
3934: PetscErrorCode VecGetArray2dRead(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3935: {
3936: PetscErrorCode ierr;
3937: PetscInt i,N;
3938: const PetscScalar *aa;
3944: VecGetLocalSize(x,&N);
3945: 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);
3946: VecGetArrayRead(x,&aa);
3948: PetscMalloc1(m,a);
3949: for (i=0; i<m; i++) (*a)[i] = (PetscScalar*) aa + i*n - nstart;
3950: *a -= mstart;
3951: return(0);
3952: }
3954: /*@C
3955: VecRestoreArray2dRead - Restores a vector after VecGetArray2dRead() has been called.
3957: Logically Collective
3959: Input Parameters:
3960: + x - the vector
3961: . m - first dimension of two dimensional array
3962: . n - second dimension of the two dimensional array
3963: . mstart - first index you will use in first coordinate direction (often 0)
3964: . nstart - first index in the second coordinate direction (often 0)
3965: - a - location of pointer to array obtained from VecGetArray2d()
3967: Level: developer
3969: Notes:
3970: For regular PETSc vectors this routine does not involve any copies. For
3971: any special vectors that do not store local vector data in a contiguous
3972: array, this routine will copy the data back into the underlying
3973: vector data structure from the array obtained with VecGetArray().
3975: This routine actually zeros out the a pointer.
3977: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3978: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3979: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3980: @*/
3981: PetscErrorCode VecRestoreArray2dRead(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3982: {
3984: void *dummy;
3990: dummy = (void*)(*a + mstart);
3991: PetscFree(dummy);
3992: VecRestoreArrayRead(x,NULL);
3993: return(0);
3994: }
3996: /*@C
3997: VecGetArray1dRead - Returns a pointer to a 1d contiguous array that contains this
3998: processor's portion of the vector data. You MUST call VecRestoreArray1dRead()
3999: when you no longer need access to the array.
4001: Logically Collective
4003: Input Parameters:
4004: + x - the vector
4005: . m - first dimension of two dimensional array
4006: - mstart - first index you will use in first coordinate direction (often 0)
4008: Output Parameter:
4009: . a - location to put pointer to the array
4011: Level: developer
4013: Notes:
4014: For a vector obtained from DMCreateLocalVector() mstart are likely
4015: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
4016: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
4018: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
4020: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
4021: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
4022: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
4023: @*/
4024: PetscErrorCode VecGetArray1dRead(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
4025: {
4027: PetscInt N;
4033: VecGetLocalSize(x,&N);
4034: if (m != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Local array size %D does not match 1d array dimensions %D",N,m);
4035: VecGetArrayRead(x,(const PetscScalar**)a);
4036: *a -= mstart;
4037: return(0);
4038: }
4040: /*@C
4041: VecRestoreArray1dRead - Restores a vector after VecGetArray1dRead() has been called.
4043: Logically Collective
4045: Input Parameters:
4046: + x - the vector
4047: . m - first dimension of two dimensional array
4048: . mstart - first index you will use in first coordinate direction (often 0)
4049: - a - location of pointer to array obtained from VecGetArray21()
4051: Level: developer
4053: Notes:
4054: For regular PETSc vectors this routine does not involve any copies. For
4055: any special vectors that do not store local vector data in a contiguous
4056: array, this routine will copy the data back into the underlying
4057: vector data structure from the array obtained with VecGetArray1dRead().
4059: This routine actually zeros out the a pointer.
4061: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
4062: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
4063: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
4064: @*/
4065: PetscErrorCode VecRestoreArray1dRead(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
4066: {
4072: VecRestoreArrayRead(x,NULL);
4073: return(0);
4074: }
4076: /*@C
4077: VecGetArray3dRead - Returns a pointer to a 3d contiguous array that contains this
4078: processor's portion of the vector data. You MUST call VecRestoreArray3dRead()
4079: when you no longer need access to the array.
4081: Logically Collective
4083: Input Parameters:
4084: + x - the vector
4085: . m - first dimension of three dimensional array
4086: . n - second dimension of three dimensional array
4087: . p - third dimension of three dimensional array
4088: . mstart - first index you will use in first coordinate direction (often 0)
4089: . nstart - first index in the second coordinate direction (often 0)
4090: - pstart - first index in the third coordinate direction (often 0)
4092: Output Parameter:
4093: . a - location to put pointer to the array
4095: Level: developer
4097: Notes:
4098: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
4099: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
4100: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
4101: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3dRead().
4103: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
4105: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
4106: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
4107: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
4108: @*/
4109: PetscErrorCode VecGetArray3dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
4110: {
4111: PetscErrorCode ierr;
4112: PetscInt i,N,j;
4113: const PetscScalar *aa;
4114: PetscScalar **b;
4120: VecGetLocalSize(x,&N);
4121: 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);
4122: VecGetArrayRead(x,&aa);
4124: PetscMalloc1(m*sizeof(PetscScalar**)+m*n,a);
4125: b = (PetscScalar**)((*a) + m);
4126: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
4127: for (i=0; i<m; i++)
4128: for (j=0; j<n; j++)
4129: b[i*n+j] = (PetscScalar *)aa + i*n*p + j*p - pstart;
4131: *a -= mstart;
4132: return(0);
4133: }
4135: /*@C
4136: VecRestoreArray3dRead - Restores a vector after VecGetArray3dRead() has been called.
4138: Logically Collective
4140: Input Parameters:
4141: + x - the vector
4142: . m - first dimension of three dimensional array
4143: . n - second dimension of the three dimensional array
4144: . p - third dimension of the three dimensional array
4145: . mstart - first index you will use in first coordinate direction (often 0)
4146: . nstart - first index in the second coordinate direction (often 0)
4147: . pstart - first index in the third coordinate direction (often 0)
4148: - a - location of pointer to array obtained from VecGetArray3dRead()
4150: Level: developer
4152: Notes:
4153: For regular PETSc vectors this routine does not involve any copies. For
4154: any special vectors that do not store local vector data in a contiguous
4155: array, this routine will copy the data back into the underlying
4156: vector data structure from the array obtained with VecGetArray().
4158: This routine actually zeros out the a pointer.
4160: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
4161: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
4162: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
4163: @*/
4164: PetscErrorCode VecRestoreArray3dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
4165: {
4167: void *dummy;
4173: dummy = (void*)(*a + mstart);
4174: PetscFree(dummy);
4175: VecRestoreArrayRead(x,NULL);
4176: return(0);
4177: }
4179: /*@C
4180: VecGetArray4dRead - Returns a pointer to a 4d contiguous array that contains this
4181: processor's portion of the vector data. You MUST call VecRestoreArray4dRead()
4182: when you no longer need access to the array.
4184: Logically Collective
4186: Input Parameters:
4187: + x - the vector
4188: . m - first dimension of four dimensional array
4189: . n - second dimension of four dimensional array
4190: . p - third dimension of four dimensional array
4191: . q - fourth dimension of four dimensional array
4192: . mstart - first index you will use in first coordinate direction (often 0)
4193: . nstart - first index in the second coordinate direction (often 0)
4194: . pstart - first index in the third coordinate direction (often 0)
4195: - qstart - first index in the fourth coordinate direction (often 0)
4197: Output Parameter:
4198: . a - location to put pointer to the array
4200: Level: beginner
4202: Notes:
4203: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
4204: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
4205: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
4206: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
4208: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
4210: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
4211: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
4212: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
4213: @*/
4214: PetscErrorCode VecGetArray4dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
4215: {
4216: PetscErrorCode ierr;
4217: PetscInt i,N,j,k;
4218: const PetscScalar *aa;
4219: PetscScalar ***b,**c;
4225: VecGetLocalSize(x,&N);
4226: 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);
4227: VecGetArrayRead(x,&aa);
4229: PetscMalloc1(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p,a);
4230: b = (PetscScalar***)((*a) + m);
4231: c = (PetscScalar**)(b + m*n);
4232: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
4233: for (i=0; i<m; i++)
4234: for (j=0; j<n; j++)
4235: b[i*n+j] = c + i*n*p + j*p - pstart;
4236: for (i=0; i<m; i++)
4237: for (j=0; j<n; j++)
4238: for (k=0; k<p; k++)
4239: c[i*n*p+j*p+k] = (PetscScalar*) aa + i*n*p*q + j*p*q + k*q - qstart;
4240: *a -= mstart;
4241: return(0);
4242: }
4244: /*@C
4245: VecRestoreArray4dRead - Restores a vector after VecGetArray3d() has been called.
4247: Logically Collective
4249: Input Parameters:
4250: + x - the vector
4251: . m - first dimension of four dimensional array
4252: . n - second dimension of the four dimensional array
4253: . p - third dimension of the four dimensional array
4254: . q - fourth dimension of the four dimensional array
4255: . mstart - first index you will use in first coordinate direction (often 0)
4256: . nstart - first index in the second coordinate direction (often 0)
4257: . pstart - first index in the third coordinate direction (often 0)
4258: . qstart - first index in the fourth coordinate direction (often 0)
4259: - a - location of pointer to array obtained from VecGetArray4dRead()
4261: Level: beginner
4263: Notes:
4264: For regular PETSc vectors this routine does not involve any copies. For
4265: any special vectors that do not store local vector data in a contiguous
4266: array, this routine will copy the data back into the underlying
4267: vector data structure from the array obtained with VecGetArray().
4269: This routine actually zeros out the a pointer.
4271: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
4272: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
4273: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
4274: @*/
4275: PetscErrorCode VecRestoreArray4dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
4276: {
4278: void *dummy;
4284: dummy = (void*)(*a + mstart);
4285: PetscFree(dummy);
4286: VecRestoreArrayRead(x,NULL);
4287: return(0);
4288: }
4290: #if defined(PETSC_USE_DEBUG)
4292: /*@
4293: VecLockGet - Gets the current lock status of a vector
4295: Logically Collective on Vec
4297: Input Parameter:
4298: . x - the vector
4300: Output Parameter:
4301: . state - greater than zero indicates the vector is locked for read; less then zero indicates the vector is
4302: locked for write; equal to zero means the vector is unlocked, that is, it is free to read or write.
4304: Level: beginner
4306: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPush(), VecLockReadPop()
4307: @*/
4308: PetscErrorCode VecLockGet(Vec x,PetscInt *state)
4309: {
4312: *state = x->lock;
4313: return(0);
4314: }
4316: /*@
4317: VecLockReadPush - Pushes a read-only lock on a vector to prevent it from writing
4319: Logically Collective on Vec
4321: Input Parameter:
4322: . x - the vector
4324: Notes:
4325: If this is set then calls to VecGetArray() or VecSetValues() or any other routines that change the vectors values will fail.
4327: The call can be nested, i.e., called multiple times on the same vector, but each VecLockReadPush(x) has to have one matching
4328: VecLockReadPop(x), which removes the latest read-only lock.
4330: Level: beginner
4332: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPop(), VecLockGet()
4333: @*/
4334: PetscErrorCode VecLockReadPush(Vec x)
4335: {
4338: 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");
4339: x->lock++;
4340: return(0);
4341: }
4343: /*@
4344: VecLockReadPop - Pops a read-only lock from a vector
4346: Logically Collective on Vec
4348: Input Parameter:
4349: . x - the vector
4351: Level: beginner
4353: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPush(), VecLockGet()
4354: @*/
4355: PetscErrorCode VecLockReadPop(Vec x)
4356: {
4359: x->lock--;
4360: if (x->lock < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Vector has been unlocked from read-only access too many times");
4361: return(0);
4362: }
4364: /*@C
4365: VecLockWriteSet_Private - Lock or unlock a vector for exclusive read/write access
4367: Logically Collective on Vec
4369: Input Parameters:
4370: + x - the vector
4371: - flg - PETSC_TRUE to lock the vector for writing; PETSC_FALSE to unlock it.
4373: Notes:
4374: The function is usefull in split-phase computations, which usually have a begin phase and an end phase.
4375: One can call VecLockWriteSet_Private(x,PETSC_TRUE) in the begin phase to lock a vector for exclusive
4376: access, and call VecLockWriteSet_Private(x,PETSC_FALSE) in the end phase to unlock the vector from exclusive
4377: access. In this way, one is ensured no other operations can access the vector in between. The code may like
4379: VecGetArray(x,&xdata); // begin phase
4380: VecLockWriteSet_Private(v,PETSC_TRUE);
4382: Other operations, which can not acceess x anymore (they can access xdata, of course)
4384: VecRestoreArray(x,&vdata); // end phase
4385: VecLockWriteSet_Private(v,PETSC_FALSE);
4387: The call can not be nested on the same vector, in other words, one can not call VecLockWriteSet_Private(x,PETSC_TRUE)
4388: again before calling VecLockWriteSet_Private(v,PETSC_FALSE).
4390: Level: beginner
4392: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPush(), VecLockReadPop(), VecLockGet()
4393: @*/
4394: PetscErrorCode VecLockWriteSet_Private(Vec x,PetscBool flg)
4395: {
4398: if (flg) {
4399: 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");
4400: 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");
4401: else x->lock = -1;
4402: } else {
4403: 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");
4404: x->lock = 0;
4405: }
4406: return(0);
4407: }
4409: /*@
4410: VecLockPush - Pushes a read-only lock on a vector to prevent it from writing
4412: Level: deprecated
4414: .seealso: VecLockReadPush()
4415: @*/
4416: PetscErrorCode VecLockPush(Vec x)
4417: {
4420: VecLockReadPush(x);
4421: return(0);
4422: }
4424: /*@
4425: VecLockPop - Pops a read-only lock from a vector
4427: Level: deprecated
4429: .seealso: VecLockReadPop()
4430: @*/
4431: PetscErrorCode VecLockPop(Vec x)
4432: {
4435: VecLockReadPop(x);
4436: return(0);
4437: }
4439: #endif