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: PetscInt VecGetSubVectorSavedStateId = -1;
18: PETSC_EXTERN PetscErrorCode VecValidValues(Vec vec,PetscInt argnum,PetscBool begin)
19: {
20: #if defined(PETSC_USE_DEBUG)
21: PetscInt n,i;
22: const PetscScalar *x;
24: #if defined(PETSC_HAVE_DEVICE)
25: if ((vec->petscnative || vec->ops->getarray) && (vec->offloadmask & PETSC_OFFLOAD_CPU)) {
26: #else
27: if (vec->petscnative || vec->ops->getarray) {
28: #endif
29: VecGetLocalSize(vec,&n);
30: VecGetArrayRead(vec,&x);
31: for (i=0; i<n; i++) {
32: if (begin) {
34: } else {
36: }
37: }
38: VecRestoreArrayRead(vec,&x);
39: }
40: #else
41: #endif
42: return 0;
43: }
45: /*@
46: VecMaxPointwiseDivide - Computes the maximum of the componentwise division max = max_i abs(x_i/y_i).
48: Logically Collective on Vec
50: Input Parameters:
51: . x, y - the vectors
53: Output Parameter:
54: . max - the result
56: Level: advanced
58: Notes:
59: x and y may be the same vector
60: if a particular y_i is zero, it is treated as 1 in the above formula
62: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMax(), VecPointwiseMin(), VecPointwiseMaxAbs()
63: @*/
64: PetscErrorCode VecMaxPointwiseDivide(Vec x,Vec y,PetscReal *max)
65: {
72: VecCheckSameSize(x,1,y,2);
73: (*x->ops->maxpointwisedivide)(x,y,max);
74: return 0;
75: }
77: /*@
78: VecDot - Computes the vector dot product.
80: Collective on Vec
82: Input Parameters:
83: . x, y - the vectors
85: Output Parameter:
86: . val - the dot product
88: Performance Issues:
89: $ per-processor memory bandwidth
90: $ interprocessor latency
91: $ work load imbalance that causes certain processes to arrive much earlier than others
93: Notes for Users of Complex Numbers:
94: For complex vectors, VecDot() computes
95: $ val = (x,y) = y^H x,
96: where y^H denotes the conjugate transpose of y. Note that this corresponds to the usual "mathematicians" complex
97: inner product where the SECOND argument gets the complex conjugate. Since the BLASdot() complex conjugates the first
98: first argument we call the BLASdot() with the arguments reversed.
100: Use VecTDot() for the indefinite form
101: $ val = (x,y) = y^T x,
102: where y^T denotes the transpose of y.
104: Level: intermediate
106: .seealso: VecMDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecDotRealPart()
107: @*/
108: PetscErrorCode VecDot(Vec x,Vec y,PetscScalar *val)
109: {
116: VecCheckSameSize(x,1,y,2);
118: PetscLogEventBegin(VEC_Dot,x,y,0,0);
119: (*x->ops->dot)(x,y,val);
120: PetscLogEventEnd(VEC_Dot,x,y,0,0);
121: return 0;
122: }
124: /*@
125: VecDotRealPart - Computes the real part of the vector dot product.
127: Collective on Vec
129: Input Parameters:
130: . x, y - the vectors
132: Output Parameter:
133: . val - the real part of the dot product;
135: Performance Issues:
136: $ per-processor memory bandwidth
137: $ interprocessor latency
138: $ work load imbalance that causes certain processes to arrive much earlier than others
140: Notes for Users of Complex Numbers:
141: See VecDot() for more details on the definition of the dot product for complex numbers
143: For real numbers this returns the same value as VecDot()
145: 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
146: the space R^{2n} (that is a vector of 2n components with the real or imaginary part of the complex numbers for components)
148: Developer Note: This is not currently optimized to compute only the real part of the dot product.
150: Level: intermediate
152: .seealso: VecMDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecDot(), VecDotNorm2()
153: @*/
154: PetscErrorCode VecDotRealPart(Vec x,Vec y,PetscReal *val)
155: {
156: PetscScalar fdot;
158: VecDot(x,y,&fdot);
159: *val = PetscRealPart(fdot);
160: return 0;
161: }
163: /*@
164: VecNorm - Computes the vector norm.
166: Collective on Vec
168: Input Parameters:
169: + x - the vector
170: - type - the type of the norm requested
172: Output Parameter:
173: . val - the norm
175: Values of NormType:
176: + NORM_1 - sum_i |x_i|
177: . NORM_2 - sqrt(sum_i |x_i|^2)
178: . NORM_INFINITY - max_i |x_i|
179: - NORM_1_AND_2 - computes efficiently both NORM_1 and NORM_2 and stores them each in an output array
181: Notes:
182: For complex numbers NORM_1 will return the traditional 1 norm of the 2 norm of the complex numbers; that is the 1
183: norm of the absolute values of the complex entries. In PETSc 3.6 and earlier releases it returned the 1 norm of
184: the 1 norm of the complex entries (what is returned by the BLAS routine asum()). Both are valid norms but most
185: people expect the former.
187: This routine stashes the computed norm value, repeated calls before the vector entries are changed are then rapid since the
188: precomputed value is immediately available. Certain vector operations such as VecSet() store the norms so the value is
189: immediately available and does not need to be explicitly computed. VecScale() updates any stashed norm values, thus calls after VecScale()
190: do not need to explicitly recompute the norm.
192: Level: intermediate
194: Performance Issues:
195: + per-processor memory bandwidth - limits the speed of the computation of local portion of the norm
196: . interprocessor latency - limits the accumulation of the result across ranks, .i.e. MPI_Allreduce() time
197: . number of ranks - the time for the result will grow with the log base 2 of the number of ranks sharing the vector
198: - work load imbalance - the rank with the largest number of vector entries will limit the speed up
200: .seealso: VecDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecNormAvailable(),
201: VecNormBegin(), VecNormEnd(), NormType()
203: @*/
204: PetscErrorCode VecNorm(Vec x,NormType type,PetscReal *val)
205: {
206: PetscBool flg;
212: /*
213: * Cached data?
214: */
215: if (type!=NORM_1_AND_2) {
216: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[type],*val,flg);
217: if (flg) return 0;
218: }
219: PetscLogEventBegin(VEC_Norm,x,0,0,0);
220: (*x->ops->norm)(x,type,val);
221: PetscLogEventEnd(VEC_Norm,x,0,0,0);
222: if (type!=NORM_1_AND_2) {
223: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[type],*val);
224: }
225: return 0;
226: }
228: /*@
229: VecNormAvailable - Returns the vector norm if it is already known.
231: Not Collective
233: Input Parameters:
234: + x - the vector
235: - type - one of NORM_1, NORM_2, NORM_INFINITY. Also available
236: NORM_1_AND_2, which computes both norms and stores them
237: in a two element array.
239: Output Parameters:
240: + available - PETSC_TRUE if the val returned is valid
241: - val - the norm
243: Notes:
244: $ NORM_1 denotes sum_i |x_i|
245: $ NORM_2 denotes sqrt(sum_i (x_i)^2)
246: $ NORM_INFINITY denotes max_i |x_i|
248: Level: intermediate
250: Performance Issues:
251: $ per-processor memory bandwidth
252: $ interprocessor latency
253: $ work load imbalance that causes certain processes to arrive much earlier than others
255: Compile Option:
256: PETSC_HAVE_SLOW_BLAS_NORM2 will cause a C (loop unrolled) version of the norm to be used, rather
257: than the BLAS. This should probably only be used when one is using the FORTRAN BLAS routines
258: (as opposed to vendor provided) because the FORTRAN BLAS NRM2() routine is very slow.
260: .seealso: VecDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecNorm()
261: VecNormBegin(), VecNormEnd()
263: @*/
264: PetscErrorCode VecNormAvailable(Vec x,NormType type,PetscBool *available,PetscReal *val)
265: {
270: *available = PETSC_FALSE;
271: if (type!=NORM_1_AND_2) {
272: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[type],*val,*available);
273: }
274: return 0;
275: }
277: /*@
278: VecNormalize - Normalizes a vector by 2-norm.
280: Collective on Vec
282: Input Parameter:
283: . x - the vector
285: Output Parameter:
286: . val - the vector norm before normalization
288: Level: intermediate
290: @*/
291: PetscErrorCode VecNormalize(Vec x,PetscReal *val)
292: {
293: PetscReal norm;
297: PetscLogEventBegin(VEC_Normalize,x,0,0,0);
298: VecNorm(x,NORM_2,&norm);
299: if (norm == 0.0) {
300: PetscInfo(x,"Vector of zero norm can not be normalized; Returning only the zero norm\n");
301: } else if (norm != 1.0) {
302: PetscScalar tmp = 1.0/norm;
303: VecScale(x,tmp);
304: }
305: if (val) *val = norm;
306: PetscLogEventEnd(VEC_Normalize,x,0,0,0);
307: return 0;
308: }
310: /*@C
311: VecMax - Determines the vector component with maximum real part and its location.
313: Collective on Vec
315: Input Parameter:
316: . x - the vector
318: Output Parameters:
319: + p - the location of val (pass NULL if you don't want this)
320: - val - the maximum component
322: Notes:
323: Returns the value PETSC_MIN_REAL and negative p if the vector is of length 0.
325: Returns the smallest index with the maximum value
326: Level: intermediate
328: .seealso: VecNorm(), VecMin()
329: @*/
330: PetscErrorCode VecMax(Vec x,PetscInt *p,PetscReal *val)
331: {
335: PetscLogEventBegin(VEC_Max,x,0,0,0);
336: (*x->ops->max)(x,p,val);
337: PetscLogEventEnd(VEC_Max,x,0,0,0);
338: return 0;
339: }
341: /*@C
342: VecMin - Determines the vector component with minimum real part and its location.
344: Collective on Vec
346: Input Parameter:
347: . x - the vector
349: Output Parameters:
350: + p - the location of val (pass NULL if you don't want this location)
351: - val - the minimum component
353: Level: intermediate
355: Notes:
356: Returns the value PETSC_MAX_REAL and negative p if the vector is of length 0.
358: This returns the smallest index with the minumum value
360: .seealso: VecMax()
361: @*/
362: PetscErrorCode VecMin(Vec x,PetscInt *p,PetscReal *val)
363: {
367: PetscLogEventBegin(VEC_Min,x,0,0,0);
368: (*x->ops->min)(x,p,val);
369: PetscLogEventEnd(VEC_Min,x,0,0,0);
370: return 0;
371: }
373: /*@
374: VecTDot - Computes an indefinite vector dot product. That is, this
375: routine does NOT use the complex conjugate.
377: Collective on Vec
379: Input Parameters:
380: . x, y - the vectors
382: Output Parameter:
383: . val - the dot product
385: Notes for Users of Complex Numbers:
386: For complex vectors, VecTDot() computes the indefinite form
387: $ val = (x,y) = y^T x,
388: where y^T denotes the transpose of y.
390: Use VecDot() for the inner product
391: $ val = (x,y) = y^H x,
392: where y^H denotes the conjugate transpose of y.
394: Level: intermediate
396: .seealso: VecDot(), VecMTDot()
397: @*/
398: PetscErrorCode VecTDot(Vec x,Vec y,PetscScalar *val)
399: {
406: VecCheckSameSize(x,1,y,2);
408: PetscLogEventBegin(VEC_TDot,x,y,0,0);
409: (*x->ops->tdot)(x,y,val);
410: PetscLogEventEnd(VEC_TDot,x,y,0,0);
411: return 0;
412: }
414: /*@
415: VecScale - Scales a vector.
417: Not collective on Vec
419: Input Parameters:
420: + x - the vector
421: - alpha - the scalar
423: Note:
424: For a vector with n components, VecScale() computes
425: $ x[i] = alpha * x[i], for i=1,...,n.
427: Level: intermediate
429: @*/
430: PetscErrorCode VecScale(Vec x, PetscScalar alpha)
431: {
432: PetscReal norms[4] = {0.0,0.0,0.0, 0.0};
433: PetscBool flgs[4];
434: PetscInt i;
439: PetscLogEventBegin(VEC_Scale,x,0,0,0);
440: if (alpha != (PetscScalar)1.0) {
441: VecSetErrorIfLocked(x,1);
442: /* get current stashed norms */
443: for (i=0; i<4; i++) {
444: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[i],norms[i],flgs[i]);
445: }
446: (*x->ops->scale)(x,alpha);
447: PetscObjectStateIncrease((PetscObject)x);
448: /* put the scaled stashed norms back into the Vec */
449: for (i=0; i<4; i++) {
450: if (flgs[i]) {
451: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[i],PetscAbsScalar(alpha)*norms[i]);
452: }
453: }
454: }
455: PetscLogEventEnd(VEC_Scale,x,0,0,0);
456: return 0;
457: }
459: /*@
460: VecSet - Sets all components of a vector to a single scalar value.
462: Logically Collective on Vec
464: Input Parameters:
465: + x - the vector
466: - alpha - the scalar
468: Output Parameter:
469: . x - the vector
471: Note:
472: For a vector of dimension n, VecSet() computes
473: $ x[i] = alpha, for i=1,...,n,
474: so that all vector entries then equal the identical
475: scalar value, alpha. Use the more general routine
476: VecSetValues() to set different vector entries.
478: You CANNOT call this after you have called VecSetValues() but before you call
479: VecAssemblyBegin/End().
481: Level: beginner
483: .seealso VecSetValues(), VecSetValuesBlocked(), VecSetRandom()
485: @*/
486: PetscErrorCode VecSet(Vec x,PetscScalar alpha)
487: {
488: PetscReal val;
494: VecSetErrorIfLocked(x,1);
496: PetscLogEventBegin(VEC_Set,x,0,0,0);
497: (*x->ops->set)(x,alpha);
498: PetscLogEventEnd(VEC_Set,x,0,0,0);
499: PetscObjectStateIncrease((PetscObject)x);
501: /* norms can be simply set (if |alpha|*N not too large) */
502: val = PetscAbsScalar(alpha);
503: if (x->map->N == 0) {
504: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_1],0.0l);
505: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],0.0);
506: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_2],0.0);
507: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_FROBENIUS],0.0);
508: } else if (val > PETSC_MAX_REAL/x->map->N) {
509: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],val);
510: } else {
511: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_1],x->map->N * val);
512: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],val);
513: val = PetscSqrtReal((PetscReal)x->map->N) * val;
514: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_2],val);
515: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_FROBENIUS],val);
516: }
517: return 0;
518: }
520: /*@
521: VecAXPY - Computes y = alpha x + y.
523: Logically Collective on Vec
525: Input Parameters:
526: + alpha - the scalar
527: - x, y - the vectors
529: Output Parameter:
530: . y - output vector
532: Level: intermediate
534: Notes:
535: x and y MUST be different vectors
536: This routine is optimized for alpha of 0.0, otherwise it calls the BLAS routine
538: $ VecAXPY(y,alpha,x) y = alpha x + y
539: $ VecAYPX(y,beta,x) y = x + beta y
540: $ VecAXPBY(y,alpha,beta,x) y = alpha x + beta y
541: $ VecWAXPY(w,alpha,x,y) w = alpha x + y
542: $ VecAXPBYPCZ(w,alpha,beta,gamma,x,y) z = alpha x + beta y + gamma z
543: $ VecMAXPY(y,nv,alpha[],x[]) y = sum alpha[i] x[i] + y
545: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPBYPCZ(), VecAXPBY()
546: @*/
547: PetscErrorCode VecAXPY(Vec y,PetscScalar alpha,Vec x)
548: {
554: VecCheckSameSize(x,3,y,1);
557: if (alpha == (PetscScalar)0.0) return 0;
558: VecSetErrorIfLocked(y,1);
560: VecLockReadPush(x);
561: PetscLogEventBegin(VEC_AXPY,x,y,0,0);
562: (*y->ops->axpy)(y,alpha,x);
563: PetscLogEventEnd(VEC_AXPY,x,y,0,0);
564: VecLockReadPop(x);
565: PetscObjectStateIncrease((PetscObject)y);
566: return 0;
567: }
569: /*@
570: VecAXPBY - Computes y = alpha x + beta y.
572: Logically Collective on Vec
574: Input Parameters:
575: + alpha,beta - the scalars
576: - x, y - the vectors
578: Output Parameter:
579: . y - output vector
581: Level: intermediate
583: Notes:
584: x and y MUST be different vectors
585: The implementation is optimized for alpha and/or beta values of 0.0 and 1.0
587: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPY(), VecAXPBYPCZ()
588: @*/
589: PetscErrorCode VecAXPBY(Vec y,PetscScalar alpha,PetscScalar beta,Vec x)
590: {
596: VecCheckSameSize(y,1,x,4);
600: if (alpha == (PetscScalar)0.0 && beta == (PetscScalar)1.0) return 0;
601: VecSetErrorIfLocked(y,1);
602: PetscLogEventBegin(VEC_AXPY,x,y,0,0);
603: (*y->ops->axpby)(y,alpha,beta,x);
604: PetscLogEventEnd(VEC_AXPY,x,y,0,0);
605: PetscObjectStateIncrease((PetscObject)y);
606: return 0;
607: }
609: /*@
610: VecAXPBYPCZ - Computes z = alpha x + beta y + gamma z
612: Logically Collective on Vec
614: Input Parameters:
615: + alpha,beta, gamma - the scalars
616: - x, y, z - the vectors
618: Output Parameter:
619: . z - output vector
621: Level: intermediate
623: Notes:
624: x, y and z must be different vectors
625: The implementation is optimized for alpha of 1.0 and gamma of 1.0 or 0.0
627: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPY(), VecAXPBY()
628: @*/
629: PetscErrorCode VecAXPBYPCZ(Vec z,PetscScalar alpha,PetscScalar beta,PetscScalar gamma,Vec x,Vec y)
630: {
639: VecCheckSameSize(x,1,y,5);
640: VecCheckSameSize(x,1,z,6);
646: if (alpha == (PetscScalar)0.0 && beta == (PetscScalar)0.0 && gamma == (PetscScalar)1.0) return 0;
647: VecSetErrorIfLocked(z,1);
649: PetscLogEventBegin(VEC_AXPBYPCZ,x,y,z,0);
650: (*y->ops->axpbypcz)(z,alpha,beta,gamma,x,y);
651: PetscLogEventEnd(VEC_AXPBYPCZ,x,y,z,0);
652: PetscObjectStateIncrease((PetscObject)z);
653: return 0;
654: }
656: /*@
657: VecAYPX - Computes y = x + beta y.
659: Logically Collective on Vec
661: Input Parameters:
662: + beta - the scalar
663: - x, y - the vectors
665: Output Parameter:
666: . y - output vector
668: Level: intermediate
670: Notes:
671: x and y MUST be different vectors
672: The implementation is optimized for beta of -1.0, 0.0, and 1.0
674: .seealso: VecMAXPY(), VecWAXPY(), VecAXPY(), VecAXPBYPCZ(), VecAXPBY()
675: @*/
676: PetscErrorCode VecAYPX(Vec y,PetscScalar beta,Vec x)
677: {
683: VecCheckSameSize(x,1,y,3);
686: VecSetErrorIfLocked(y,1);
688: PetscLogEventBegin(VEC_AYPX,x,y,0,0);
689: (*y->ops->aypx)(y,beta,x);
690: PetscLogEventEnd(VEC_AYPX,x,y,0,0);
691: PetscObjectStateIncrease((PetscObject)y);
692: return 0;
693: }
695: /*@
696: VecWAXPY - Computes w = alpha x + y.
698: Logically Collective on Vec
700: Input Parameters:
701: + alpha - the scalar
702: - x, y - the vectors
704: Output Parameter:
705: . w - the result
707: Level: intermediate
709: Notes:
710: w cannot be either x or y, but x and y can be the same
711: The implementation is optimzed for alpha of -1.0, 0.0, and 1.0
713: .seealso: VecAXPY(), VecAYPX(), VecAXPBY(), VecMAXPY(), VecAXPBYPCZ()
714: @*/
715: PetscErrorCode VecWAXPY(Vec w,PetscScalar alpha,Vec x,Vec y)
716: {
725: VecCheckSameSize(x,3,y,4);
726: VecCheckSameSize(x,3,w,1);
730: VecSetErrorIfLocked(w,1);
732: PetscLogEventBegin(VEC_WAXPY,x,y,w,0);
733: (*w->ops->waxpy)(w,alpha,x,y);
734: PetscLogEventEnd(VEC_WAXPY,x,y,w,0);
735: PetscObjectStateIncrease((PetscObject)w);
736: return 0;
737: }
739: /*@C
740: VecSetValues - Inserts or adds values into certain locations of a vector.
742: Not Collective
744: Input Parameters:
745: + x - vector to insert in
746: . ni - number of elements to add
747: . ix - indices where to add
748: . y - array of values
749: - iora - either INSERT_VALUES or ADD_VALUES, where
750: ADD_VALUES adds values to any existing entries, and
751: INSERT_VALUES replaces existing entries with new values
753: Notes:
754: VecSetValues() sets x[ix[i]] = y[i], for i=0,...,ni-1.
756: Calls to VecSetValues() with the INSERT_VALUES and ADD_VALUES
757: options cannot be mixed without intervening calls to the assembly
758: routines.
760: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
761: MUST be called after all calls to VecSetValues() have been completed.
763: VecSetValues() uses 0-based indices in Fortran as well as in C.
765: If you call VecSetOption(x, VEC_IGNORE_NEGATIVE_INDICES,PETSC_TRUE),
766: negative indices may be passed in ix. These rows are
767: simply ignored. This allows easily inserting element load matrices
768: with homogeneous Dirchlet boundary conditions that you don't want represented
769: in the vector.
771: Level: beginner
773: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValuesLocal(),
774: VecSetValue(), VecSetValuesBlocked(), InsertMode, INSERT_VALUES, ADD_VALUES, VecGetValues()
775: @*/
776: PetscErrorCode VecSetValues(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
777: {
780: if (!ni) return 0;
785: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
786: (*x->ops->setvalues)(x,ni,ix,y,iora);
787: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
788: PetscObjectStateIncrease((PetscObject)x);
789: return 0;
790: }
792: /*@C
793: VecGetValues - Gets values from certain locations of a vector. Currently
794: can only get values on the same processor
796: Not Collective
798: Input Parameters:
799: + x - vector to get values from
800: . ni - number of elements to get
801: - ix - indices where to get them from (in global 1d numbering)
803: Output Parameter:
804: . y - array of values
806: Notes:
807: The user provides the allocated array y; it is NOT allocated in this routine
809: VecGetValues() gets y[i] = x[ix[i]], for i=0,...,ni-1.
811: VecAssemblyBegin() and VecAssemblyEnd() MUST be called before calling this
813: VecGetValues() uses 0-based indices in Fortran as well as in C.
815: If you call VecSetOption(x, VEC_IGNORE_NEGATIVE_INDICES,PETSC_TRUE),
816: negative indices may be passed in ix. These rows are
817: simply ignored.
819: Level: beginner
821: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues()
822: @*/
823: PetscErrorCode VecGetValues(Vec x,PetscInt ni,const PetscInt ix[],PetscScalar y[])
824: {
826: if (!ni) return 0;
830: (*x->ops->getvalues)(x,ni,ix,y);
831: return 0;
832: }
834: /*@C
835: VecSetValuesBlocked - Inserts or adds blocks of values into certain locations of a vector.
837: Not Collective
839: Input Parameters:
840: + x - vector to insert in
841: . ni - number of blocks to add
842: . ix - indices where to add in block count, rather than element count
843: . y - array of values
844: - iora - either INSERT_VALUES or ADD_VALUES, where
845: ADD_VALUES adds values to any existing entries, and
846: INSERT_VALUES replaces existing entries with new values
848: Notes:
849: VecSetValuesBlocked() sets x[bs*ix[i]+j] = y[bs*i+j],
850: for j=0,...,bs-1, for i=0,...,ni-1. where bs was set with VecSetBlockSize().
852: Calls to VecSetValuesBlocked() with the INSERT_VALUES and ADD_VALUES
853: options cannot be mixed without intervening calls to the assembly
854: routines.
856: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
857: MUST be called after all calls to VecSetValuesBlocked() have been completed.
859: VecSetValuesBlocked() uses 0-based indices in Fortran as well as in C.
861: Negative indices may be passed in ix, these rows are
862: simply ignored. This allows easily inserting element load matrices
863: with homogeneous Dirchlet boundary conditions that you don't want represented
864: in the vector.
866: Level: intermediate
868: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValuesBlockedLocal(),
869: VecSetValues()
870: @*/
871: PetscErrorCode VecSetValuesBlocked(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
872: {
875: if (!ni) return 0;
880: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
881: (*x->ops->setvaluesblocked)(x,ni,ix,y,iora);
882: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
883: PetscObjectStateIncrease((PetscObject)x);
884: return 0;
885: }
887: /*@C
888: VecSetValuesLocal - Inserts or adds values into certain locations of a vector,
889: using a local ordering of the nodes.
891: Not Collective
893: Input Parameters:
894: + x - vector to insert in
895: . ni - number of elements to add
896: . ix - indices where to add
897: . y - array of values
898: - iora - either INSERT_VALUES or ADD_VALUES, where
899: ADD_VALUES adds values to any existing entries, and
900: INSERT_VALUES replaces existing entries with new values
902: Level: intermediate
904: Notes:
905: VecSetValuesLocal() sets x[ix[i]] = y[i], for i=0,...,ni-1.
907: Calls to VecSetValues() with the INSERT_VALUES and ADD_VALUES
908: options cannot be mixed without intervening calls to the assembly
909: routines.
911: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
912: MUST be called after all calls to VecSetValuesLocal() have been completed.
914: VecSetValuesLocal() uses 0-based indices in Fortran as well as in C.
916: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetLocalToGlobalMapping(),
917: VecSetValuesBlockedLocal()
918: @*/
919: PetscErrorCode VecSetValuesLocal(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
920: {
921: PetscInt lixp[128],*lix = lixp;
925: if (!ni) return 0;
930: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
931: if (!x->ops->setvalueslocal) {
932: if (x->map->mapping) {
933: if (ni > 128) {
934: PetscMalloc1(ni,&lix);
935: }
936: ISLocalToGlobalMappingApply(x->map->mapping,ni,(PetscInt*)ix,lix);
937: (*x->ops->setvalues)(x,ni,lix,y,iora);
938: if (ni > 128) {
939: PetscFree(lix);
940: }
941: } else {
942: (*x->ops->setvalues)(x,ni,ix,y,iora);
943: }
944: } else {
945: (*x->ops->setvalueslocal)(x,ni,ix,y,iora);
946: }
947: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
948: PetscObjectStateIncrease((PetscObject)x);
949: return 0;
950: }
952: /*@
953: VecSetValuesBlockedLocal - Inserts or adds values into certain locations of a vector,
954: using a local ordering of the nodes.
956: Not Collective
958: Input Parameters:
959: + x - vector to insert in
960: . ni - number of blocks to add
961: . ix - indices where to add in block count, not element count
962: . y - array of values
963: - iora - either INSERT_VALUES or ADD_VALUES, where
964: ADD_VALUES adds values to any existing entries, and
965: INSERT_VALUES replaces existing entries with new values
967: Level: intermediate
969: Notes:
970: VecSetValuesBlockedLocal() sets x[bs*ix[i]+j] = y[bs*i+j],
971: for j=0,..bs-1, for i=0,...,ni-1, where bs has been set with VecSetBlockSize().
973: Calls to VecSetValuesBlockedLocal() with the INSERT_VALUES and ADD_VALUES
974: options cannot be mixed without intervening calls to the assembly
975: routines.
977: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
978: MUST be called after all calls to VecSetValuesBlockedLocal() have been completed.
980: VecSetValuesBlockedLocal() uses 0-based indices in Fortran as well as in C.
982: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetValuesBlocked(),
983: VecSetLocalToGlobalMapping()
984: @*/
985: PetscErrorCode VecSetValuesBlockedLocal(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
986: {
987: PetscInt lixp[128],*lix = lixp;
991: if (!ni) return 0;
995: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
996: if (x->map->mapping) {
997: if (ni > 128) {
998: PetscMalloc1(ni,&lix);
999: }
1000: ISLocalToGlobalMappingApplyBlock(x->map->mapping,ni,(PetscInt*)ix,lix);
1001: (*x->ops->setvaluesblocked)(x,ni,lix,y,iora);
1002: if (ni > 128) {
1003: PetscFree(lix);
1004: }
1005: } else {
1006: (*x->ops->setvaluesblocked)(x,ni,ix,y,iora);
1007: }
1008: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
1009: PetscObjectStateIncrease((PetscObject)x);
1010: return 0;
1011: }
1013: /*@
1014: VecMTDot - Computes indefinite vector multiple dot products.
1015: That is, it does NOT use the complex conjugate.
1017: Collective on Vec
1019: Input Parameters:
1020: + x - one vector
1021: . nv - number of vectors
1022: - y - array of vectors. Note that vectors are pointers
1024: Output Parameter:
1025: . val - array of the dot products
1027: Notes for Users of Complex Numbers:
1028: For complex vectors, VecMTDot() computes the indefinite form
1029: $ val = (x,y) = y^T x,
1030: where y^T denotes the transpose of y.
1032: Use VecMDot() for the inner product
1033: $ val = (x,y) = y^H x,
1034: where y^H denotes the conjugate transpose of y.
1036: Level: intermediate
1038: .seealso: VecMDot(), VecTDot()
1039: @*/
1040: PetscErrorCode VecMTDot(Vec x,PetscInt nv,const Vec y[],PetscScalar val[])
1041: {
1044: if (!nv) return 0;
1051: VecCheckSameSize(x,1,*y,3);
1053: PetscLogEventBegin(VEC_MTDot,x,*y,0,0);
1054: (*x->ops->mtdot)(x,nv,y,val);
1055: PetscLogEventEnd(VEC_MTDot,x,*y,0,0);
1056: return 0;
1057: }
1059: /*@
1060: VecMDot - Computes vector multiple dot products.
1062: Collective on Vec
1064: Input Parameters:
1065: + x - one vector
1066: . nv - number of vectors
1067: - y - array of vectors.
1069: Output Parameter:
1070: . val - array of the dot products (does not allocate the array)
1072: Notes for Users of Complex Numbers:
1073: For complex vectors, VecMDot() computes
1074: $ val = (x,y) = y^H x,
1075: where y^H denotes the conjugate transpose of y.
1077: Use VecMTDot() for the indefinite form
1078: $ val = (x,y) = y^T x,
1079: where y^T denotes the transpose of y.
1081: Level: intermediate
1083: .seealso: VecMTDot(), VecDot()
1084: @*/
1085: PetscErrorCode VecMDot(Vec x,PetscInt nv,const Vec y[],PetscScalar val[])
1086: {
1089: if (!nv) return 0;
1097: VecCheckSameSize(x,1,*y,3);
1099: PetscLogEventBegin(VEC_MDot,x,*y,0,0);
1100: (*x->ops->mdot)(x,nv,y,val);
1101: PetscLogEventEnd(VEC_MDot,x,*y,0,0);
1102: return 0;
1103: }
1105: /*@
1106: VecMAXPY - Computes y = y + sum alpha[i] x[i]
1108: Logically Collective on Vec
1110: Input Parameters:
1111: + nv - number of scalars and x-vectors
1112: . alpha - array of scalars
1113: . y - one vector
1114: - x - array of vectors
1116: Level: intermediate
1118: Notes:
1119: y cannot be any of the x vectors
1121: .seealso: VecAYPX(), VecWAXPY(), VecAXPY(), VecAXPBYPCZ(), VecAXPBY()
1122: @*/
1123: PetscErrorCode VecMAXPY(Vec y,PetscInt nv,const PetscScalar alpha[],Vec x[])
1124: {
1125: PetscInt i;
1126: PetscBool nonzero;
1130: if (!nv) return 0;
1138: VecCheckSameSize(y,1,*x,4);
1140: for (i=0, nonzero = PETSC_FALSE; i<nv && !nonzero; i++) nonzero = (PetscBool)(nonzero || alpha[i] != (PetscScalar)0.0);
1141: if (!nonzero) return 0;
1142: VecSetErrorIfLocked(y,1);
1143: PetscLogEventBegin(VEC_MAXPY,*x,y,0,0);
1144: (*y->ops->maxpy)(y,nv,alpha,x);
1145: PetscLogEventEnd(VEC_MAXPY,*x,y,0,0);
1146: PetscObjectStateIncrease((PetscObject)y);
1147: return 0;
1148: }
1150: /*@
1151: VecConcatenate - Creates a new vector that is a vertical concatenation of all the given array of vectors
1152: in the order they appear in the array. The concatenated vector resides on the same
1153: communicator and is the same type as the source vectors.
1155: Collective on X
1157: Input Parameters:
1158: + nx - number of vectors to be concatenated
1159: - X - array containing the vectors to be concatenated in the order of concatenation
1161: Output Parameters:
1162: + Y - concatenated vector
1163: - x_is - array of index sets corresponding to the concatenated components of Y (NULL if not needed)
1165: Notes:
1166: Concatenation is similar to the functionality of a VecNest object; they both represent combination of
1167: different vector spaces. However, concatenated vectors do not store any information about their
1168: sub-vectors and own their own data. Consequently, this function provides index sets to enable the
1169: manipulation of data in the concatenated vector that corresponds to the original components at creation.
1171: This is a useful tool for outer loop algorithms, particularly constrained optimizers, where the solver
1172: has to operate on combined vector spaces and cannot utilize VecNest objects due to incompatibility with
1173: bound projections.
1175: Level: advanced
1177: .seealso: VECNEST, VECSCATTER, VecScatterCreate()
1178: @*/
1179: PetscErrorCode VecConcatenate(PetscInt nx, const Vec X[], Vec *Y, IS *x_is[])
1180: {
1181: MPI_Comm comm;
1182: VecType vec_type;
1183: Vec Ytmp, Xtmp;
1184: IS *is_tmp;
1185: PetscInt i, shift=0, Xnl, Xng, Xbegin;
1192: if ((*X)->ops->concatenate) {
1193: /* use the dedicated concatenation function if available */
1194: (*(*X)->ops->concatenate)(nx,X,Y,x_is);
1195: } else {
1196: /* loop over vectors and start creating IS */
1197: comm = PetscObjectComm((PetscObject)(*X));
1198: VecGetType(*X, &vec_type);
1199: PetscMalloc1(nx, &is_tmp);
1200: for (i=0; i<nx; i++) {
1201: VecGetSize(X[i], &Xng);
1202: VecGetLocalSize(X[i], &Xnl);
1203: VecGetOwnershipRange(X[i], &Xbegin, NULL);
1204: ISCreateStride(comm, Xnl, shift + Xbegin, 1, &is_tmp[i]);
1205: shift += Xng;
1206: }
1207: /* create the concatenated vector */
1208: VecCreate(comm, &Ytmp);
1209: VecSetType(Ytmp, vec_type);
1210: VecSetSizes(Ytmp, PETSC_DECIDE, shift);
1211: VecSetUp(Ytmp);
1212: /* copy data from X array to Y and return */
1213: for (i=0; i<nx; i++) {
1214: VecGetSubVector(Ytmp, is_tmp[i], &Xtmp);
1215: VecCopy(X[i], Xtmp);
1216: VecRestoreSubVector(Ytmp, is_tmp[i], &Xtmp);
1217: }
1218: *Y = Ytmp;
1219: if (x_is) {
1220: *x_is = is_tmp;
1221: } else {
1222: for (i=0; i<nx; i++) {
1223: ISDestroy(&is_tmp[i]);
1224: }
1225: PetscFree(is_tmp);
1226: }
1227: }
1228: return 0;
1229: }
1231: /* A helper function for VecGetSubVector to check if we can implement it with no-copy (i.e. the subvector shares
1232: memory with the original vector), and the block size of the subvector.
1234: Input Parameters:
1235: + X - the original vector
1236: - is - the index set of the subvector
1238: Output Parameters:
1239: + contig - PETSC_TRUE if the index set refers to contiguous entries on this process, else PETSC_FALSE
1240: . start - start of contiguous block, as an offset from the start of the ownership range of the original vector
1241: - blocksize - the block size of the subvector
1243: */
1244: PetscErrorCode VecGetSubVectorContiguityAndBS_Private(Vec X,IS is,PetscBool *contig,PetscInt *start,PetscInt *blocksize)
1245: {
1246: PetscInt gstart,gend,lstart;
1247: PetscBool red[2] = {PETSC_TRUE/*contiguous*/,PETSC_TRUE/*validVBS*/};
1248: PetscInt n,N,ibs,vbs,bs = -1;
1250: ISGetLocalSize(is,&n);
1251: ISGetSize(is,&N);
1252: ISGetBlockSize(is,&ibs);
1253: VecGetBlockSize(X,&vbs);
1254: VecGetOwnershipRange(X,&gstart,&gend);
1255: ISContiguousLocal(is,gstart,gend,&lstart,&red[0]);
1256: /* block size is given by IS if ibs > 1; otherwise, check the vector */
1257: if (ibs > 1) {
1258: MPIU_Allreduce(MPI_IN_PLACE,red,1,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)is));
1259: bs = ibs;
1260: } else {
1261: if (n%vbs || vbs == 1) red[1] = PETSC_FALSE; /* this process invalidate the collectiveness of block size */
1262: MPIU_Allreduce(MPI_IN_PLACE,red,2,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)is));
1263: if (red[0] && red[1]) bs = vbs; /* all processes have a valid block size and the access will be contiguous */
1264: }
1266: *contig = red[0];
1267: *start = lstart;
1268: *blocksize = bs;
1269: return 0;
1270: }
1272: /* A helper function for VecGetSubVector, to be used when we have to build a standalone subvector through VecScatter
1274: Input Parameters:
1275: + X - the original vector
1276: . is - the index set of the subvector
1277: - bs - the block size of the subvector, gotten from VecGetSubVectorContiguityAndBS_Private()
1279: Output Parameters:
1280: . Z - the subvector, which will compose the VecScatter context on output
1281: */
1282: PetscErrorCode VecGetSubVectorThroughVecScatter_Private(Vec X,IS is,PetscInt bs,Vec *Z)
1283: {
1284: PetscInt n,N;
1285: VecScatter vscat;
1286: Vec Y;
1288: ISGetLocalSize(is,&n);
1289: ISGetSize(is,&N);
1290: VecCreate(PetscObjectComm((PetscObject)is),&Y);
1291: VecSetSizes(Y,n,N);
1292: VecSetBlockSize(Y,bs);
1293: VecSetType(Y,((PetscObject)X)->type_name);
1294: VecScatterCreate(X,is,Y,NULL,&vscat);
1295: VecScatterBegin(vscat,X,Y,INSERT_VALUES,SCATTER_FORWARD);
1296: VecScatterEnd(vscat,X,Y,INSERT_VALUES,SCATTER_FORWARD);
1297: PetscObjectCompose((PetscObject)Y,"VecGetSubVector_Scatter",(PetscObject)vscat);
1298: VecScatterDestroy(&vscat);
1299: *Z = Y;
1300: return 0;
1301: }
1303: /*@
1304: VecGetSubVector - Gets a vector representing part of another vector
1306: Collective on X and IS
1308: Input Parameters:
1309: + X - vector from which to extract a subvector
1310: - is - index set representing portion of X to extract
1312: Output Parameter:
1313: . Y - subvector corresponding to is
1315: Level: advanced
1317: Notes:
1318: The subvector Y should be returned with VecRestoreSubVector().
1319: X and is must be defined on the same communicator
1321: This function may return a subvector without making a copy, therefore it is not safe to use the original vector while
1322: modifying the subvector. Other non-overlapping subvectors can still be obtained from X using this function.
1323: 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.
1325: .seealso: MatCreateSubMatrix()
1326: @*/
1327: PetscErrorCode VecGetSubVector(Vec X,IS is,Vec *Y)
1328: {
1329: Vec Z;
1335: if (X->ops->getsubvector) {
1336: (*X->ops->getsubvector)(X,is,&Z);
1337: } else { /* Default implementation currently does no caching */
1338: PetscBool contig;
1339: PetscInt n,N,start,bs;
1341: ISGetLocalSize(is,&n);
1342: ISGetSize(is,&N);
1343: VecGetSubVectorContiguityAndBS_Private(X,is,&contig,&start,&bs);
1344: if (contig) { /* We can do a no-copy implementation */
1345: const PetscScalar *x;
1346: PetscInt state = 0;
1347: PetscBool isstd,iscuda,iship;
1349: PetscObjectTypeCompareAny((PetscObject)X,&isstd,VECSEQ,VECMPI,VECSTANDARD,"");
1350: PetscObjectTypeCompareAny((PetscObject)X,&iscuda,VECSEQCUDA,VECMPICUDA,"");
1351: PetscObjectTypeCompareAny((PetscObject)X,&iship,VECSEQHIP,VECMPIHIP,"");
1352: if (iscuda) {
1353: #if defined(PETSC_HAVE_CUDA)
1354: const PetscScalar *x_d;
1355: PetscMPIInt size;
1356: PetscOffloadMask flg;
1358: VecCUDAGetArrays_Private(X,&x,&x_d,&flg);
1361: if (x) x += start;
1362: if (x_d) x_d += start;
1363: MPI_Comm_size(PetscObjectComm((PetscObject)X),&size);
1364: if (size == 1) {
1365: VecCreateSeqCUDAWithArrays(PetscObjectComm((PetscObject)X),bs,n,x,x_d,&Z);
1366: } else {
1367: VecCreateMPICUDAWithArrays(PetscObjectComm((PetscObject)X),bs,n,N,x,x_d,&Z);
1368: }
1369: Z->offloadmask = flg;
1370: #endif
1371: } else if (iship) {
1372: #if defined(PETSC_HAVE_HIP)
1373: const PetscScalar *x_d;
1374: PetscMPIInt size;
1375: PetscOffloadMask flg;
1377: VecHIPGetArrays_Private(X,&x,&x_d,&flg);
1380: if (x) x += start;
1381: if (x_d) x_d += start;
1382: MPI_Comm_size(PetscObjectComm((PetscObject)X),&size);
1383: if (size == 1) {
1384: VecCreateSeqHIPWithArrays(PetscObjectComm((PetscObject)X),bs,n,x,x_d,&Z);
1385: } else {
1386: VecCreateMPIHIPWithArrays(PetscObjectComm((PetscObject)X),bs,n,N,x,x_d,&Z);
1387: }
1388: Z->offloadmask = flg;
1389: #endif
1390: } else if (isstd) {
1391: PetscMPIInt size;
1393: MPI_Comm_size(PetscObjectComm((PetscObject)X),&size);
1394: VecGetArrayRead(X,&x);
1395: if (x) x += start;
1396: if (size == 1) {
1397: VecCreateSeqWithArray(PetscObjectComm((PetscObject)X),bs,n,x,&Z);
1398: } else {
1399: VecCreateMPIWithArray(PetscObjectComm((PetscObject)X),bs,n,N,x,&Z);
1400: }
1401: VecRestoreArrayRead(X,&x);
1402: } else { /* default implementation: use place array */
1403: VecGetArrayRead(X,&x);
1404: VecCreate(PetscObjectComm((PetscObject)X),&Z);
1405: VecSetType(Z,((PetscObject)X)->type_name);
1406: VecSetSizes(Z,n,N);
1407: VecSetBlockSize(Z,bs);
1408: VecPlaceArray(Z,x ? x+start : NULL);
1409: VecRestoreArrayRead(X,&x);
1410: }
1412: /* this is relevant only in debug mode */
1413: VecLockGet(X,&state);
1414: if (state) {
1415: VecLockReadPush(Z);
1416: }
1417: Z->ops->placearray = NULL;
1418: Z->ops->replacearray = NULL;
1419: } else { /* Have to create a scatter and do a copy */
1420: VecGetSubVectorThroughVecScatter_Private(X,is,bs,&Z);
1421: }
1422: }
1423: /* Record the state when the subvector was gotten so we know whether its values need to be put back */
1424: if (VecGetSubVectorSavedStateId < 0) PetscObjectComposedDataRegister(&VecGetSubVectorSavedStateId);
1425: PetscObjectComposedDataSetInt((PetscObject)Z,VecGetSubVectorSavedStateId,1);
1426: *Y = Z;
1427: return 0;
1428: }
1430: /*@
1431: VecRestoreSubVector - Restores a subvector extracted using VecGetSubVector()
1433: Collective on IS
1435: Input Parameters:
1436: + X - vector from which subvector was obtained
1437: . is - index set representing the subset of X
1438: - Y - subvector being restored
1440: Level: advanced
1442: .seealso: VecGetSubVector()
1443: @*/
1444: PetscErrorCode VecRestoreSubVector(Vec X,IS is,Vec *Y)
1445: {
1446: PETSC_UNUSED PetscObjectState dummystate = 0;
1447: PetscBool unchanged;
1455: if (X->ops->restoresubvector) {
1456: (*X->ops->restoresubvector)(X,is,Y);
1457: } else {
1458: PetscObjectComposedDataGetInt((PetscObject)*Y,VecGetSubVectorSavedStateId,dummystate,unchanged);
1459: if (!unchanged) { /* If Y's state has not changed since VecGetSubVector(), we only need to destroy Y */
1460: VecScatter scatter;
1461: PetscInt state;
1463: VecLockGet(X,&state);
1466: PetscObjectQuery((PetscObject)*Y,"VecGetSubVector_Scatter",(PetscObject*)&scatter);
1467: if (scatter) {
1468: VecScatterBegin(scatter,*Y,X,INSERT_VALUES,SCATTER_REVERSE);
1469: VecScatterEnd(scatter,*Y,X,INSERT_VALUES,SCATTER_REVERSE);
1470: } else {
1471: PetscBool iscuda,iship;
1472: PetscObjectTypeCompareAny((PetscObject)X,&iscuda,VECSEQCUDA,VECMPICUDA,"");
1473: PetscObjectTypeCompareAny((PetscObject)X,&iship,VECSEQHIP,VECMPIHIP,"");
1475: if (iscuda) {
1476: #if defined(PETSC_HAVE_CUDA)
1477: PetscOffloadMask ymask = (*Y)->offloadmask;
1479: /* The offloadmask of X dictates where to move memory
1480: If X GPU data is valid, then move Y data on GPU if needed
1481: Otherwise, move back to the CPU */
1482: switch (X->offloadmask) {
1483: case PETSC_OFFLOAD_BOTH:
1484: if (ymask == PETSC_OFFLOAD_CPU) {
1485: VecCUDAResetArray(*Y);
1486: } else if (ymask == PETSC_OFFLOAD_GPU) {
1487: X->offloadmask = PETSC_OFFLOAD_GPU;
1488: }
1489: break;
1490: case PETSC_OFFLOAD_GPU:
1491: if (ymask == PETSC_OFFLOAD_CPU) {
1492: VecCUDAResetArray(*Y);
1493: }
1494: break;
1495: case PETSC_OFFLOAD_CPU:
1496: if (ymask == PETSC_OFFLOAD_GPU) {
1497: VecResetArray(*Y);
1498: }
1499: break;
1500: case PETSC_OFFLOAD_UNALLOCATED:
1501: case PETSC_OFFLOAD_KOKKOS:
1502: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"This should not happen");
1503: }
1504: #endif
1505: } else if (iship) {
1506: #if defined(PETSC_HAVE_HIP)
1507: PetscOffloadMask ymask = (*Y)->offloadmask;
1509: /* The offloadmask of X dictates where to move memory
1510: If X GPU data is valid, then move Y data on GPU if needed
1511: Otherwise, move back to the CPU */
1512: switch (X->offloadmask) {
1513: case PETSC_OFFLOAD_BOTH:
1514: if (ymask == PETSC_OFFLOAD_CPU) {
1515: VecHIPResetArray(*Y);
1516: } else if (ymask == PETSC_OFFLOAD_GPU) {
1517: X->offloadmask = PETSC_OFFLOAD_GPU;
1518: }
1519: break;
1520: case PETSC_OFFLOAD_GPU:
1521: if (ymask == PETSC_OFFLOAD_CPU) {
1522: VecHIPResetArray(*Y);
1523: }
1524: break;
1525: case PETSC_OFFLOAD_CPU:
1526: if (ymask == PETSC_OFFLOAD_GPU) {
1527: VecResetArray(*Y);
1528: }
1529: break;
1530: case PETSC_OFFLOAD_UNALLOCATED:
1531: case PETSC_OFFLOAD_KOKKOS:
1532: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"This should not happen");
1533: }
1534: #endif
1535: } else {
1536: /* If OpenCL vecs updated the device memory, this triggers a copy on the CPU */
1537: VecResetArray(*Y);
1538: }
1539: PetscObjectStateIncrease((PetscObject)X);
1540: }
1541: }
1542: }
1543: VecDestroy(Y);
1544: return 0;
1545: }
1547: /*@
1548: VecGetLocalVectorRead - Maps the local portion of a vector into a
1549: vector. You must call VecRestoreLocalVectorRead() when the local
1550: vector is no longer needed.
1552: Not collective.
1554: Input parameter:
1555: . v - The vector for which the local vector is desired.
1557: Output parameter:
1558: . w - Upon exit this contains the local vector.
1560: Level: beginner
1562: Notes:
1563: This function is similar to VecGetArrayRead() which maps the local
1564: portion into a raw pointer. VecGetLocalVectorRead() is usually
1565: almost as efficient as VecGetArrayRead() but in certain circumstances
1566: VecGetLocalVectorRead() can be much more efficient than
1567: VecGetArrayRead(). This is because the construction of a contiguous
1568: array representing the vector data required by VecGetArrayRead() can
1569: be an expensive operation for certain vector types. For example, for
1570: GPU vectors VecGetArrayRead() requires that the data between device
1571: and host is synchronized.
1573: Unlike VecGetLocalVector(), this routine is not collective and
1574: preserves cached information.
1576: .seealso: VecRestoreLocalVectorRead(), VecGetLocalVector(), VecGetArrayRead(), VecGetArray()
1577: @*/
1578: PetscErrorCode VecGetLocalVectorRead(Vec v,Vec w)
1579: {
1580: PetscScalar *a;
1584: VecCheckSameLocalSize(v,1,w,2);
1585: if (v->ops->getlocalvectorread) {
1586: (*v->ops->getlocalvectorread)(v,w);
1587: } else {
1588: VecGetArrayRead(v,(const PetscScalar**)&a);
1589: VecPlaceArray(w,a);
1590: }
1591: PetscObjectStateIncrease((PetscObject)w);
1592: VecLockReadPush(v);
1593: VecLockReadPush(w);
1594: return 0;
1595: }
1597: /*@
1598: VecRestoreLocalVectorRead - Unmaps the local portion of a vector
1599: previously mapped into a vector using VecGetLocalVectorRead().
1601: Not collective.
1603: Input parameter:
1604: + v - The local portion of this vector was previously mapped into w using VecGetLocalVectorRead().
1605: - w - The vector into which the local portion of v was mapped.
1607: Level: beginner
1609: .seealso: VecGetLocalVectorRead(), VecGetLocalVector(), VecGetArrayRead(), VecGetArray()
1610: @*/
1611: PetscErrorCode VecRestoreLocalVectorRead(Vec v,Vec w)
1612: {
1613: PetscScalar *a;
1617: if (v->ops->restorelocalvectorread) {
1618: (*v->ops->restorelocalvectorread)(v,w);
1619: } else {
1620: VecGetArrayRead(w,(const PetscScalar**)&a);
1621: VecRestoreArrayRead(v,(const PetscScalar**)&a);
1622: VecResetArray(w);
1623: }
1624: VecLockReadPop(v);
1625: VecLockReadPop(w);
1626: PetscObjectStateIncrease((PetscObject)w);
1627: return 0;
1628: }
1630: /*@
1631: VecGetLocalVector - Maps the local portion of a vector into a
1632: vector.
1634: Collective on v, not collective on w.
1636: Input parameter:
1637: . v - The vector for which the local vector is desired.
1639: Output parameter:
1640: . w - Upon exit this contains the local vector.
1642: Level: beginner
1644: Notes:
1645: This function is similar to VecGetArray() which maps the local
1646: portion into a raw pointer. VecGetLocalVector() is usually about as
1647: efficient as VecGetArray() but in certain circumstances
1648: VecGetLocalVector() can be much more efficient than VecGetArray().
1649: This is because the construction of a contiguous array representing
1650: the vector data required by VecGetArray() can be an expensive
1651: operation for certain vector types. For example, for GPU vectors
1652: VecGetArray() requires that the data between device and host is
1653: synchronized.
1655: .seealso: VecRestoreLocalVector(), VecGetLocalVectorRead(), VecGetArrayRead(), VecGetArray()
1656: @*/
1657: PetscErrorCode VecGetLocalVector(Vec v,Vec w)
1658: {
1659: PetscScalar *a;
1663: VecCheckSameLocalSize(v,1,w,2);
1664: if (v->ops->getlocalvector) {
1665: (*v->ops->getlocalvector)(v,w);
1666: } else {
1667: VecGetArray(v,&a);
1668: VecPlaceArray(w,a);
1669: }
1670: PetscObjectStateIncrease((PetscObject)w);
1671: return 0;
1672: }
1674: /*@
1675: VecRestoreLocalVector - Unmaps the local portion of a vector
1676: previously mapped into a vector using VecGetLocalVector().
1678: Logically collective.
1680: Input parameter:
1681: + v - The local portion of this vector was previously mapped into w using VecGetLocalVector().
1682: - w - The vector into which the local portion of v was mapped.
1684: Level: beginner
1686: .seealso: VecGetLocalVector(), VecGetLocalVectorRead(), VecRestoreLocalVectorRead(), LocalVectorRead(), VecGetArrayRead(), VecGetArray()
1687: @*/
1688: PetscErrorCode VecRestoreLocalVector(Vec v,Vec w)
1689: {
1690: PetscScalar *a;
1694: if (v->ops->restorelocalvector) {
1695: (*v->ops->restorelocalvector)(v,w);
1696: } else {
1697: VecGetArray(w,&a);
1698: VecRestoreArray(v,&a);
1699: VecResetArray(w);
1700: }
1701: PetscObjectStateIncrease((PetscObject)w);
1702: PetscObjectStateIncrease((PetscObject)v);
1703: return 0;
1704: }
1706: /*@C
1707: VecGetArray - Returns a pointer to a contiguous array that contains this
1708: processor's portion of the vector data. For the standard PETSc
1709: vectors, VecGetArray() returns a pointer to the local data array and
1710: does not use any copies. If the underlying vector data is not stored
1711: in a contiguous array this routine will copy the data to a contiguous
1712: array and return a pointer to that. You MUST call VecRestoreArray()
1713: when you no longer need access to the array.
1715: Logically Collective on Vec
1717: Input Parameter:
1718: . x - the vector
1720: Output Parameter:
1721: . a - location to put pointer to the array
1723: Fortran Note:
1724: This routine is used differently from Fortran 77
1725: $ Vec x
1726: $ PetscScalar x_array(1)
1727: $ PetscOffset i_x
1728: $ PetscErrorCode ierr
1729: $ call VecGetArray(x,x_array,i_x,ierr)
1730: $
1731: $ Access first local entry in vector with
1732: $ value = x_array(i_x + 1)
1733: $
1734: $ ...... other code
1735: $ call VecRestoreArray(x,x_array,i_x,ierr)
1736: For Fortran 90 see VecGetArrayF90()
1738: See the Fortran chapter of the users manual and
1739: petsc/src/snes/tutorials/ex5f.F for details.
1741: Level: beginner
1743: .seealso: VecRestoreArray(), VecGetArrayRead(), VecGetArrays(), VecGetArrayF90(), VecGetArrayReadF90(), VecPlaceArray(), VecGetArray2d(),
1744: VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayWrite(), VecRestoreArrayWrite()
1745: @*/
1746: PetscErrorCode VecGetArray(Vec x,PetscScalar **a)
1747: {
1749: VecSetErrorIfLocked(x,1);
1750: if (x->ops->getarray) { /* The if-else order matters! VECNEST, VECCUDA etc should have ops->getarray while VECCUDA etc are petscnative */
1751: (*x->ops->getarray)(x,a);
1752: } else if (x->petscnative) { /* VECSTANDARD */
1753: *a = *((PetscScalar**)x->data);
1754: } else SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Cannot get array for vector type \"%s\"",((PetscObject)x)->type_name);
1755: return 0;
1756: }
1758: /*@C
1759: VecRestoreArray - Restores a vector after VecGetArray() has been called.
1761: Logically Collective on Vec
1763: Input Parameters:
1764: + x - the vector
1765: - a - location of pointer to array obtained from VecGetArray()
1767: Level: beginner
1769: .seealso: VecGetArray(), VecRestoreArrayRead(), VecRestoreArrays(), VecRestoreArrayF90(), VecRestoreArrayReadF90(), VecPlaceArray(), VecRestoreArray2d(),
1770: VecGetArrayPair(), VecRestoreArrayPair()
1771: @*/
1772: PetscErrorCode VecRestoreArray(Vec x,PetscScalar **a)
1773: {
1775: if (x->ops->restorearray) { /* VECNEST, VECCUDA etc */
1776: (*x->ops->restorearray)(x,a);
1777: } else if (x->petscnative) { /* VECSTANDARD */
1778: /* nothing */
1779: } else SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Cannot restore array for vector type \"%s\"",((PetscObject)x)->type_name);
1780: if (a) *a = NULL;
1781: PetscObjectStateIncrease((PetscObject)x);
1782: return 0;
1783: }
1784: /*@C
1785: VecGetArrayRead - Get read-only pointer to contiguous array containing this processor's portion of the vector data.
1787: Not Collective
1789: Input Parameter:
1790: . x - the vector
1792: Output Parameter:
1793: . a - the array
1795: Level: beginner
1797: Notes:
1798: The array must be returned using a matching call to VecRestoreArrayRead().
1800: Unlike VecGetArray(), this routine is not collective and preserves cached information like vector norms.
1802: Standard PETSc vectors use contiguous storage so that this routine does not perform a copy. Other vector
1803: implementations may require a copy, but must such implementations should cache the contiguous representation so that
1804: only one copy is performed when this routine is called multiple times in sequence.
1806: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair()
1807: @*/
1808: PetscErrorCode VecGetArrayRead(Vec x,const PetscScalar **a)
1809: {
1811: if (x->ops->getarray) { /* VECNEST, VECCUDA, VECKOKKOS etc */
1812: (*x->ops->getarray)(x,(PetscScalar**)a);
1813: } else if (x->petscnative) { /* VECSTANDARD */
1814: *a = *((PetscScalar**)x->data);
1815: } else SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Cannot get array read for vector type \"%s\"",((PetscObject)x)->type_name);
1816: return 0;
1817: }
1819: /*@C
1820: VecRestoreArrayRead - Restore array obtained with VecGetArrayRead()
1822: Not Collective
1824: Input Parameters:
1825: + vec - the vector
1826: - array - the array
1828: Level: beginner
1830: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair()
1831: @*/
1832: PetscErrorCode VecRestoreArrayRead(Vec x,const PetscScalar **a)
1833: {
1835: if (x->petscnative) { /* VECSTANDARD, VECCUDA, VECKOKKOS etc */
1836: /* nothing */
1837: } else if (x->ops->restorearrayread) { /* VECNEST */
1838: (*x->ops->restorearrayread)(x,a);
1839: } else { /* No one? */
1840: (*x->ops->restorearray)(x,(PetscScalar**)a);
1841: }
1842: if (a) *a = NULL;
1843: return 0;
1844: }
1846: /*@C
1847: VecGetArrayWrite - Returns a pointer to a contiguous array that WILL contains this
1848: processor's portion of the vector data. The values in this array are NOT valid, the routine calling this
1849: routine is responsible for putting values into the array; any values it does not set will be invalid
1851: Logically Collective on Vec
1853: Input Parameter:
1854: . x - the vector
1856: Output Parameter:
1857: . a - location to put pointer to the array
1859: Level: intermediate
1861: This is for vectors associate with GPUs, the vector is not copied up before giving access. If you need correct
1862: values in the array use VecGetArray()
1864: .seealso: VecRestoreArray(), VecGetArrayRead(), VecGetArrays(), VecGetArrayF90(), VecGetArrayReadF90(), VecPlaceArray(), VecGetArray2d(),
1865: VecGetArrayPair(), VecRestoreArrayPair(), VecGetArray(), VecRestoreArrayWrite()
1866: @*/
1867: PetscErrorCode VecGetArrayWrite(Vec x,PetscScalar **a)
1868: {
1870: VecSetErrorIfLocked(x,1);
1871: if (x->ops->getarraywrite) {
1872: (*x->ops->getarraywrite)(x,a);
1873: } else {
1874: VecGetArray(x,a);
1875: }
1876: return 0;
1877: }
1879: /*@C
1880: VecRestoreArrayWrite - Restores a vector after VecGetArrayWrite() has been called.
1882: Logically Collective on Vec
1884: Input Parameters:
1885: + x - the vector
1886: - a - location of pointer to array obtained from VecGetArray()
1888: Level: beginner
1890: .seealso: VecGetArray(), VecRestoreArrayRead(), VecRestoreArrays(), VecRestoreArrayF90(), VecRestoreArrayReadF90(), VecPlaceArray(), VecRestoreArray2d(),
1891: VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayWrite()
1892: @*/
1893: PetscErrorCode VecRestoreArrayWrite(Vec x,PetscScalar **a)
1894: {
1896: if (x->ops->restorearraywrite) {
1897: (*x->ops->restorearraywrite)(x,a);
1898: } else if (x->ops->restorearray) {
1899: (*x->ops->restorearray)(x,a);
1900: }
1901: if (a) *a = NULL;
1902: PetscObjectStateIncrease((PetscObject)x);
1903: return 0;
1904: }
1906: /*@C
1907: VecGetArrays - Returns a pointer to the arrays in a set of vectors
1908: that were created by a call to VecDuplicateVecs(). You MUST call
1909: VecRestoreArrays() when you no longer need access to the array.
1911: Logically Collective on Vec
1913: Input Parameters:
1914: + x - the vectors
1915: - n - the number of vectors
1917: Output Parameter:
1918: . a - location to put pointer to the array
1920: Fortran Note:
1921: This routine is not supported in Fortran.
1923: Level: intermediate
1925: .seealso: VecGetArray(), VecRestoreArrays()
1926: @*/
1927: PetscErrorCode VecGetArrays(const Vec x[],PetscInt n,PetscScalar **a[])
1928: {
1929: PetscInt i;
1930: PetscScalar **q;
1936: PetscMalloc1(n,&q);
1937: for (i=0; i<n; ++i) {
1938: VecGetArray(x[i],&q[i]);
1939: }
1940: *a = q;
1941: return 0;
1942: }
1944: /*@C
1945: VecRestoreArrays - Restores a group of vectors after VecGetArrays()
1946: has been called.
1948: Logically Collective on Vec
1950: Input Parameters:
1951: + x - the vector
1952: . n - the number of vectors
1953: - a - location of pointer to arrays obtained from VecGetArrays()
1955: Notes:
1956: For regular PETSc vectors this routine does not involve any copies. For
1957: any special vectors that do not store local vector data in a contiguous
1958: array, this routine will copy the data back into the underlying
1959: vector data structure from the arrays obtained with VecGetArrays().
1961: Fortran Note:
1962: This routine is not supported in Fortran.
1964: Level: intermediate
1966: .seealso: VecGetArrays(), VecRestoreArray()
1967: @*/
1968: PetscErrorCode VecRestoreArrays(const Vec x[],PetscInt n,PetscScalar **a[])
1969: {
1970: PetscInt i;
1971: PetscScalar **q = *a;
1977: for (i=0; i<n; ++i) {
1978: VecRestoreArray(x[i],&q[i]);
1979: }
1980: PetscFree(q);
1981: return 0;
1982: }
1984: /*@C
1985: VecGetArrayAndMemType - Like VecGetArray(), but if this is a standard device vector (e.g., VECCUDA), the returned pointer will be a device
1986: pointer to the device memory that contains this processor's portion of the vector data. Device data is guaranteed to have the latest value.
1987: Otherwise, when this is a host vector (e.g., VECMPI), this routine functions the same as VecGetArray() and returns a host pointer.
1989: For VECKOKKOS, if Kokkos is configured without device (e.g., use serial or openmp), per this function, the vector works like VECSEQ/VECMPI;
1990: otherwise, it works like VECCUDA or VECHIP etc.
1992: Logically Collective on Vec
1994: Input Parameter:
1995: . x - the vector
1997: Output Parameters:
1998: + a - location to put pointer to the array
1999: - mtype - memory type of the array
2001: Level: beginner
2003: .seealso: VecRestoreArrayAndMemType(), VecGetArrayReadAndMemType(), VecGetArrayWriteAndMemType(), VecRestoreArray(), VecGetArrayRead(), VecGetArrays(), VecGetArrayF90(), VecGetArrayReadF90(),
2004: VecPlaceArray(), VecGetArray2d(), VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayWrite(), VecRestoreArrayWrite()
2005: @*/
2006: PetscErrorCode VecGetArrayAndMemType(Vec x,PetscScalar **a,PetscMemType *mtype)
2007: {
2008: PetscMemType omtype;
2012: VecSetErrorIfLocked(x,1);
2013: if (x->ops->getarrayandmemtype) { /* VECCUDA, VECKOKKOS etc */
2014: (*x->ops->getarrayandmemtype)(x,a,&omtype);
2015: } else { /* VECSTANDARD, VECNEST, VECVIENNACL */
2016: VecGetArray(x,a);
2017: omtype = PETSC_MEMTYPE_HOST;
2018: }
2019: if (mtype) *mtype = omtype;
2020: return 0;
2021: }
2023: /*@C
2024: VecRestoreArrayAndMemType - Restores a vector after VecGetArrayAndMemType() has been called.
2026: Logically Collective on Vec
2028: Input Parameters:
2029: + x - the vector
2030: - a - location of pointer to array obtained from VecGetArrayAndMemType()
2032: Level: beginner
2034: .seealso: VecGetArrayAndMemType(), VecGetArray(), VecRestoreArrayRead(), VecRestoreArrays(), VecRestoreArrayF90(), VecRestoreArrayReadF90(),
2035: VecPlaceArray(), VecRestoreArray2d(), VecGetArrayPair(), VecRestoreArrayPair()
2036: @*/
2037: PetscErrorCode VecRestoreArrayAndMemType(Vec x,PetscScalar **a)
2038: {
2041: if (x->ops->restorearrayandmemtype) { /* VECCUDA, VECKOKKOS etc */
2042: (*x->ops->restorearrayandmemtype)(x,a);
2043: } else if (x->ops->restorearray) { /* VECNEST, VECVIENNACL */
2044: (*x->ops->restorearray)(x,a);
2045: } /* VECSTANDARD does nothing */
2046: if (a) *a = NULL;
2047: PetscObjectStateIncrease((PetscObject)x);
2048: return 0;
2049: }
2051: /*@C
2052: VecGetArrayReadAndMemType - Like VecGetArrayRead(), but if the input vector is a device vector, it will return a read-only device pointer. The returned pointer is guarenteed to point to up-to-date data. For host vectors, it functions as VecGetArrayRead().
2054: Not Collective
2056: Input Parameter:
2057: . x - the vector
2059: Output Parameters:
2060: + a - the array
2061: - mtype - memory type of the array
2063: Level: beginner
2065: Notes:
2066: The array must be returned using a matching call to VecRestoreArrayReadAndMemType().
2068: .seealso: VecRestoreArrayReadAndMemType(), VecGetArrayAndMemType(), VecGetArrayWriteAndMemType(), VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayAndMemType()
2069: @*/
2070: PetscErrorCode VecGetArrayReadAndMemType(Vec x,const PetscScalar **a,PetscMemType *mtype)
2071: {
2072: PetscMemType omtype;
2076: #if defined(PETSC_USE_DEBUG)
2078: #endif
2080: if (x->ops->getarrayandmemtype) { /* VECCUDA, VECKOKKOS etc, though they are also petscnative */
2081: (*x->ops->getarrayandmemtype)(x,(PetscScalar**)a,&omtype);
2082: } else if (x->ops->getarray) { /* VECNEST, VECVIENNACL */
2083: (*x->ops->getarray)(x,(PetscScalar**)a);
2084: omtype = PETSC_MEMTYPE_HOST;
2085: } else if (x->petscnative) { /* VECSTANDARD */
2086: *a = *((PetscScalar**)x->data);
2087: omtype = PETSC_MEMTYPE_HOST;
2088: } else SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Cannot get array read in place for vector type \"%s\"",((PetscObject)x)->type_name);
2089: if (mtype) *mtype = omtype;
2090: return 0;
2091: }
2093: /*@C
2094: VecRestoreArrayReadAndMemType - Restore array obtained with VecGetArrayReadAndMemType()
2096: Not Collective
2098: Input Parameters:
2099: + vec - the vector
2100: - array - the array
2102: Level: beginner
2104: .seealso: VecGetArrayReadAndMemType(), VecRestoreArrayAndMemType(), VecRestoreArrayWriteAndMemType(), VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair()
2105: @*/
2106: PetscErrorCode VecRestoreArrayReadAndMemType(Vec x,const PetscScalar **a)
2107: {
2110: if (x->petscnative) { /* VECSTANDARD, VECCUDA, VECKOKKOS, VECVIENNACL etc */
2111: /* nothing */
2112: } else if (x->ops->restorearrayread) { /* VECNEST */
2113: (*x->ops->restorearrayread)(x,a);
2114: } else SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Cannot restore array read in place for vector type \"%s\"",((PetscObject)x)->type_name);
2115: if (a) *a = NULL;
2116: return 0;
2117: }
2119: /*@C
2120: VecGetArrayWriteAndMemType - Like VecGetArrayWrite(), but if this is a device vector it will aways return
2121: a device pointer to the device memory that contains this processor's portion of the vector data.
2123: Not Collective
2125: Input Parameter:
2126: . x - the vector
2128: Output Parameters:
2129: + a - the array
2130: - mtype - memory type of the array
2132: Level: beginner
2134: Notes:
2135: The array must be returned using a matching call to VecRestoreArrayWriteAndMemType(), where it will label the device memory as most recent.
2137: .seealso: VecRestoreArrayWriteAndMemType(), VecGetArrayReadAndMemType(), VecGetArrayAndMemType(), VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair(),
2138: @*/
2139: PetscErrorCode VecGetArrayWriteAndMemType(Vec x,PetscScalar **a,PetscMemType *mtype)
2140: {
2141: PetscMemType omtype;
2145: if (x->ops->getarraywriteandmemtype) { /* VECCUDA, VECHIP, VECKOKKOS etc, though they are also petscnative */
2146: (*x->ops->getarrayandmemtype)(x,a,&omtype);
2147: } else if (x->ops->getarraywrite) { /* VECNEST, VECVIENNACL */
2148: (*x->ops->getarraywrite)(x,a);
2149: omtype = PETSC_MEMTYPE_HOST;
2150: } else if (x->petscnative) { /* VECSTANDARD */
2151: *a = *((PetscScalar**)x->data);
2152: omtype = PETSC_MEMTYPE_HOST;
2153: } else SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Cannot get array read in place for vector type \"%s\"",((PetscObject)x)->type_name);
2154: if (mtype) *mtype = omtype;
2155: return 0;
2156: }
2158: /*@C
2159: VecRestoreArrayWriteAndMemType - Restore array obtained with VecGetArrayWriteAndMemType()
2161: Not Collective
2163: Input Parameters:
2164: + vec - the vector
2165: - array - the array
2167: Level: beginner
2169: .seealso: VecGetArrayWriteAndMemType(), VecRestoreArrayAndMemType(), VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair()
2170: @*/
2171: PetscErrorCode VecRestoreArrayWriteAndMemType(Vec x,PetscScalar **a)
2172: {
2173: VecRestoreArrayAndMemType(x,a);
2174: return 0;
2175: }
2177: /*@
2178: VecPlaceArray - Allows one to replace the array in a vector with an
2179: array provided by the user. This is useful to avoid copying an array
2180: into a vector.
2182: Not Collective
2184: Input Parameters:
2185: + vec - the vector
2186: - array - the array
2188: Notes:
2189: You can return to the original array with a call to `VecResetArray()`. `vec` does not take
2190: ownership of `array` in any way. The user must free `array` themselves but be careful not to
2191: do so before the vector has either been destroyed, had its original array restored with
2192: `VecResetArray()` or permanently replaced with `VecReplaceArray()`.
2194: Level: developer
2196: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecResetArray()
2198: @*/
2199: PetscErrorCode VecPlaceArray(Vec vec,const PetscScalar array[])
2200: {
2204: if (vec->ops->placearray) {
2205: (*vec->ops->placearray)(vec,array);
2206: } else SETERRQ(PetscObjectComm((PetscObject)vec),PETSC_ERR_SUP,"Cannot place array in this type of vector");
2207: PetscObjectStateIncrease((PetscObject)vec);
2208: return 0;
2209: }
2211: /*@C
2212: VecReplaceArray - Allows one to replace the array in a vector with an
2213: array provided by the user. This is useful to avoid copying an array
2214: into a vector.
2216: Not Collective
2218: Input Parameters:
2219: + vec - the vector
2220: - array - the array
2222: Notes:
2223: This permanently replaces the array and frees the memory associated
2224: with the old array.
2226: The memory passed in MUST be obtained with PetscMalloc() and CANNOT be
2227: freed by the user. It will be freed when the vector is destroyed.
2229: Not supported from Fortran
2231: Level: developer
2233: .seealso: VecGetArray(), VecRestoreArray(), VecPlaceArray(), VecResetArray()
2235: @*/
2236: PetscErrorCode VecReplaceArray(Vec vec,const PetscScalar array[])
2237: {
2240: if (vec->ops->replacearray) {
2241: (*vec->ops->replacearray)(vec,array);
2242: } else SETERRQ(PetscObjectComm((PetscObject)vec),PETSC_ERR_SUP,"Cannot replace array in this type of vector");
2243: PetscObjectStateIncrease((PetscObject)vec);
2244: return 0;
2245: }
2247: /*@C
2248: VecCUDAGetArray - Provides access to the CUDA buffer inside a vector.
2250: This function has semantics similar to VecGetArray(): the pointer
2251: returned by this function points to a consistent view of the vector
2252: data. This may involve a copy operation of data from the host to the
2253: device if the data on the device is out of date. If the device
2254: memory hasn't been allocated previously it will be allocated as part
2255: of this function call. VecCUDAGetArray() assumes that
2256: the user will modify the vector data. This is similar to
2257: intent(inout) in fortran.
2259: The CUDA device pointer has to be released by calling
2260: VecCUDARestoreArray(). Upon restoring the vector data
2261: the data on the host will be marked as out of date. A subsequent
2262: access of the host data will thus incur a data transfer from the
2263: device to the host.
2265: Input Parameter:
2266: . v - the vector
2268: Output Parameter:
2269: . a - the CUDA device pointer
2271: Fortran note:
2272: This function is not currently available from Fortran.
2274: Level: intermediate
2276: .seealso: VecCUDARestoreArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2277: @*/
2278: PETSC_EXTERN PetscErrorCode VecCUDAGetArray(Vec v, PetscScalar **a)
2279: {
2281: #if defined(PETSC_HAVE_CUDA)
2282: {
2283: VecCUDACopyToGPU(v);
2284: *a = ((Vec_CUDA*)v->spptr)->GPUarray;
2285: }
2286: #endif
2287: return 0;
2288: }
2290: /*@C
2291: VecCUDARestoreArray - Restore a CUDA device pointer previously acquired with VecCUDAGetArray().
2293: This marks the host data as out of date. Subsequent access to the
2294: vector data on the host side with for instance VecGetArray() incurs a
2295: data transfer.
2297: Input Parameters:
2298: + v - the vector
2299: - a - the CUDA device pointer. This pointer is invalid after
2300: VecCUDARestoreArray() returns.
2302: Fortran note:
2303: This function is not currently available from Fortran.
2305: Level: intermediate
2307: .seealso: VecCUDAGetArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2308: @*/
2309: PETSC_EXTERN PetscErrorCode VecCUDARestoreArray(Vec v, PetscScalar **a)
2310: {
2312: #if defined(PETSC_HAVE_CUDA)
2313: v->offloadmask = PETSC_OFFLOAD_GPU;
2314: #endif
2315: PetscObjectStateIncrease((PetscObject)v);
2316: return 0;
2317: }
2319: /*@C
2320: VecCUDAGetArrayRead - Provides read access to the CUDA buffer inside a vector.
2322: This function is analogous to VecGetArrayRead(): The pointer
2323: returned by this function points to a consistent view of the vector
2324: data. This may involve a copy operation of data from the host to the
2325: device if the data on the device is out of date. If the device
2326: memory hasn't been allocated previously it will be allocated as part
2327: of this function call. VecCUDAGetArrayRead() assumes that the
2328: user will not modify the vector data. This is analgogous to
2329: intent(in) in Fortran.
2331: The CUDA device pointer has to be released by calling
2332: VecCUDARestoreArrayRead(). If the data on the host side was
2333: previously up to date it will remain so, i.e. data on both the device
2334: and the host is up to date. Accessing data on the host side does not
2335: incur a device to host data transfer.
2337: Input Parameter:
2338: . v - the vector
2340: Output Parameter:
2341: . a - the CUDA pointer.
2343: Fortran note:
2344: This function is not currently available from Fortran.
2346: Level: intermediate
2348: .seealso: VecCUDARestoreArrayRead(), VecCUDAGetArray(), VecCUDAGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2349: @*/
2350: PETSC_EXTERN PetscErrorCode VecCUDAGetArrayRead(Vec v,const PetscScalar** a)
2351: {
2352: VecCUDAGetArray(v,(PetscScalar**)a);
2353: return 0;
2354: }
2356: /*@C
2357: VecCUDARestoreArrayRead - Restore a CUDA device pointer previously acquired with VecCUDAGetArrayRead().
2359: If the data on the host side was previously up to date it will remain
2360: so, i.e. data on both the device and the host is up to date.
2361: Accessing data on the host side e.g. with VecGetArray() does not
2362: incur a device to host data transfer.
2364: Input Parameters:
2365: + v - the vector
2366: - a - the CUDA device pointer. This pointer is invalid after
2367: VecCUDARestoreArrayRead() returns.
2369: Fortran note:
2370: This function is not currently available from Fortran.
2372: Level: intermediate
2374: .seealso: VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecCUDAGetArray(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2375: @*/
2376: PETSC_EXTERN PetscErrorCode VecCUDARestoreArrayRead(Vec v, const PetscScalar **a)
2377: {
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: {
2414: #if defined(PETSC_HAVE_CUDA)
2415: {
2416: VecCUDAAllocateCheck(v);
2417: *a = ((Vec_CUDA*)v->spptr)->GPUarray;
2418: }
2419: #endif
2420: return 0;
2421: }
2423: /*@C
2424: VecCUDARestoreArrayWrite - Restore a CUDA device pointer previously acquired with VecCUDAGetArrayWrite().
2426: Data on the host will be marked as out of date. Subsequent access of
2427: the data on the host side e.g. with VecGetArray() will incur a device
2428: to host data transfer.
2430: Input Parameters:
2431: + v - the vector
2432: - a - the CUDA device pointer. This pointer is invalid after
2433: VecCUDARestoreArrayWrite() returns.
2435: Fortran note:
2436: This function is not currently available from Fortran.
2438: Level: intermediate
2440: .seealso: VecCUDAGetArrayWrite(), VecCUDAGetArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2441: @*/
2442: PETSC_EXTERN PetscErrorCode VecCUDARestoreArrayWrite(Vec v, PetscScalar **a)
2443: {
2445: #if defined(PETSC_HAVE_CUDA)
2446: v->offloadmask = PETSC_OFFLOAD_GPU;
2447: if (a) *a = NULL;
2448: #endif
2449: PetscObjectStateIncrease((PetscObject)v);
2450: return 0;
2451: }
2453: /*@C
2454: VecCUDAPlaceArray - Allows one to replace the GPU array in a vector with a
2455: GPU array provided by the user. This is useful to avoid copying an
2456: array into a vector.
2458: Not Collective
2460: Input Parameters:
2461: + vec - the vector
2462: - array - the GPU array
2464: Notes:
2465: You can return to the original GPU array with a call to VecCUDAResetArray()
2466: It is not possible to use VecCUDAPlaceArray() and VecPlaceArray() at the
2467: same time on the same vector.
2469: `vec` does not take ownership of `array` in any way. The user must free `array` themselves
2470: but be careful not to do so before the vector has either been destroyed, had its original
2471: array restored with `VecCUDAResetArray()` or permanently replaced with
2472: `VecCUDAReplaceArray()`.
2474: Level: developer
2476: .seealso: VecPlaceArray(), VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecResetArray(), VecCUDAResetArray(), VecCUDAReplaceArray()
2478: @*/
2479: PetscErrorCode VecCUDAPlaceArray(Vec vin,const PetscScalar a[])
2480: {
2482: #if defined(PETSC_HAVE_CUDA)
2483: VecCUDACopyToGPU(vin);
2485: ((Vec_Seq*)vin->data)->unplacedarray = (PetscScalar *) ((Vec_CUDA*)vin->spptr)->GPUarray; /* save previous GPU array so reset can bring it back */
2486: ((Vec_CUDA*)vin->spptr)->GPUarray = (PetscScalar*)a;
2487: vin->offloadmask = PETSC_OFFLOAD_GPU;
2488: #endif
2489: PetscObjectStateIncrease((PetscObject)vin);
2490: return 0;
2491: }
2493: /*@C
2494: VecCUDAReplaceArray - Allows one to replace the GPU array in a vector
2495: with a GPU array provided by the user. This is useful to avoid copying
2496: a GPU array into a vector.
2498: Not Collective
2500: Input Parameters:
2501: + vec - the vector
2502: - array - the GPU array
2504: Notes:
2505: This permanently replaces the GPU array and frees the memory associated
2506: with the old GPU array.
2508: The memory passed in CANNOT be freed by the user. It will be freed
2509: when the vector is destroyed.
2511: Not supported from Fortran
2513: Level: developer
2515: .seealso: VecGetArray(), VecRestoreArray(), VecPlaceArray(), VecResetArray(), VecCUDAResetArray(), VecCUDAPlaceArray(), VecReplaceArray()
2517: @*/
2518: PetscErrorCode VecCUDAReplaceArray(Vec vin,const PetscScalar a[])
2519: {
2520: #if defined(PETSC_HAVE_CUDA)
2521: #endif
2524: #if defined(PETSC_HAVE_CUDA)
2525: if (((Vec_CUDA*)vin->spptr)->GPUarray_allocated) {
2526: cudaFree(((Vec_CUDA*)vin->spptr)->GPUarray_allocated);
2527: }
2528: ((Vec_CUDA*)vin->spptr)->GPUarray_allocated = ((Vec_CUDA*)vin->spptr)->GPUarray = (PetscScalar*)a;
2529: vin->offloadmask = PETSC_OFFLOAD_GPU;
2530: #endif
2531: PetscObjectStateIncrease((PetscObject)vin);
2532: return 0;
2533: }
2535: /*@C
2536: VecCUDAResetArray - Resets a vector to use its default memory. Call this
2537: after the use of VecCUDAPlaceArray().
2539: Not Collective
2541: Input Parameters:
2542: . vec - the vector
2544: Level: developer
2546: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecPlaceArray(), VecResetArray(), VecCUDAPlaceArray(), VecCUDAReplaceArray()
2548: @*/
2549: PetscErrorCode VecCUDAResetArray(Vec vin)
2550: {
2552: #if defined(PETSC_HAVE_CUDA)
2553: VecCUDACopyToGPU(vin);
2554: ((Vec_CUDA*)vin->spptr)->GPUarray = (PetscScalar *) ((Vec_Seq*)vin->data)->unplacedarray;
2555: ((Vec_Seq*)vin->data)->unplacedarray = 0;
2556: vin->offloadmask = PETSC_OFFLOAD_GPU;
2557: #endif
2558: PetscObjectStateIncrease((PetscObject)vin);
2559: return 0;
2560: }
2562: /*@C
2563: VecHIPGetArray - Provides access to the HIP buffer inside a vector.
2565: This function has semantics similar to VecGetArray(): the pointer
2566: returned by this function points to a consistent view of the vector
2567: data. This may involve a copy operation of data from the host to the
2568: device if the data on the device is out of date. If the device
2569: memory hasn't been allocated previously it will be allocated as part
2570: of this function call. VecHIPGetArray() assumes that
2571: the user will modify the vector data. This is similar to
2572: intent(inout) in fortran.
2574: The HIP device pointer has to be released by calling
2575: VecHIPRestoreArray(). Upon restoring the vector data
2576: the data on the host will be marked as out of date. A subsequent
2577: access of the host data will thus incur a data transfer from the
2578: device to the host.
2580: Input Parameter:
2581: . v - the vector
2583: Output Parameter:
2584: . a - the HIP device pointer
2586: Fortran note:
2587: This function is not currently available from Fortran.
2589: Level: intermediate
2591: .seealso: VecHIPRestoreArray(), VecHIPGetArrayRead(), VecHIPGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2592: @*/
2593: PETSC_EXTERN PetscErrorCode VecHIPGetArray(Vec v, PetscScalar **a)
2594: {
2596: #if defined(PETSC_HAVE_HIP)
2597: *a = 0;
2598: VecHIPCopyToGPU(v);
2599: *a = ((Vec_HIP*)v->spptr)->GPUarray;
2600: #endif
2601: return 0;
2602: }
2604: /*@C
2605: VecHIPRestoreArray - Restore a HIP device pointer previously acquired with VecHIPGetArray().
2607: This marks the host data as out of date. Subsequent access to the
2608: vector data on the host side with for instance VecGetArray() incurs a
2609: data transfer.
2611: Input Parameters:
2612: + v - the vector
2613: - a - the HIP device pointer. This pointer is invalid after
2614: VecHIPRestoreArray() returns.
2616: Fortran note:
2617: This function is not currently available from Fortran.
2619: Level: intermediate
2621: .seealso: VecHIPGetArray(), VecHIPGetArrayRead(), VecHIPGetArrayWrite(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2622: @*/
2623: PETSC_EXTERN PetscErrorCode VecHIPRestoreArray(Vec v, PetscScalar **a)
2624: {
2626: #if defined(PETSC_HAVE_HIP)
2627: v->offloadmask = PETSC_OFFLOAD_GPU;
2628: #endif
2630: PetscObjectStateIncrease((PetscObject)v);
2631: return 0;
2632: }
2634: /*@C
2635: VecHIPGetArrayRead - Provides read access to the HIP buffer inside a vector.
2637: This function is analogous to VecGetArrayRead(): The pointer
2638: returned by this function points to a consistent view of the vector
2639: data. This may involve a copy operation of data from the host to the
2640: device if the data on the device is out of date. If the device
2641: memory hasn't been allocated previously it will be allocated as part
2642: of this function call. VecHIPGetArrayRead() assumes that the
2643: user will not modify the vector data. This is analgogous to
2644: intent(in) in Fortran.
2646: The HIP device pointer has to be released by calling
2647: VecHIPRestoreArrayRead(). If the data on the host side was
2648: previously up to date it will remain so, i.e. data on both the device
2649: and the host is up to date. Accessing data on the host side does not
2650: incur a device to host data transfer.
2652: Input Parameter:
2653: . v - the vector
2655: Output Parameter:
2656: . a - the HIP pointer.
2658: Fortran note:
2659: This function is not currently available from Fortran.
2661: Level: intermediate
2663: .seealso: VecHIPRestoreArrayRead(), VecHIPGetArray(), VecHIPGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2664: @*/
2665: PETSC_EXTERN PetscErrorCode VecHIPGetArrayRead(Vec v, const PetscScalar **a)
2666: {
2668: #if defined(PETSC_HAVE_HIP)
2669: *a = 0;
2670: VecHIPCopyToGPU(v);
2671: *a = ((Vec_HIP*)v->spptr)->GPUarray;
2672: #endif
2673: return 0;
2674: }
2676: /*@C
2677: VecHIPRestoreArrayRead - Restore a HIP device pointer previously acquired with VecHIPGetArrayRead().
2679: If the data on the host side was previously up to date it will remain
2680: so, i.e. data on both the device and the host is up to date.
2681: Accessing data on the host side e.g. with VecGetArray() does not
2682: incur a device to host data transfer.
2684: Input Parameters:
2685: + v - the vector
2686: - a - the HIP device pointer. This pointer is invalid after
2687: VecHIPRestoreArrayRead() returns.
2689: Fortran note:
2690: This function is not currently available from Fortran.
2692: Level: intermediate
2694: .seealso: VecHIPGetArrayRead(), VecHIPGetArrayWrite(), VecHIPGetArray(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2695: @*/
2696: PETSC_EXTERN PetscErrorCode VecHIPRestoreArrayRead(Vec v, const PetscScalar **a)
2697: {
2699: *a = NULL;
2700: return 0;
2701: }
2703: /*@C
2704: VecHIPGetArrayWrite - Provides write access to the HIP buffer inside a vector.
2706: The data pointed to by the device pointer is uninitialized. The user
2707: may not read from this data. Furthermore, the entire array needs to
2708: be filled by the user to obtain well-defined behaviour. The device
2709: memory will be allocated by this function if it hasn't been allocated
2710: previously. This is analogous to intent(out) in Fortran.
2712: The device pointer needs to be released with
2713: VecHIPRestoreArrayWrite(). When the pointer is released the
2714: host data of the vector is marked as out of data. Subsequent access
2715: of the host data with e.g. VecGetArray() incurs a device to host data
2716: transfer.
2718: Input Parameter:
2719: . v - the vector
2721: Output Parameter:
2722: . a - the HIP pointer
2724: Fortran note:
2725: This function is not currently available from Fortran.
2727: Level: advanced
2729: .seealso: VecHIPRestoreArrayWrite(), VecHIPGetArray(), VecHIPGetArrayRead(), VecHIPGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2730: @*/
2731: PETSC_EXTERN PetscErrorCode VecHIPGetArrayWrite(Vec v, PetscScalar **a)
2732: {
2734: #if defined(PETSC_HAVE_HIP)
2735: *a = 0;
2736: VecHIPAllocateCheck(v);
2737: *a = ((Vec_HIP*)v->spptr)->GPUarray;
2738: #endif
2739: return 0;
2740: }
2742: /*@C
2743: VecHIPRestoreArrayWrite - Restore a HIP device pointer previously acquired with VecHIPGetArrayWrite().
2745: Data on the host will be marked as out of date. Subsequent access of
2746: the data on the host side e.g. with VecGetArray() will incur a device
2747: to host data transfer.
2749: Input Parameters:
2750: + v - the vector
2751: - a - the HIP device pointer. This pointer is invalid after
2752: VecHIPRestoreArrayWrite() returns.
2754: Fortran note:
2755: This function is not currently available from Fortran.
2757: Level: intermediate
2759: .seealso: VecHIPGetArrayWrite(), VecHIPGetArray(), VecHIPGetArrayRead(), VecHIPGetArrayWrite(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2760: @*/
2761: PETSC_EXTERN PetscErrorCode VecHIPRestoreArrayWrite(Vec v, PetscScalar **a)
2762: {
2764: #if defined(PETSC_HAVE_HIP)
2765: v->offloadmask = PETSC_OFFLOAD_GPU;
2766: #endif
2768: PetscObjectStateIncrease((PetscObject)v);
2769: return 0;
2770: }
2772: /*@C
2773: VecHIPPlaceArray - Allows one to replace the GPU array in a vector with a
2774: GPU array provided by the user. This is useful to avoid copying an
2775: array into a vector.
2777: Not Collective
2779: Input Parameters:
2780: + vec - the vector
2781: - array - the GPU array
2783: Notes:
2784: You can return to the original GPU array with a call to VecHIPResetArray()
2785: It is not possible to use VecHIPPlaceArray() and VecPlaceArray() at the
2786: same time on the same vector.
2788: `vec` does not take ownership of `array` in any way. The user must free `array` themselves
2789: but be careful not to do so before the vector has either been destroyed, had its original
2790: array restored with `VecHIPResetArray()` or permanently replaced with
2791: `VecHIPReplaceArray()`.
2793: Level: developer
2795: .seealso: VecPlaceArray(), VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecResetArray(), VecHIPResetArray(), VecHIPReplaceArray()
2797: @*/
2798: PetscErrorCode VecHIPPlaceArray(Vec vin,const PetscScalar a[])
2799: {
2801: #if defined(PETSC_HAVE_HIP)
2802: VecHIPCopyToGPU(vin);
2804: ((Vec_Seq*)vin->data)->unplacedarray = (PetscScalar *) ((Vec_HIP*)vin->spptr)->GPUarray; /* save previous GPU array so reset can bring it back */
2805: ((Vec_HIP*)vin->spptr)->GPUarray = (PetscScalar*)a;
2806: vin->offloadmask = PETSC_OFFLOAD_GPU;
2807: #endif
2808: PetscObjectStateIncrease((PetscObject)vin);
2809: return 0;
2810: }
2812: /*@C
2813: VecHIPReplaceArray - Allows one to replace the GPU array in a vector
2814: with a GPU array provided by the user. This is useful to avoid copying
2815: a GPU array into a vector.
2817: Not Collective
2819: Input Parameters:
2820: + vec - the vector
2821: - array - the GPU array
2823: Notes:
2824: This permanently replaces the GPU array and frees the memory associated
2825: with the old GPU array.
2827: The memory passed in CANNOT be freed by the user. It will be freed
2828: when the vector is destroyed.
2830: Not supported from Fortran
2832: Level: developer
2834: .seealso: VecGetArray(), VecRestoreArray(), VecPlaceArray(), VecResetArray(), VecHIPResetArray(), VecHIPPlaceArray(), VecReplaceArray()
2836: @*/
2837: PetscErrorCode VecHIPReplaceArray(Vec vin,const PetscScalar a[])
2838: {
2840: #if defined(PETSC_HAVE_HIP)
2841: hipFree(((Vec_HIP*)vin->spptr)->GPUarray);
2842: ((Vec_HIP*)vin->spptr)->GPUarray = (PetscScalar*)a;
2843: vin->offloadmask = PETSC_OFFLOAD_GPU;
2844: #endif
2845: PetscObjectStateIncrease((PetscObject)vin);
2846: return 0;
2847: }
2849: /*@C
2850: VecHIPResetArray - Resets a vector to use its default memory. Call this
2851: after the use of VecHIPPlaceArray().
2853: Not Collective
2855: Input Parameters:
2856: . vec - the vector
2858: Level: developer
2860: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecPlaceArray(), VecResetArray(), VecHIPPlaceArray(), VecHIPReplaceArray()
2862: @*/
2863: PetscErrorCode VecHIPResetArray(Vec vin)
2864: {
2866: #if defined(PETSC_HAVE_HIP)
2867: VecHIPCopyToGPU(vin);
2868: ((Vec_HIP*)vin->spptr)->GPUarray = (PetscScalar *) ((Vec_Seq*)vin->data)->unplacedarray;
2869: ((Vec_Seq*)vin->data)->unplacedarray = 0;
2870: vin->offloadmask = PETSC_OFFLOAD_GPU;
2871: #endif
2872: PetscObjectStateIncrease((PetscObject)vin);
2873: return 0;
2874: }
2876: /*MC
2877: VecDuplicateVecsF90 - Creates several vectors of the same type as an existing vector
2878: and makes them accessible via a Fortran90 pointer.
2880: Synopsis:
2881: VecDuplicateVecsF90(Vec x,PetscInt n,{Vec, pointer :: y(:)},integer ierr)
2883: Collective on Vec
2885: Input Parameters:
2886: + x - a vector to mimic
2887: - n - the number of vectors to obtain
2889: Output Parameters:
2890: + y - Fortran90 pointer to the array of vectors
2891: - ierr - error code
2893: Example of Usage:
2894: .vb
2895: #include <petsc/finclude/petscvec.h>
2896: use petscvec
2898: Vec x
2899: Vec, pointer :: y(:)
2900: ....
2901: call VecDuplicateVecsF90(x,2,y,ierr)
2902: call VecSet(y(2),alpha,ierr)
2903: call VecSet(y(2),alpha,ierr)
2904: ....
2905: call VecDestroyVecsF90(2,y,ierr)
2906: .ve
2908: Notes:
2909: Not yet supported for all F90 compilers
2911: Use VecDestroyVecsF90() to free the space.
2913: Level: beginner
2915: .seealso: VecDestroyVecsF90(), VecDuplicateVecs()
2917: M*/
2919: /*MC
2920: VecRestoreArrayF90 - Restores a vector to a usable state after a call to
2921: VecGetArrayF90().
2923: Synopsis:
2924: VecRestoreArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2926: Logically Collective on Vec
2928: Input Parameters:
2929: + x - vector
2930: - xx_v - the Fortran90 pointer to the array
2932: Output Parameter:
2933: . ierr - error code
2935: Example of Usage:
2936: .vb
2937: #include <petsc/finclude/petscvec.h>
2938: use petscvec
2940: PetscScalar, pointer :: xx_v(:)
2941: ....
2942: call VecGetArrayF90(x,xx_v,ierr)
2943: xx_v(3) = a
2944: call VecRestoreArrayF90(x,xx_v,ierr)
2945: .ve
2947: Level: beginner
2949: .seealso: VecGetArrayF90(), VecGetArray(), VecRestoreArray(), VecRestoreArrayReadF90()
2951: M*/
2953: /*MC
2954: VecDestroyVecsF90 - Frees a block of vectors obtained with VecDuplicateVecsF90().
2956: Synopsis:
2957: VecDestroyVecsF90(PetscInt n,{Vec, pointer :: x(:)},PetscErrorCode ierr)
2959: Collective on Vec
2961: Input Parameters:
2962: + n - the number of vectors previously obtained
2963: - x - pointer to array of vector pointers
2965: Output Parameter:
2966: . ierr - error code
2968: Notes:
2969: Not yet supported for all F90 compilers
2971: Level: beginner
2973: .seealso: VecDestroyVecs(), VecDuplicateVecsF90()
2975: M*/
2977: /*MC
2978: VecGetArrayF90 - Accesses a vector array from Fortran90. For default PETSc
2979: vectors, VecGetArrayF90() returns a pointer to the local data array. Otherwise,
2980: this routine is implementation dependent. You MUST call VecRestoreArrayF90()
2981: when you no longer need access to the array.
2983: Synopsis:
2984: VecGetArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2986: Logically Collective on Vec
2988: Input Parameter:
2989: . x - vector
2991: Output Parameters:
2992: + xx_v - the Fortran90 pointer to the array
2993: - ierr - error code
2995: Example of Usage:
2996: .vb
2997: #include <petsc/finclude/petscvec.h>
2998: use petscvec
3000: PetscScalar, pointer :: xx_v(:)
3001: ....
3002: call VecGetArrayF90(x,xx_v,ierr)
3003: xx_v(3) = a
3004: call VecRestoreArrayF90(x,xx_v,ierr)
3005: .ve
3007: If you ONLY intend to read entries from the array and not change any entries you should use VecGetArrayReadF90().
3009: Level: beginner
3011: .seealso: VecRestoreArrayF90(), VecGetArray(), VecRestoreArray(), VecGetArrayReadF90()
3013: M*/
3015: /*MC
3016: VecGetArrayReadF90 - Accesses a read only 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 VecRestoreArrayReadF90()
3019: when you no longer need access to the array.
3021: Synopsis:
3022: VecGetArrayReadF90(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 VecGetArrayReadF90(x,xx_v,ierr)
3041: a = xx_v(3)
3042: call VecRestoreArrayReadF90(x,xx_v,ierr)
3043: .ve
3045: If you intend to write entries into the array you must use VecGetArrayF90().
3047: Level: beginner
3049: .seealso: VecRestoreArrayReadF90(), VecGetArray(), VecRestoreArray(), VecGetArrayRead(), VecRestoreArrayRead(), VecGetArrayF90()
3051: M*/
3053: /*MC
3054: VecRestoreArrayReadF90 - Restores a readonly vector to a usable state after a call to
3055: VecGetArrayReadF90().
3057: Synopsis:
3058: VecRestoreArrayReadF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
3060: Logically Collective on Vec
3062: Input Parameters:
3063: + x - vector
3064: - xx_v - the Fortran90 pointer to the array
3066: Output Parameter:
3067: . ierr - error code
3069: Example of Usage:
3070: .vb
3071: #include <petsc/finclude/petscvec.h>
3072: use petscvec
3074: PetscScalar, pointer :: xx_v(:)
3075: ....
3076: call VecGetArrayReadF90(x,xx_v,ierr)
3077: a = xx_v(3)
3078: call VecRestoreArrayReadF90(x,xx_v,ierr)
3079: .ve
3081: Level: beginner
3083: .seealso: VecGetArrayReadF90(), VecGetArray(), VecRestoreArray(), VecGetArrayRead(), VecRestoreArrayRead(), VecRestoreArrayF90()
3085: M*/
3087: /*@C
3088: VecGetArray2d - Returns a pointer to a 2d contiguous array that contains this
3089: processor's portion of the vector data. You MUST call VecRestoreArray2d()
3090: when you no longer need access to the array.
3092: Logically Collective
3094: Input Parameters:
3095: + x - the vector
3096: . m - first dimension of two dimensional array
3097: . n - second dimension of two dimensional array
3098: . mstart - first index you will use in first coordinate direction (often 0)
3099: - nstart - first index in the second coordinate direction (often 0)
3101: Output Parameter:
3102: . a - location to put pointer to the array
3104: Level: developer
3106: Notes:
3107: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
3108: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3109: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3110: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
3112: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3114: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3115: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3116: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3117: @*/
3118: PetscErrorCode VecGetArray2d(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3119: {
3120: PetscInt i,N;
3121: PetscScalar *aa;
3126: VecGetLocalSize(x,&N);
3128: VecGetArray(x,&aa);
3130: PetscMalloc1(m,a);
3131: for (i=0; i<m; i++) (*a)[i] = aa + i*n - nstart;
3132: *a -= mstart;
3133: return 0;
3134: }
3136: /*@C
3137: VecGetArray2dWrite - Returns a pointer to a 2d contiguous array that will contain this
3138: processor's portion of the vector data. You MUST call VecRestoreArray2dWrite()
3139: when you no longer need access to the array.
3141: Logically Collective
3143: Input Parameters:
3144: + x - the vector
3145: . m - first dimension of two dimensional array
3146: . n - second dimension of two dimensional array
3147: . mstart - first index you will use in first coordinate direction (often 0)
3148: - nstart - first index in the second coordinate direction (often 0)
3150: Output Parameter:
3151: . a - location to put pointer to the array
3153: Level: developer
3155: Notes:
3156: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
3157: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3158: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3159: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
3161: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3163: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3164: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3165: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3166: @*/
3167: PetscErrorCode VecGetArray2dWrite(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3168: {
3169: PetscInt i,N;
3170: PetscScalar *aa;
3175: VecGetLocalSize(x,&N);
3177: VecGetArrayWrite(x,&aa);
3179: PetscMalloc1(m,a);
3180: for (i=0; i<m; i++) (*a)[i] = aa + i*n - nstart;
3181: *a -= mstart;
3182: return 0;
3183: }
3185: /*@C
3186: VecRestoreArray2d - Restores a vector after VecGetArray2d() has been called.
3188: Logically Collective
3190: Input Parameters:
3191: + x - the vector
3192: . m - first dimension of two dimensional array
3193: . n - second dimension of the two dimensional array
3194: . mstart - first index you will use in first coordinate direction (often 0)
3195: . nstart - first index in the second coordinate direction (often 0)
3196: - a - location of pointer to array obtained from VecGetArray2d()
3198: Level: developer
3200: Notes:
3201: For regular PETSc vectors this routine does not involve any copies. For
3202: any special vectors that do not store local vector data in a contiguous
3203: array, this routine will copy the data back into the underlying
3204: vector data structure from the array obtained with VecGetArray().
3206: This routine actually zeros out the a pointer.
3208: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3209: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3210: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3211: @*/
3212: PetscErrorCode VecRestoreArray2d(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3213: {
3214: void *dummy;
3219: dummy = (void*)(*a + mstart);
3220: PetscFree(dummy);
3221: VecRestoreArray(x,NULL);
3222: return 0;
3223: }
3225: /*@C
3226: VecRestoreArray2dWrite - Restores a vector after VecGetArray2dWrite() has been called.
3228: Logically Collective
3230: Input Parameters:
3231: + x - the vector
3232: . m - first dimension of two dimensional array
3233: . n - second dimension of the two dimensional array
3234: . mstart - first index you will use in first coordinate direction (often 0)
3235: . nstart - first index in the second coordinate direction (often 0)
3236: - a - location of pointer to array obtained from VecGetArray2d()
3238: Level: developer
3240: Notes:
3241: For regular PETSc vectors this routine does not involve any copies. For
3242: any special vectors that do not store local vector data in a contiguous
3243: array, this routine will copy the data back into the underlying
3244: vector data structure from the array obtained with VecGetArray().
3246: This routine actually zeros out the a pointer.
3248: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3249: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3250: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3251: @*/
3252: PetscErrorCode VecRestoreArray2dWrite(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3253: {
3254: void *dummy;
3259: dummy = (void*)(*a + mstart);
3260: PetscFree(dummy);
3261: VecRestoreArrayWrite(x,NULL);
3262: return 0;
3263: }
3265: /*@C
3266: VecGetArray1d - Returns a pointer to a 1d contiguous array that contains this
3267: processor's portion of the vector data. You MUST call VecRestoreArray1d()
3268: when you no longer need access to the array.
3270: Logically Collective
3272: Input Parameters:
3273: + x - the vector
3274: . m - first dimension of two dimensional array
3275: - mstart - first index you will use in first coordinate direction (often 0)
3277: Output Parameter:
3278: . a - location to put pointer to the array
3280: Level: developer
3282: Notes:
3283: For a vector obtained from DMCreateLocalVector() mstart are likely
3284: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3285: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
3287: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3289: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3290: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3291: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3292: @*/
3293: PetscErrorCode VecGetArray1d(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
3294: {
3295: PetscInt N;
3300: VecGetLocalSize(x,&N);
3302: VecGetArray(x,a);
3303: *a -= mstart;
3304: return 0;
3305: }
3307: /*@C
3308: VecGetArray1dWrite - Returns a pointer to a 1d contiguous array that will contain this
3309: processor's portion of the vector data. You MUST call VecRestoreArray1dWrite()
3310: when you no longer need access to the array.
3312: Logically Collective
3314: Input Parameters:
3315: + x - the vector
3316: . m - first dimension of two dimensional array
3317: - mstart - first index you will use in first coordinate direction (often 0)
3319: Output Parameter:
3320: . a - location to put pointer to the array
3322: Level: developer
3324: Notes:
3325: For a vector obtained from DMCreateLocalVector() mstart are likely
3326: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3327: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
3329: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3331: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3332: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3333: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3334: @*/
3335: PetscErrorCode VecGetArray1dWrite(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
3336: {
3337: PetscInt N;
3342: VecGetLocalSize(x,&N);
3344: VecGetArrayWrite(x,a);
3345: *a -= mstart;
3346: return 0;
3347: }
3349: /*@C
3350: VecRestoreArray1d - Restores a vector after VecGetArray1d() has been called.
3352: Logically Collective
3354: Input Parameters:
3355: + x - the vector
3356: . m - first dimension of two dimensional array
3357: . mstart - first index you will use in first coordinate direction (often 0)
3358: - a - location of pointer to array obtained from VecGetArray21()
3360: Level: developer
3362: Notes:
3363: For regular PETSc vectors this routine does not involve any copies. For
3364: any special vectors that do not store local vector data in a contiguous
3365: array, this routine will copy the data back into the underlying
3366: vector data structure from the array obtained with VecGetArray1d().
3368: This routine actually zeros out the a pointer.
3370: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3371: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3372: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
3373: @*/
3374: PetscErrorCode VecRestoreArray1d(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
3375: {
3378: VecRestoreArray(x,NULL);
3379: return 0;
3380: }
3382: /*@C
3383: VecRestoreArray1dWrite - Restores a vector after VecGetArray1dWrite() has been called.
3385: Logically Collective
3387: Input Parameters:
3388: + x - the vector
3389: . m - first dimension of two dimensional array
3390: . mstart - first index you will use in first coordinate direction (often 0)
3391: - a - location of pointer to array obtained from VecGetArray21()
3393: Level: developer
3395: Notes:
3396: For regular PETSc vectors this routine does not involve any copies. For
3397: any special vectors that do not store local vector data in a contiguous
3398: array, this routine will copy the data back into the underlying
3399: vector data structure from the array obtained with VecGetArray1d().
3401: This routine actually zeros out the a pointer.
3403: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3404: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3405: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
3406: @*/
3407: PetscErrorCode VecRestoreArray1dWrite(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
3408: {
3411: VecRestoreArrayWrite(x,NULL);
3412: return 0;
3413: }
3415: /*@C
3416: VecGetArray3d - Returns a pointer to a 3d contiguous array that contains this
3417: processor's portion of the vector data. You MUST call VecRestoreArray3d()
3418: when you no longer need access to the array.
3420: Logically Collective
3422: Input Parameters:
3423: + x - the vector
3424: . m - first dimension of three dimensional array
3425: . n - second dimension of three dimensional array
3426: . p - third dimension of three dimensional array
3427: . mstart - first index you will use in first coordinate direction (often 0)
3428: . nstart - first index in the second coordinate direction (often 0)
3429: - pstart - first index in the third coordinate direction (often 0)
3431: Output Parameter:
3432: . a - location to put pointer to the array
3434: Level: developer
3436: Notes:
3437: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3438: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3439: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3440: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3442: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3444: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3445: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3446: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3447: @*/
3448: PetscErrorCode VecGetArray3d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3449: {
3450: PetscInt i,N,j;
3451: PetscScalar *aa,**b;
3456: VecGetLocalSize(x,&N);
3458: VecGetArray(x,&aa);
3460: PetscMalloc(m*sizeof(PetscScalar**)+m*n*sizeof(PetscScalar*),a);
3461: b = (PetscScalar**)((*a) + m);
3462: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3463: for (i=0; i<m; i++)
3464: for (j=0; j<n; j++)
3465: b[i*n+j] = aa + i*n*p + j*p - pstart;
3466: *a -= mstart;
3467: return 0;
3468: }
3470: /*@C
3471: VecGetArray3dWrite - Returns a pointer to a 3d contiguous array that will contain this
3472: processor's portion of the vector data. You MUST call VecRestoreArray3dWrite()
3473: when you no longer need access to the array.
3475: Logically Collective
3477: Input Parameters:
3478: + x - the vector
3479: . m - first dimension of three dimensional array
3480: . n - second dimension of three dimensional array
3481: . p - third dimension of three dimensional array
3482: . mstart - first index you will use in first coordinate direction (often 0)
3483: . nstart - first index in the second coordinate direction (often 0)
3484: - pstart - first index in the third coordinate direction (often 0)
3486: Output Parameter:
3487: . a - location to put pointer to the array
3489: Level: developer
3491: Notes:
3492: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3493: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3494: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3495: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3497: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3499: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3500: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3501: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3502: @*/
3503: PetscErrorCode VecGetArray3dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3504: {
3505: PetscInt i,N,j;
3506: PetscScalar *aa,**b;
3511: VecGetLocalSize(x,&N);
3513: VecGetArrayWrite(x,&aa);
3515: PetscMalloc(m*sizeof(PetscScalar**)+m*n*sizeof(PetscScalar*),a);
3516: b = (PetscScalar**)((*a) + m);
3517: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3518: for (i=0; i<m; i++)
3519: for (j=0; j<n; j++)
3520: b[i*n+j] = aa + i*n*p + j*p - pstart;
3522: *a -= mstart;
3523: return 0;
3524: }
3526: /*@C
3527: VecRestoreArray3d - Restores a vector after VecGetArray3d() has been called.
3529: Logically Collective
3531: Input Parameters:
3532: + x - the vector
3533: . m - first dimension of three dimensional array
3534: . n - second dimension of the three dimensional array
3535: . p - third dimension of the three dimensional array
3536: . mstart - first index you will use in first coordinate direction (often 0)
3537: . nstart - first index in the second coordinate direction (often 0)
3538: . pstart - first index in the third coordinate direction (often 0)
3539: - a - location of pointer to array obtained from VecGetArray3d()
3541: Level: developer
3543: Notes:
3544: For regular PETSc vectors this routine does not involve any copies. For
3545: any special vectors that do not store local vector data in a contiguous
3546: array, this routine will copy the data back into the underlying
3547: vector data structure from the array obtained with VecGetArray().
3549: This routine actually zeros out the a pointer.
3551: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3552: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3553: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3554: @*/
3555: PetscErrorCode VecRestoreArray3d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3556: {
3557: void *dummy;
3562: dummy = (void*)(*a + mstart);
3563: PetscFree(dummy);
3564: VecRestoreArray(x,NULL);
3565: return 0;
3566: }
3568: /*@C
3569: VecRestoreArray3dWrite - Restores a vector after VecGetArray3dWrite() has been called.
3571: Logically Collective
3573: Input Parameters:
3574: + x - the vector
3575: . m - first dimension of three dimensional array
3576: . n - second dimension of the three dimensional array
3577: . p - third dimension of the three dimensional array
3578: . mstart - first index you will use in first coordinate direction (often 0)
3579: . nstart - first index in the second coordinate direction (often 0)
3580: . pstart - first index in the third coordinate direction (often 0)
3581: - a - location of pointer to array obtained from VecGetArray3d()
3583: Level: developer
3585: Notes:
3586: For regular PETSc vectors this routine does not involve any copies. For
3587: any special vectors that do not store local vector data in a contiguous
3588: array, this routine will copy the data back into the underlying
3589: vector data structure from the array obtained with VecGetArray().
3591: This routine actually zeros out the a pointer.
3593: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3594: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3595: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3596: @*/
3597: PetscErrorCode VecRestoreArray3dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3598: {
3599: void *dummy;
3604: dummy = (void*)(*a + mstart);
3605: PetscFree(dummy);
3606: VecRestoreArrayWrite(x,NULL);
3607: return 0;
3608: }
3610: /*@C
3611: VecGetArray4d - Returns a pointer to a 4d contiguous array that contains this
3612: processor's portion of the vector data. You MUST call VecRestoreArray4d()
3613: when you no longer need access to the array.
3615: Logically Collective
3617: Input Parameters:
3618: + x - the vector
3619: . m - first dimension of four dimensional array
3620: . n - second dimension of four dimensional array
3621: . p - third dimension of four dimensional array
3622: . q - fourth dimension of four dimensional array
3623: . mstart - first index you will use in first coordinate direction (often 0)
3624: . nstart - first index in the second coordinate direction (often 0)
3625: . pstart - first index in the third coordinate direction (often 0)
3626: - qstart - first index in the fourth coordinate direction (often 0)
3628: Output Parameter:
3629: . a - location to put pointer to the array
3631: Level: beginner
3633: Notes:
3634: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3635: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3636: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3637: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3639: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3641: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3642: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3643: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3644: @*/
3645: PetscErrorCode VecGetArray4d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3646: {
3647: PetscInt i,N,j,k;
3648: PetscScalar *aa,***b,**c;
3653: VecGetLocalSize(x,&N);
3655: VecGetArray(x,&aa);
3657: PetscMalloc(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p*sizeof(PetscScalar*),a);
3658: b = (PetscScalar***)((*a) + m);
3659: c = (PetscScalar**)(b + m*n);
3660: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3661: for (i=0; i<m; i++)
3662: for (j=0; j<n; j++)
3663: b[i*n+j] = c + i*n*p + j*p - pstart;
3664: for (i=0; i<m; i++)
3665: for (j=0; j<n; j++)
3666: for (k=0; k<p; k++)
3667: c[i*n*p+j*p+k] = aa + i*n*p*q + j*p*q + k*q - qstart;
3668: *a -= mstart;
3669: return 0;
3670: }
3672: /*@C
3673: VecGetArray4dWrite - Returns a pointer to a 4d contiguous array that will contain this
3674: processor's portion of the vector data. You MUST call VecRestoreArray4dWrite()
3675: when you no longer need access to the array.
3677: Logically Collective
3679: Input Parameters:
3680: + x - the vector
3681: . m - first dimension of four dimensional array
3682: . n - second dimension of four dimensional array
3683: . p - third dimension of four dimensional array
3684: . q - fourth dimension of four dimensional array
3685: . mstart - first index you will use in first coordinate direction (often 0)
3686: . nstart - first index in the second coordinate direction (often 0)
3687: . pstart - first index in the third coordinate direction (often 0)
3688: - qstart - first index in the fourth coordinate direction (often 0)
3690: Output Parameter:
3691: . a - location to put pointer to the array
3693: Level: beginner
3695: Notes:
3696: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3697: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3698: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3699: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3701: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3703: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3704: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3705: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3706: @*/
3707: PetscErrorCode VecGetArray4dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3708: {
3709: PetscInt i,N,j,k;
3710: PetscScalar *aa,***b,**c;
3715: VecGetLocalSize(x,&N);
3717: VecGetArrayWrite(x,&aa);
3719: PetscMalloc(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p*sizeof(PetscScalar*),a);
3720: b = (PetscScalar***)((*a) + m);
3721: c = (PetscScalar**)(b + m*n);
3722: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3723: for (i=0; i<m; i++)
3724: for (j=0; j<n; j++)
3725: b[i*n+j] = c + i*n*p + j*p - pstart;
3726: for (i=0; i<m; i++)
3727: for (j=0; j<n; j++)
3728: for (k=0; k<p; k++)
3729: c[i*n*p+j*p+k] = aa + i*n*p*q + j*p*q + k*q - qstart;
3730: *a -= mstart;
3731: return 0;
3732: }
3734: /*@C
3735: VecRestoreArray4d - Restores a vector after VecGetArray3d() has been called.
3737: Logically Collective
3739: Input Parameters:
3740: + x - the vector
3741: . m - first dimension of four dimensional array
3742: . n - second dimension of the four dimensional array
3743: . p - third dimension of the four dimensional array
3744: . q - fourth dimension of the four dimensional array
3745: . mstart - first index you will use in first coordinate direction (often 0)
3746: . nstart - first index in the second coordinate direction (often 0)
3747: . pstart - first index in the third coordinate direction (often 0)
3748: . qstart - first index in the fourth coordinate direction (often 0)
3749: - a - location of pointer to array obtained from VecGetArray4d()
3751: Level: beginner
3753: Notes:
3754: For regular PETSc vectors this routine does not involve any copies. For
3755: any special vectors that do not store local vector data in a contiguous
3756: array, this routine will copy the data back into the underlying
3757: vector data structure from the array obtained with VecGetArray().
3759: This routine actually zeros out the a pointer.
3761: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3762: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3763: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3764: @*/
3765: PetscErrorCode VecRestoreArray4d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3766: {
3767: void *dummy;
3772: dummy = (void*)(*a + mstart);
3773: PetscFree(dummy);
3774: VecRestoreArray(x,NULL);
3775: return 0;
3776: }
3778: /*@C
3779: VecRestoreArray4dWrite - Restores a vector after VecGetArray3dWrite() has been called.
3781: Logically Collective
3783: Input Parameters:
3784: + x - the vector
3785: . m - first dimension of four dimensional array
3786: . n - second dimension of the four dimensional array
3787: . p - third dimension of the four dimensional array
3788: . q - fourth dimension of the four dimensional array
3789: . mstart - first index you will use in first coordinate direction (often 0)
3790: . nstart - first index in the second coordinate direction (often 0)
3791: . pstart - first index in the third coordinate direction (often 0)
3792: . qstart - first index in the fourth coordinate direction (often 0)
3793: - a - location of pointer to array obtained from VecGetArray4d()
3795: Level: beginner
3797: Notes:
3798: For regular PETSc vectors this routine does not involve any copies. For
3799: any special vectors that do not store local vector data in a contiguous
3800: array, this routine will copy the data back into the underlying
3801: vector data structure from the array obtained with VecGetArray().
3803: This routine actually zeros out the a pointer.
3805: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3806: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3807: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3808: @*/
3809: PetscErrorCode VecRestoreArray4dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3810: {
3811: void *dummy;
3816: dummy = (void*)(*a + mstart);
3817: PetscFree(dummy);
3818: VecRestoreArrayWrite(x,NULL);
3819: return 0;
3820: }
3822: /*@C
3823: VecGetArray2dRead - Returns a pointer to a 2d contiguous array that contains this
3824: processor's portion of the vector data. You MUST call VecRestoreArray2dRead()
3825: when you no longer need access to the array.
3827: Logically Collective
3829: Input Parameters:
3830: + x - the vector
3831: . m - first dimension of two dimensional array
3832: . n - second dimension of two dimensional array
3833: . mstart - first index you will use in first coordinate direction (often 0)
3834: - nstart - first index in the second coordinate direction (often 0)
3836: Output Parameter:
3837: . a - location to put pointer to the array
3839: Level: developer
3841: Notes:
3842: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
3843: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3844: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3845: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
3847: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3849: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3850: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3851: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3852: @*/
3853: PetscErrorCode VecGetArray2dRead(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3854: {
3855: PetscInt i,N;
3856: const PetscScalar *aa;
3861: VecGetLocalSize(x,&N);
3863: VecGetArrayRead(x,&aa);
3865: PetscMalloc1(m,a);
3866: for (i=0; i<m; i++) (*a)[i] = (PetscScalar*) aa + i*n - nstart;
3867: *a -= mstart;
3868: return 0;
3869: }
3871: /*@C
3872: VecRestoreArray2dRead - Restores a vector after VecGetArray2dRead() has been called.
3874: Logically Collective
3876: Input Parameters:
3877: + x - the vector
3878: . m - first dimension of two dimensional array
3879: . n - second dimension of the two dimensional array
3880: . mstart - first index you will use in first coordinate direction (often 0)
3881: . nstart - first index in the second coordinate direction (often 0)
3882: - a - location of pointer to array obtained from VecGetArray2d()
3884: Level: developer
3886: Notes:
3887: For regular PETSc vectors this routine does not involve any copies. For
3888: any special vectors that do not store local vector data in a contiguous
3889: array, this routine will copy the data back into the underlying
3890: vector data structure from the array obtained with VecGetArray().
3892: This routine actually zeros out the a pointer.
3894: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3895: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3896: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3897: @*/
3898: PetscErrorCode VecRestoreArray2dRead(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3899: {
3900: void *dummy;
3905: dummy = (void*)(*a + mstart);
3906: PetscFree(dummy);
3907: VecRestoreArrayRead(x,NULL);
3908: return 0;
3909: }
3911: /*@C
3912: VecGetArray1dRead - Returns a pointer to a 1d contiguous array that contains this
3913: processor's portion of the vector data. You MUST call VecRestoreArray1dRead()
3914: when you no longer need access to the array.
3916: Logically Collective
3918: Input Parameters:
3919: + x - the vector
3920: . m - first dimension of two dimensional array
3921: - mstart - first index you will use in first coordinate direction (often 0)
3923: Output Parameter:
3924: . a - location to put pointer to the array
3926: Level: developer
3928: Notes:
3929: For a vector obtained from DMCreateLocalVector() mstart are likely
3930: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3931: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
3933: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3935: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3936: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3937: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3938: @*/
3939: PetscErrorCode VecGetArray1dRead(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
3940: {
3941: PetscInt N;
3946: VecGetLocalSize(x,&N);
3948: VecGetArrayRead(x,(const PetscScalar**)a);
3949: *a -= mstart;
3950: return 0;
3951: }
3953: /*@C
3954: VecRestoreArray1dRead - Restores a vector after VecGetArray1dRead() has been called.
3956: Logically Collective
3958: Input Parameters:
3959: + x - the vector
3960: . m - first dimension of two dimensional array
3961: . mstart - first index you will use in first coordinate direction (often 0)
3962: - a - location of pointer to array obtained from VecGetArray21()
3964: Level: developer
3966: Notes:
3967: For regular PETSc vectors this routine does not involve any copies. For
3968: any special vectors that do not store local vector data in a contiguous
3969: array, this routine will copy the data back into the underlying
3970: vector data structure from the array obtained with VecGetArray1dRead().
3972: This routine actually zeros out the a pointer.
3974: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3975: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3976: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
3977: @*/
3978: PetscErrorCode VecRestoreArray1dRead(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
3979: {
3982: VecRestoreArrayRead(x,NULL);
3983: return 0;
3984: }
3986: /*@C
3987: VecGetArray3dRead - Returns a pointer to a 3d contiguous array that contains this
3988: processor's portion of the vector data. You MUST call VecRestoreArray3dRead()
3989: when you no longer need access to the array.
3991: Logically Collective
3993: Input Parameters:
3994: + x - the vector
3995: . m - first dimension of three dimensional array
3996: . n - second dimension of three dimensional array
3997: . p - third dimension of three dimensional array
3998: . mstart - first index you will use in first coordinate direction (often 0)
3999: . nstart - first index in the second coordinate direction (often 0)
4000: - pstart - first index in the third coordinate direction (often 0)
4002: Output Parameter:
4003: . a - location to put pointer to the array
4005: Level: developer
4007: Notes:
4008: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
4009: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
4010: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
4011: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3dRead().
4013: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
4015: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
4016: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
4017: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
4018: @*/
4019: PetscErrorCode VecGetArray3dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
4020: {
4021: PetscInt i,N,j;
4022: const PetscScalar *aa;
4023: PetscScalar **b;
4028: VecGetLocalSize(x,&N);
4030: VecGetArrayRead(x,&aa);
4032: PetscMalloc(m*sizeof(PetscScalar**)+m*n*sizeof(PetscScalar*),a);
4033: b = (PetscScalar**)((*a) + m);
4034: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
4035: for (i=0; i<m; i++)
4036: for (j=0; j<n; j++)
4037: b[i*n+j] = (PetscScalar *)aa + i*n*p + j*p - pstart;
4038: *a -= mstart;
4039: return 0;
4040: }
4042: /*@C
4043: VecRestoreArray3dRead - Restores a vector after VecGetArray3dRead() has been called.
4045: Logically Collective
4047: Input Parameters:
4048: + x - the vector
4049: . m - first dimension of three dimensional array
4050: . n - second dimension of the three dimensional array
4051: . p - third dimension of the three dimensional array
4052: . mstart - first index you will use in first coordinate direction (often 0)
4053: . nstart - first index in the second coordinate direction (often 0)
4054: . pstart - first index in the third coordinate direction (often 0)
4055: - a - location of pointer to array obtained from VecGetArray3dRead()
4057: Level: developer
4059: Notes:
4060: For regular PETSc vectors this routine does not involve any copies. For
4061: any special vectors that do not store local vector data in a contiguous
4062: array, this routine will copy the data back into the underlying
4063: vector data structure from the array obtained with VecGetArray().
4065: This routine actually zeros out the a pointer.
4067: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
4068: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
4069: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
4070: @*/
4071: PetscErrorCode VecRestoreArray3dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
4072: {
4073: void *dummy;
4078: dummy = (void*)(*a + mstart);
4079: PetscFree(dummy);
4080: VecRestoreArrayRead(x,NULL);
4081: return 0;
4082: }
4084: /*@C
4085: VecGetArray4dRead - Returns a pointer to a 4d contiguous array that contains this
4086: processor's portion of the vector data. You MUST call VecRestoreArray4dRead()
4087: when you no longer need access to the array.
4089: Logically Collective
4091: Input Parameters:
4092: + x - the vector
4093: . m - first dimension of four dimensional array
4094: . n - second dimension of four dimensional array
4095: . p - third dimension of four dimensional array
4096: . q - fourth dimension of four dimensional array
4097: . mstart - first index you will use in first coordinate direction (often 0)
4098: . nstart - first index in the second coordinate direction (often 0)
4099: . pstart - first index in the third coordinate direction (often 0)
4100: - qstart - first index in the fourth coordinate direction (often 0)
4102: Output Parameter:
4103: . a - location to put pointer to the array
4105: Level: beginner
4107: Notes:
4108: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
4109: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
4110: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
4111: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
4113: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
4115: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
4116: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
4117: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
4118: @*/
4119: PetscErrorCode VecGetArray4dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
4120: {
4121: PetscInt i,N,j,k;
4122: const PetscScalar *aa;
4123: PetscScalar ***b,**c;
4128: VecGetLocalSize(x,&N);
4130: VecGetArrayRead(x,&aa);
4132: PetscMalloc(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p*sizeof(PetscScalar*),a);
4133: b = (PetscScalar***)((*a) + m);
4134: c = (PetscScalar**)(b + m*n);
4135: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
4136: for (i=0; i<m; i++)
4137: for (j=0; j<n; j++)
4138: b[i*n+j] = c + i*n*p + j*p - pstart;
4139: for (i=0; i<m; i++)
4140: for (j=0; j<n; j++)
4141: for (k=0; k<p; k++)
4142: c[i*n*p+j*p+k] = (PetscScalar*) aa + i*n*p*q + j*p*q + k*q - qstart;
4143: *a -= mstart;
4144: return 0;
4145: }
4147: /*@C
4148: VecRestoreArray4dRead - Restores a vector after VecGetArray3d() has been called.
4150: Logically Collective
4152: Input Parameters:
4153: + x - the vector
4154: . m - first dimension of four dimensional array
4155: . n - second dimension of the four dimensional array
4156: . p - third dimension of the four dimensional array
4157: . q - fourth dimension of the four dimensional array
4158: . mstart - first index you will use in first coordinate direction (often 0)
4159: . nstart - first index in the second coordinate direction (often 0)
4160: . pstart - first index in the third coordinate direction (often 0)
4161: . qstart - first index in the fourth coordinate direction (often 0)
4162: - a - location of pointer to array obtained from VecGetArray4dRead()
4164: Level: beginner
4166: Notes:
4167: For regular PETSc vectors this routine does not involve any copies. For
4168: any special vectors that do not store local vector data in a contiguous
4169: array, this routine will copy the data back into the underlying
4170: vector data structure from the array obtained with VecGetArray().
4172: This routine actually zeros out the a pointer.
4174: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
4175: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
4176: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
4177: @*/
4178: PetscErrorCode VecRestoreArray4dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
4179: {
4180: void *dummy;
4185: dummy = (void*)(*a + mstart);
4186: PetscFree(dummy);
4187: VecRestoreArrayRead(x,NULL);
4188: return 0;
4189: }
4191: #if defined(PETSC_USE_DEBUG)
4193: /*@
4194: VecLockGet - Gets the current lock status of a vector
4196: Logically Collective on Vec
4198: Input Parameter:
4199: . x - the vector
4201: Output Parameter:
4202: . state - greater than zero indicates the vector is locked for read; less then zero indicates the vector is
4203: locked for write; equal to zero means the vector is unlocked, that is, it is free to read or write.
4205: Level: beginner
4207: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPush(), VecLockReadPop()
4208: @*/
4209: PetscErrorCode VecLockGet(Vec x,PetscInt *state)
4210: {
4212: *state = x->lock;
4213: return 0;
4214: }
4216: /*@
4217: VecLockReadPush - Pushes a read-only lock on a vector to prevent it from writing
4219: Logically Collective on Vec
4221: Input Parameter:
4222: . x - the vector
4224: Notes:
4225: If this is set then calls to VecGetArray() or VecSetValues() or any other routines that change the vectors values will fail.
4227: The call can be nested, i.e., called multiple times on the same vector, but each VecLockReadPush(x) has to have one matching
4228: VecLockReadPop(x), which removes the latest read-only lock.
4230: Level: beginner
4232: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPop(), VecLockGet()
4233: @*/
4234: PetscErrorCode VecLockReadPush(Vec x)
4235: {
4238: x->lock++;
4239: return 0;
4240: }
4242: /*@
4243: VecLockReadPop - Pops a read-only lock from a vector
4245: Logically Collective on Vec
4247: Input Parameter:
4248: . x - the vector
4250: Level: beginner
4252: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPush(), VecLockGet()
4253: @*/
4254: PetscErrorCode VecLockReadPop(Vec x)
4255: {
4257: x->lock--;
4259: return 0;
4260: }
4262: /*@C
4263: VecLockWriteSet_Private - Lock or unlock a vector for exclusive read/write access
4265: Logically Collective on Vec
4267: Input Parameters:
4268: + x - the vector
4269: - flg - PETSC_TRUE to lock the vector for writing; PETSC_FALSE to unlock it.
4271: Notes:
4272: The function is usefull in split-phase computations, which usually have a begin phase and an end phase.
4273: One can call VecLockWriteSet_Private(x,PETSC_TRUE) in the begin phase to lock a vector for exclusive
4274: access, and call VecLockWriteSet_Private(x,PETSC_FALSE) in the end phase to unlock the vector from exclusive
4275: access. In this way, one is ensured no other operations can access the vector in between. The code may like
4277: VecGetArray(x,&xdata); // begin phase
4278: VecLockWriteSet_Private(v,PETSC_TRUE);
4280: Other operations, which can not acceess x anymore (they can access xdata, of course)
4282: VecRestoreArray(x,&vdata); // end phase
4283: VecLockWriteSet_Private(v,PETSC_FALSE);
4285: The call can not be nested on the same vector, in other words, one can not call VecLockWriteSet_Private(x,PETSC_TRUE)
4286: again before calling VecLockWriteSet_Private(v,PETSC_FALSE).
4288: Level: beginner
4290: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPush(), VecLockReadPop(), VecLockGet()
4291: @*/
4292: PetscErrorCode VecLockWriteSet_Private(Vec x,PetscBool flg)
4293: {
4295: if (flg) {
4298: else x->lock = -1;
4299: } else {
4301: x->lock = 0;
4302: }
4303: return 0;
4304: }
4306: /*@
4307: VecLockPush - Pushes a read-only lock on a vector to prevent it from writing
4309: Level: deprecated
4311: .seealso: VecLockReadPush()
4312: @*/
4313: PetscErrorCode VecLockPush(Vec x)
4314: {
4315: VecLockReadPush(x);
4316: return 0;
4317: }
4319: /*@
4320: VecLockPop - Pops a read-only lock from a vector
4322: Level: deprecated
4324: .seealso: VecLockReadPop()
4325: @*/
4326: PetscErrorCode VecLockPop(Vec x)
4327: {
4328: VecLockReadPop(x);
4329: return 0;
4330: }
4332: #endif