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: #if PetscDefined(USE_DEBUG)
19: // this is a no-op '0' macro in optimized builds
20: PetscErrorCode VecValidValues_Internal(Vec vec, PetscInt argnum, PetscBool begin)
21: {
22: if (vec->petscnative || vec->ops->getarray) {
23: PetscInt n;
24: const PetscScalar *x;
25: PetscOffloadMask mask;
27: VecGetOffloadMask(vec, &mask);
28: if (!PetscOffloadHost(mask)) return 0;
29: VecGetLocalSize(vec, &n);
30: VecGetArrayRead(vec, &x);
31: for (PetscInt i = 0; i < n; i++) {
32: if (begin) {
34: } else {
36: }
37: }
38: VecRestoreArrayRead(vec, &x);
39: }
40: return 0;
41: }
42: #endif
44: /*@
45: VecMaxPointwiseDivide - Computes the maximum of the componentwise division max = max_i abs(x_i/y_i).
47: Logically Collective on x
49: Input Parameters:
50: . x, y - the vectors
52: Output Parameter:
53: . max - the result
55: Level: advanced
57: Notes:
58: 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: `Vec`, `VecPointwiseDivide()`, `VecPointwiseMult()`, `VecPointwiseMax()`, `VecPointwiseMin()`, `VecPointwiseMaxAbs()`
63: @*/
64: PetscErrorCode VecMaxPointwiseDivide(Vec x, Vec y, PetscReal *max)
65: {
72: VecCheckSameSize(x, 1, y, 2);
73: VecLockReadPush(x);
74: VecLockReadPush(y);
75: PetscUseTypeMethod(x, maxpointwisedivide, y, max);
76: VecLockReadPop(x);
77: VecLockReadPop(y);
78: return 0;
79: }
81: /*@
82: VecDot - Computes the vector dot product.
84: Collective on x
86: Input Parameters:
87: . x, y - the vectors
89: Output Parameter:
90: . val - the dot product
92: Performance Issues:
93: .vb
94: per-processor memory bandwidth
95: interprocessor latency
96: work load imbalance that causes certain processes to arrive much earlier than others
97: .ve
99: Notes for Users of Complex Numbers:
100: For complex vectors, `VecDot()` computes
101: $ val = (x,y) = y^H x,
102: where y^H denotes the conjugate transpose of y. Note that this corresponds to the usual "mathematicians" complex
103: inner product where the SECOND argument gets the complex conjugate. Since the `BLASdot()` complex conjugates the first
104: first argument we call the `BLASdot()` with the arguments reversed.
106: Use `VecTDot()` for the indefinite form
107: $ val = (x,y) = y^T x,
108: where y^T denotes the transpose of y.
110: Level: intermediate
112: .seealso: `Vec`, `VecMDot()`, `VecTDot()`, `VecNorm()`, `VecDotBegin()`, `VecDotEnd()`, `VecDotRealPart()`
113: @*/
114: PetscErrorCode VecDot(Vec x, Vec y, PetscScalar *val)
115: {
122: VecCheckSameSize(x, 1, y, 2);
124: VecLockReadPush(x);
125: VecLockReadPush(y);
126: PetscLogEventBegin(VEC_Dot, x, y, 0, 0);
127: PetscUseTypeMethod(x, dot, y, val);
128: PetscLogEventEnd(VEC_Dot, x, y, 0, 0);
129: VecLockReadPop(x);
130: VecLockReadPop(y);
131: return 0;
132: }
134: /*@
135: VecDotRealPart - Computes the real part of the vector dot product.
137: Collective on x
139: Input Parameters:
140: . x, y - the vectors
142: Output Parameter:
143: . val - the real part of the dot product;
145: Level: intermediate
147: Performance Issues:
148: .vb
149: per-processor memory bandwidth
150: interprocessor latency
151: work load imbalance that causes certain processes to arrive much earlier than others
152: .ve
154: Notes for Users of Complex Numbers:
155: See `VecDot()` for more details on the definition of the dot product for complex numbers
157: For real numbers this returns the same value as `VecDot()`
159: 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
160: the space R^{2n} (that is a vector of 2n components with the real or imaginary part of the complex numbers for components)
162: Developer Note:
163: This is not currently optimized to compute only the real part of the dot product.
165: .seealso: `Vec`, `VecMDot()`, `VecTDot()`, `VecNorm()`, `VecDotBegin()`, `VecDotEnd()`, `VecDot()`, `VecDotNorm2()`
166: @*/
167: PetscErrorCode VecDotRealPart(Vec x, Vec y, PetscReal *val)
168: {
169: PetscScalar fdot;
171: VecDot(x, y, &fdot);
172: *val = PetscRealPart(fdot);
173: return 0;
174: }
176: /*@
177: VecNorm - Computes the vector norm.
179: Collective on Vec
181: Input Parameters:
182: + x - the vector
183: - type - the type of the norm requested
185: Output Parameter:
186: . val - the norm
188: Values of NormType:
189: + NORM_1 - sum_i |x_i|
190: . NORM_2 - sqrt(sum_i |x_i|^2)
191: . NORM_INFINITY - max_i |x_i|
192: - NORM_1_AND_2 - computes efficiently both NORM_1 and NORM_2 and stores them each in an output array
194: Notes:
195: For complex numbers NORM_1 will return the traditional 1 norm of the 2 norm of the complex numbers; that is the 1
196: norm of the absolute values of the complex entries. In PETSc 3.6 and earlier releases it returned the 1 norm of
197: the 1 norm of the complex entries (what is returned by the BLAS routine asum()). Both are valid norms but most
198: people expect the former.
200: This routine stashes the computed norm value, repeated calls before the vector entries are changed are then rapid since the
201: precomputed value is immediately available. Certain vector operations such as VecSet() store the norms so the value is
202: immediately available and does not need to be explicitly computed. VecScale() updates any stashed norm values, thus calls after VecScale()
203: do not need to explicitly recompute the norm.
205: Level: intermediate
207: Performance Issues:
208: + per-processor memory bandwidth - limits the speed of the computation of local portion of the norm
209: . interprocessor latency - limits the accumulation of the result across ranks, .i.e. MPI_Allreduce() time
210: . number of ranks - the time for the result will grow with the log base 2 of the number of ranks sharing the vector
211: - work load imbalance - the rank with the largest number of vector entries will limit the speed up
213: .seealso: `VecDot()`, `VecTDot()`, `VecNorm()`, `VecDotBegin()`, `VecDotEnd()`, `VecNormAvailable()`,
214: `VecNormBegin()`, `VecNormEnd()`, `NormType()`
216: @*/
217: PetscErrorCode VecNorm(Vec x, NormType type, PetscReal *val)
218: {
223: /* Cached data? */
224: if (type != NORM_1_AND_2) {
225: PetscBool flg;
227: PetscObjectComposedDataGetReal((PetscObject)x, NormIds[type], *val, flg);
228: if (flg) return 0;
229: }
231: VecLockReadPush(x);
232: PetscLogEventBegin(VEC_Norm, x, 0, 0, 0);
233: PetscUseTypeMethod(x, norm, type, val);
234: PetscLogEventEnd(VEC_Norm, x, 0, 0, 0);
235: VecLockReadPop(x);
237: if (type != NORM_1_AND_2) PetscObjectComposedDataSetReal((PetscObject)x, NormIds[type], *val);
238: return 0;
239: }
241: /*@
242: VecNormAvailable - Returns the vector norm if it is already known.
244: Not Collective
246: Input Parameters:
247: + x - the vector
248: - type - one of NORM_1, NORM_2, NORM_INFINITY. Also available
249: NORM_1_AND_2, which computes both norms and stores them
250: in a two element array.
252: Output Parameters:
253: + available - PETSC_TRUE if the val returned is valid
254: - val - the norm
256: Notes:
257: $ NORM_1 denotes sum_i |x_i|
258: $ NORM_2 denotes sqrt(sum_i (x_i)^2)
259: $ NORM_INFINITY denotes max_i |x_i|
261: Level: intermediate
263: Performance Issues:
264: $ per-processor memory bandwidth
265: $ interprocessor latency
266: $ work load imbalance that causes certain processes to arrive much earlier than others
268: Compile Option:
269: PETSC_HAVE_SLOW_BLAS_NORM2 will cause a C (loop unrolled) version of the norm to be used, rather
270: than the BLAS. This should probably only be used when one is using the FORTRAN BLAS routines
271: (as opposed to vendor provided) because the FORTRAN BLAS NRM2() routine is very slow.
273: .seealso: `VecDot()`, `VecTDot()`, `VecNorm()`, `VecDotBegin()`, `VecDotEnd()`, `VecNorm()`
274: `VecNormBegin()`, `VecNormEnd()`
276: @*/
277: PetscErrorCode VecNormAvailable(Vec x, NormType type, PetscBool *available, PetscReal *val)
278: {
284: if (type == NORM_1_AND_2) {
285: *available = PETSC_FALSE;
286: } else {
287: PetscObjectComposedDataGetReal((PetscObject)x, NormIds[type], *val, *available);
288: }
289: return 0;
290: }
292: /*@
293: VecNormalize - Normalizes a vector by 2-norm.
295: Collective on Vec
297: Input Parameter:
298: . x - the vector
300: Output Parameter:
301: . val - the vector norm before normalization. May be `NULL` if the value is not needed.
303: Level: intermediate
305: @*/
306: PetscErrorCode VecNormalize(Vec x, PetscReal *val)
307: {
308: PetscReal norm;
312: VecSetErrorIfLocked(x, 1);
314: PetscLogEventBegin(VEC_Normalize, x, 0, 0, 0);
315: VecNorm(x, NORM_2, &norm);
316: if (norm == 0.0) {
317: PetscInfo(x, "Vector of zero norm can not be normalized; Returning only the zero norm\n");
318: } else if (norm != 1.0) {
319: VecScale(x, 1.0 / norm);
320: }
321: PetscLogEventEnd(VEC_Normalize, x, 0, 0, 0);
322: if (val) *val = norm;
323: return 0;
324: }
326: /*@C
327: VecMax - Determines the vector component with maximum real part and its location.
329: Collective on Vec
331: Input Parameter:
332: . x - the vector
334: Output Parameters:
335: + p - the location of val (pass NULL if you don't want this)
336: - val - the maximum component
338: Notes:
339: Returns the value PETSC_MIN_REAL and negative p if the vector is of length 0.
341: Returns the smallest index with the maximum value
342: Level: intermediate
344: .seealso: `VecNorm()`, `VecMin()`
345: @*/
346: PetscErrorCode VecMax(Vec x, PetscInt *p, PetscReal *val)
347: {
352: VecLockReadPush(x);
353: PetscLogEventBegin(VEC_Max, x, 0, 0, 0);
354: PetscUseTypeMethod(x, max, p, val);
355: PetscLogEventEnd(VEC_Max, x, 0, 0, 0);
356: VecLockReadPop(x);
357: return 0;
358: }
360: /*@C
361: VecMin - Determines the vector component with minimum real part and its location.
363: Collective on Vec
365: Input Parameter:
366: . x - the vector
368: Output Parameters:
369: + p - the location of val (pass NULL if you don't want this location)
370: - val - the minimum component
372: Level: intermediate
374: Notes:
375: Returns the value PETSC_MAX_REAL and negative p if the vector is of length 0.
377: This returns the smallest index with the minumum value
379: .seealso: `VecMax()`
380: @*/
381: PetscErrorCode VecMin(Vec x, PetscInt *p, PetscReal *val)
382: {
387: VecLockReadPush(x);
388: PetscLogEventBegin(VEC_Min, x, 0, 0, 0);
389: PetscUseTypeMethod(x, min, p, val);
390: PetscLogEventEnd(VEC_Min, x, 0, 0, 0);
391: VecLockReadPop(x);
392: return 0;
393: }
395: /*@
396: VecTDot - Computes an indefinite vector dot product. That is, this
397: routine does NOT use the complex conjugate.
399: Collective on Vec
401: Input Parameters:
402: . x, y - the vectors
404: Output Parameter:
405: . val - the dot product
407: Notes for Users of Complex Numbers:
408: For complex vectors, VecTDot() computes the indefinite form
409: $ val = (x,y) = y^T x,
410: where y^T denotes the transpose of y.
412: Use VecDot() for the inner product
413: $ val = (x,y) = y^H x,
414: where y^H denotes the conjugate transpose of y.
416: Level: intermediate
418: .seealso: `VecDot()`, `VecMTDot()`
419: @*/
420: PetscErrorCode VecTDot(Vec x, Vec y, PetscScalar *val)
421: {
428: VecCheckSameSize(x, 1, y, 2);
430: VecLockReadPush(x);
431: VecLockReadPush(y);
432: PetscLogEventBegin(VEC_TDot, x, y, 0, 0);
433: PetscUseTypeMethod(x, tdot, y, val);
434: PetscLogEventEnd(VEC_TDot, x, y, 0, 0);
435: VecLockReadPop(x);
436: VecLockReadPop(y);
437: return 0;
438: }
440: /*@
441: VecScale - Scales a vector.
443: Not collective on Vec
445: Input Parameters:
446: + x - the vector
447: - alpha - the scalar
449: Note:
450: For a vector with n components, VecScale() computes
451: $ x[i] = alpha * x[i], for i=1,...,n.
453: Level: intermediate
455: @*/
456: PetscErrorCode VecScale(Vec x, PetscScalar alpha)
457: {
458: PetscReal norms[4];
459: PetscBool flgs[4];
464: VecSetErrorIfLocked(x, 1);
465: if (alpha == (PetscScalar)1.0) return 0;
467: /* get current stashed norms */
468: for (PetscInt i = 0; i < 4; i++) PetscObjectComposedDataGetReal((PetscObject)x, NormIds[i], norms[i], flgs[i]);
470: PetscLogEventBegin(VEC_Scale, x, 0, 0, 0);
471: PetscUseTypeMethod(x, scale, alpha);
472: PetscLogEventEnd(VEC_Scale, x, 0, 0, 0);
474: PetscObjectStateIncrease((PetscObject)x);
475: /* put the scaled stashed norms back into the Vec */
476: for (PetscInt i = 0; i < 4; i++) {
477: if (flgs[i]) PetscObjectComposedDataSetReal((PetscObject)x, NormIds[i], PetscAbsScalar(alpha) * norms[i]);
478: }
479: return 0;
480: }
482: /*@
483: VecSet - Sets all components of a vector to a single scalar value.
485: Logically Collective on x
487: Input Parameters:
488: + x - the vector
489: - alpha - the scalar
491: Output Parameter:
492: . x - the vector
494: Level: beginner
496: Note:
497: For a vector of dimension n, `VecSet()` computes
498: $ x[i] = alpha, for i=1,...,n,
499: so that all vector entries then equal the identical
500: scalar value, alpha. Use the more general routine
501: `VecSetValues()` to set different vector entries.
503: You CANNOT call this after you have called `VecSetValues()` but before you call
504: `VecAssemblyBegin()`
506: .seealso: `VecSetValues()`, `VecSetValuesBlocked()`, `VecSetRandom()`
507: @*/
508: PetscErrorCode VecSet(Vec x, PetscScalar alpha)
509: {
514: VecSetErrorIfLocked(x, 1);
516: PetscLogEventBegin(VEC_Set, x, 0, 0, 0);
517: PetscUseTypeMethod(x, set, alpha);
518: PetscLogEventEnd(VEC_Set, x, 0, 0, 0);
519: PetscObjectStateIncrease((PetscObject)x);
521: /* norms can be simply set (if |alpha|*N not too large) */
523: {
524: PetscReal val = PetscAbsScalar(alpha);
525: const PetscInt N = x->map->N;
527: if (N == 0) {
528: PetscObjectComposedDataSetReal((PetscObject)x, NormIds[NORM_1], 0.0l);
529: PetscObjectComposedDataSetReal((PetscObject)x, NormIds[NORM_INFINITY], 0.0);
530: PetscObjectComposedDataSetReal((PetscObject)x, NormIds[NORM_2], 0.0);
531: PetscObjectComposedDataSetReal((PetscObject)x, NormIds[NORM_FROBENIUS], 0.0);
532: } else if (val > PETSC_MAX_REAL / N) {
533: PetscObjectComposedDataSetReal((PetscObject)x, NormIds[NORM_INFINITY], val);
534: } else {
535: PetscObjectComposedDataSetReal((PetscObject)x, NormIds[NORM_1], N * val);
536: PetscObjectComposedDataSetReal((PetscObject)x, NormIds[NORM_INFINITY], val);
537: val = PetscSqrtReal((PetscReal)N) * val;
538: PetscObjectComposedDataSetReal((PetscObject)x, NormIds[NORM_2], val);
539: PetscObjectComposedDataSetReal((PetscObject)x, NormIds[NORM_FROBENIUS], val);
540: }
541: }
542: return 0;
543: }
545: /*@
546: VecAXPY - Computes y = alpha x + y.
548: Logically Collective on Vec
550: Input Parameters:
551: + alpha - the scalar
552: - x, y - the vectors
554: Output Parameter:
555: . y - output vector
557: Level: intermediate
559: Notes:
560: x and y MUST be different vectors
561: This routine is optimized for alpha of 0.0, otherwise it calls the BLAS routine
563: $ VecAXPY(y,alpha,x) y = alpha x + y
564: $ VecAYPX(y,beta,x) y = x + beta y
565: $ VecAXPBY(y,alpha,beta,x) y = alpha x + beta y
566: $ VecWAXPY(w,alpha,x,y) w = alpha x + y
567: $ VecAXPBYPCZ(w,alpha,beta,gamma,x,y) z = alpha x + beta y + gamma z
568: $ VecMAXPY(y,nv,alpha[],x[]) y = sum alpha[i] x[i] + y
570: .seealso: `VecAYPX()`, `VecMAXPY()`, `VecWAXPY()`, `VecAXPBYPCZ()`, `VecAXPBY()`
571: @*/
572: PetscErrorCode VecAXPY(Vec y, PetscScalar alpha, Vec x)
573: {
579: VecCheckSameSize(x, 3, y, 1);
582: if (alpha == (PetscScalar)0.0) return 0;
584: VecSetErrorIfLocked(y, 1);
585: VecLockReadPush(x);
586: PetscLogEventBegin(VEC_AXPY, x, y, 0, 0);
587: PetscUseTypeMethod(y, axpy, alpha, x);
588: PetscLogEventEnd(VEC_AXPY, x, y, 0, 0);
589: VecLockReadPop(x);
590: PetscObjectStateIncrease((PetscObject)y);
591: return 0;
592: }
594: /*@
595: VecAYPX - Computes y = x + beta y.
597: Logically Collective on Vec
599: Input Parameters:
600: + beta - the scalar
601: - x, y - the vectors
603: Output Parameter:
604: . y - output vector
606: Level: intermediate
608: Notes:
609: x and y MUST be different vectors
610: The implementation is optimized for beta of -1.0, 0.0, and 1.0
612: .seealso: `VecMAXPY()`, `VecWAXPY()`, `VecAXPY()`, `VecAXPBYPCZ()`, `VecAXPBY()`
613: @*/
614: PetscErrorCode VecAYPX(Vec y, PetscScalar beta, Vec x)
615: {
621: VecCheckSameSize(x, 1, y, 3);
624: VecSetErrorIfLocked(y, 1);
625: VecLockReadPush(x);
626: if (beta == (PetscScalar)0.0) {
627: VecCopy(x, y);
628: } else {
629: PetscLogEventBegin(VEC_AYPX, x, y, 0, 0);
630: PetscUseTypeMethod(y, aypx, beta, x);
631: PetscLogEventEnd(VEC_AYPX, x, y, 0, 0);
632: PetscObjectStateIncrease((PetscObject)y);
633: }
634: VecLockReadPop(x);
635: return 0;
636: }
638: /*@
639: VecAXPBY - Computes y = alpha x + beta y.
641: Logically Collective on Vec
643: Input Parameters:
644: + alpha,beta - the scalars
645: - x, y - the vectors
647: Output Parameter:
648: . y - output vector
650: Level: intermediate
652: Notes:
653: x and y MUST be different vectors
654: The implementation is optimized for alpha and/or beta values of 0.0 and 1.0
656: .seealso: `VecAYPX()`, `VecMAXPY()`, `VecWAXPY()`, `VecAXPY()`, `VecAXPBYPCZ()`
657: @*/
658: PetscErrorCode VecAXPBY(Vec y, PetscScalar alpha, PetscScalar beta, Vec x)
659: {
665: VecCheckSameSize(y, 1, x, 4);
669: if (alpha == (PetscScalar)0.0 && beta == (PetscScalar)1.0) return 0;
671: VecSetErrorIfLocked(y, 1);
672: VecLockReadPush(x);
673: PetscLogEventBegin(VEC_AXPY, y, x, 0, 0);
674: PetscUseTypeMethod(y, axpby, alpha, beta, x);
675: PetscLogEventEnd(VEC_AXPY, y, x, 0, 0);
676: PetscObjectStateIncrease((PetscObject)y);
677: VecLockReadPop(x);
678: return 0;
679: }
681: /*@
682: VecAXPBYPCZ - Computes z = alpha x + beta y + gamma z
684: Logically Collective on Vec
686: Input Parameters:
687: + alpha,beta, gamma - the scalars
688: - x, y, z - the vectors
690: Output Parameter:
691: . z - output vector
693: Level: intermediate
695: Notes:
696: x, y and z must be different vectors
697: The implementation is optimized for alpha of 1.0 and gamma of 1.0 or 0.0
699: .seealso: `VecAYPX()`, `VecMAXPY()`, `VecWAXPY()`, `VecAXPY()`, `VecAXPBY()`
700: @*/
701: PetscErrorCode VecAXPBYPCZ(Vec z, PetscScalar alpha, PetscScalar beta, PetscScalar gamma, Vec x, Vec y)
702: {
711: VecCheckSameSize(x, 1, y, 5);
712: VecCheckSameSize(x, 1, z, 6);
718: if (alpha == (PetscScalar)0.0 && beta == (PetscScalar)0.0 && gamma == (PetscScalar)1.0) return 0;
720: VecSetErrorIfLocked(z, 1);
721: VecLockReadPush(x);
722: VecLockReadPush(y);
723: PetscLogEventBegin(VEC_AXPBYPCZ, x, y, z, 0);
724: PetscUseTypeMethod(z, axpbypcz, alpha, beta, gamma, x, y);
725: PetscLogEventEnd(VEC_AXPBYPCZ, x, y, z, 0);
726: PetscObjectStateIncrease((PetscObject)z);
727: VecLockReadPop(x);
728: VecLockReadPop(y);
729: return 0;
730: }
732: /*@
733: VecWAXPY - Computes w = alpha x + y.
735: Logically Collective on Vec
737: Input Parameters:
738: + alpha - the scalar
739: - x, y - the vectors
741: Output Parameter:
742: . w - the result
744: Level: intermediate
746: Notes:
747: w cannot be either x or y, but x and y can be the same
748: The implementation is optimzed for alpha of -1.0, 0.0, and 1.0
750: .seealso: `VecAXPY()`, `VecAYPX()`, `VecAXPBY()`, `VecMAXPY()`, `VecAXPBYPCZ()`
751: @*/
752: PetscErrorCode VecWAXPY(Vec w, PetscScalar alpha, Vec x, Vec y)
753: {
762: VecCheckSameSize(x, 3, y, 4);
763: VecCheckSameSize(x, 3, w, 1);
767: VecSetErrorIfLocked(w, 1);
769: VecLockReadPush(x);
770: VecLockReadPush(y);
771: if (alpha == (PetscScalar)0.0) {
772: VecCopy(y, w);
773: } else {
774: PetscLogEventBegin(VEC_WAXPY, x, y, w, 0);
775: PetscUseTypeMethod(w, waxpy, alpha, x, y);
776: PetscLogEventEnd(VEC_WAXPY, x, y, w, 0);
777: PetscObjectStateIncrease((PetscObject)w);
778: }
779: VecLockReadPop(x);
780: VecLockReadPop(y);
781: return 0;
782: }
784: /*@C
785: VecSetValues - Inserts or adds values into certain locations of a vector.
787: Not Collective
789: Input Parameters:
790: + x - vector to insert in
791: . ni - number of elements to add
792: . ix - indices where to add
793: . y - array of values
794: - iora - either INSERT_VALUES or ADD_VALUES, where
795: ADD_VALUES adds values to any existing entries, and
796: INSERT_VALUES replaces existing entries with new values
798: Notes:
799: VecSetValues() sets x[ix[i]] = y[i], for i=0,...,ni-1.
801: Calls to VecSetValues() with the INSERT_VALUES and ADD_VALUES
802: options cannot be mixed without intervening calls to the assembly
803: routines.
805: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
806: MUST be called after all calls to VecSetValues() have been completed.
808: VecSetValues() uses 0-based indices in Fortran as well as in C.
810: If you call VecSetOption(x, VEC_IGNORE_NEGATIVE_INDICES,PETSC_TRUE),
811: negative indices may be passed in ix. These rows are
812: simply ignored. This allows easily inserting element load matrices
813: with homogeneous Dirchlet boundary conditions that you don't want represented
814: in the vector.
816: Level: beginner
818: .seealso: `VecAssemblyBegin()`, `VecAssemblyEnd()`, `VecSetValuesLocal()`,
819: `VecSetValue()`, `VecSetValuesBlocked()`, `InsertMode`, `INSERT_VALUES`, `ADD_VALUES`, `VecGetValues()`
820: @*/
821: PetscErrorCode VecSetValues(Vec x, PetscInt ni, const PetscInt ix[], const PetscScalar y[], InsertMode iora)
822: {
825: if (!ni) return 0;
830: PetscLogEventBegin(VEC_SetValues, x, 0, 0, 0);
831: PetscUseTypeMethod(x, setvalues, ni, ix, y, iora);
832: PetscLogEventEnd(VEC_SetValues, x, 0, 0, 0);
833: PetscObjectStateIncrease((PetscObject)x);
834: return 0;
835: }
837: /*@C
838: VecGetValues - Gets values from certain locations of a vector. Currently
839: can only get values on the same processor
841: Not Collective
843: Input Parameters:
844: + x - vector to get values from
845: . ni - number of elements to get
846: - ix - indices where to get them from (in global 1d numbering)
848: Output Parameter:
849: . y - array of values
851: Notes:
852: The user provides the allocated array y; it is NOT allocated in this routine
854: VecGetValues() gets y[i] = x[ix[i]], for i=0,...,ni-1.
856: VecAssemblyBegin() and VecAssemblyEnd() MUST be called before calling this
858: VecGetValues() uses 0-based indices in Fortran as well as in C.
860: If you call VecSetOption(x, VEC_IGNORE_NEGATIVE_INDICES,PETSC_TRUE),
861: negative indices may be passed in ix. These rows are
862: simply ignored.
864: Level: beginner
866: .seealso: `VecAssemblyBegin()`, `VecAssemblyEnd()`, `VecSetValues()`
867: @*/
868: PetscErrorCode VecGetValues(Vec x, PetscInt ni, const PetscInt ix[], PetscScalar y[])
869: {
871: if (!ni) return 0;
875: PetscUseTypeMethod(x, getvalues, ni, ix, y);
876: return 0;
877: }
879: /*@C
880: VecSetValuesBlocked - Inserts or adds blocks of values into certain locations of a vector.
882: Not Collective
884: Input Parameters:
885: + x - vector to insert in
886: . ni - number of blocks to add
887: . ix - indices where to add in block count, rather than element count
888: . y - array of values
889: - iora - either INSERT_VALUES or ADD_VALUES, where
890: ADD_VALUES adds values to any existing entries, and
891: INSERT_VALUES replaces existing entries with new values
893: Notes:
894: VecSetValuesBlocked() sets x[bs*ix[i]+j] = y[bs*i+j],
895: for j=0,...,bs-1, for i=0,...,ni-1. where bs was set with VecSetBlockSize().
897: Calls to VecSetValuesBlocked() with the INSERT_VALUES and ADD_VALUES
898: options cannot be mixed without intervening calls to the assembly
899: routines.
901: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
902: MUST be called after all calls to VecSetValuesBlocked() have been completed.
904: VecSetValuesBlocked() uses 0-based indices in Fortran as well as in C.
906: Negative indices may be passed in ix, these rows are
907: simply ignored. This allows easily inserting element load matrices
908: with homogeneous Dirchlet boundary conditions that you don't want represented
909: in the vector.
911: Level: intermediate
913: .seealso: `VecAssemblyBegin()`, `VecAssemblyEnd()`, `VecSetValuesBlockedLocal()`,
914: `VecSetValues()`
915: @*/
916: PetscErrorCode VecSetValuesBlocked(Vec x, PetscInt ni, const PetscInt ix[], const PetscScalar y[], InsertMode iora)
917: {
920: if (!ni) return 0;
925: PetscLogEventBegin(VEC_SetValues, x, 0, 0, 0);
926: PetscUseTypeMethod(x, setvaluesblocked, ni, ix, y, iora);
927: PetscLogEventEnd(VEC_SetValues, x, 0, 0, 0);
928: PetscObjectStateIncrease((PetscObject)x);
929: return 0;
930: }
932: /*@C
933: VecSetValuesLocal - Inserts or adds values into certain locations of a vector,
934: using a local ordering of the nodes.
936: Not Collective
938: Input Parameters:
939: + x - vector to insert in
940: . ni - number of elements to add
941: . ix - indices where to add
942: . y - array of values
943: - iora - either INSERT_VALUES or ADD_VALUES, where
944: ADD_VALUES adds values to any existing entries, and
945: INSERT_VALUES replaces existing entries with new values
947: Level: intermediate
949: Notes:
950: VecSetValuesLocal() sets x[ix[i]] = y[i], for i=0,...,ni-1.
952: Calls to VecSetValues() with the INSERT_VALUES and ADD_VALUES
953: options cannot be mixed without intervening calls to the assembly
954: routines.
956: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
957: MUST be called after all calls to VecSetValuesLocal() have been completed.
959: VecSetValuesLocal() uses 0-based indices in Fortran as well as in C.
961: .seealso: `VecAssemblyBegin()`, `VecAssemblyEnd()`, `VecSetValues()`, `VecSetLocalToGlobalMapping()`,
962: `VecSetValuesBlockedLocal()`
963: @*/
964: PetscErrorCode VecSetValuesLocal(Vec x, PetscInt ni, const PetscInt ix[], const PetscScalar y[], InsertMode iora)
965: {
966: PetscInt lixp[128], *lix = lixp;
970: if (!ni) return 0;
975: PetscLogEventBegin(VEC_SetValues, x, 0, 0, 0);
976: if (!x->ops->setvalueslocal) {
977: if (x->map->mapping) {
978: if (ni > 128) PetscMalloc1(ni, &lix);
979: ISLocalToGlobalMappingApply(x->map->mapping, ni, (PetscInt *)ix, lix);
980: PetscUseTypeMethod(x, setvalues, ni, lix, y, iora);
981: if (ni > 128) PetscFree(lix);
982: } else PetscUseTypeMethod(x, setvalues, ni, ix, y, iora);
983: } else PetscUseTypeMethod(x, setvalueslocal, ni, ix, y, iora);
984: PetscLogEventEnd(VEC_SetValues, x, 0, 0, 0);
985: PetscObjectStateIncrease((PetscObject)x);
986: return 0;
987: }
989: /*@
990: VecSetValuesBlockedLocal - Inserts or adds values into certain locations of a vector,
991: using a local ordering of the nodes.
993: Not Collective
995: Input Parameters:
996: + x - vector to insert in
997: . ni - number of blocks to add
998: . ix - indices where to add in block count, not element count
999: . y - array of values
1000: - iora - either INSERT_VALUES or ADD_VALUES, where
1001: ADD_VALUES adds values to any existing entries, and
1002: INSERT_VALUES replaces existing entries with new values
1004: Level: intermediate
1006: Notes:
1007: VecSetValuesBlockedLocal() sets x[bs*ix[i]+j] = y[bs*i+j],
1008: for j=0,..bs-1, for i=0,...,ni-1, where bs has been set with VecSetBlockSize().
1010: Calls to VecSetValuesBlockedLocal() with the INSERT_VALUES and ADD_VALUES
1011: options cannot be mixed without intervening calls to the assembly
1012: routines.
1014: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
1015: MUST be called after all calls to VecSetValuesBlockedLocal() have been completed.
1017: VecSetValuesBlockedLocal() uses 0-based indices in Fortran as well as in C.
1019: .seealso: `VecAssemblyBegin()`, `VecAssemblyEnd()`, `VecSetValues()`, `VecSetValuesBlocked()`,
1020: `VecSetLocalToGlobalMapping()`
1021: @*/
1022: PetscErrorCode VecSetValuesBlockedLocal(Vec x, PetscInt ni, const PetscInt ix[], const PetscScalar y[], InsertMode iora)
1023: {
1024: PetscInt lixp[128], *lix = lixp;
1028: if (!ni) return 0;
1032: PetscLogEventBegin(VEC_SetValues, x, 0, 0, 0);
1033: if (x->map->mapping) {
1034: if (ni > 128) PetscMalloc1(ni, &lix);
1035: ISLocalToGlobalMappingApplyBlock(x->map->mapping, ni, (PetscInt *)ix, lix);
1036: PetscUseTypeMethod(x, setvaluesblocked, ni, lix, y, iora);
1037: if (ni > 128) PetscFree(lix);
1038: } else {
1039: PetscUseTypeMethod(x, setvaluesblocked, ni, ix, y, iora);
1040: }
1041: PetscLogEventEnd(VEC_SetValues, x, 0, 0, 0);
1042: PetscObjectStateIncrease((PetscObject)x);
1043: return 0;
1044: }
1046: /*@
1047: VecMTDot - Computes indefinite vector multiple dot products.
1048: That is, it does NOT use the complex conjugate.
1050: Collective on Vec
1052: Input Parameters:
1053: + x - one vector
1054: . nv - number of vectors
1055: - y - array of vectors. Note that vectors are pointers
1057: Output Parameter:
1058: . val - array of the dot products
1060: Notes for Users of Complex Numbers:
1061: For complex vectors, VecMTDot() computes the indefinite form
1062: $ val = (x,y) = y^T x,
1063: where y^T denotes the transpose of y.
1065: Use VecMDot() for the inner product
1066: $ val = (x,y) = y^H x,
1067: where y^H denotes the conjugate transpose of y.
1069: Level: intermediate
1071: .seealso: `VecMDot()`, `VecTDot()`
1072: @*/
1073: PetscErrorCode VecMTDot(Vec x, PetscInt nv, const Vec y[], PetscScalar val[])
1074: {
1078: if (!nv) return 0;
1080: for (PetscInt i = 0; i < nv; ++i) {
1084: VecCheckSameSize(x, 1, y[i], 3);
1085: VecLockReadPush(y[i]);
1086: }
1089: VecLockReadPush(x);
1090: PetscLogEventBegin(VEC_MTDot, x, *y, 0, 0);
1091: PetscUseTypeMethod(x, mtdot, nv, y, val);
1092: PetscLogEventEnd(VEC_MTDot, x, *y, 0, 0);
1093: VecLockReadPop(x);
1094: for (PetscInt i = 0; i < nv; ++i) VecLockReadPop(y[i]);
1095: return 0;
1096: }
1098: /*@
1099: VecMDot - Computes vector multiple dot products.
1101: Collective on Vec
1103: Input Parameters:
1104: + x - one vector
1105: . nv - number of vectors
1106: - y - array of vectors.
1108: Output Parameter:
1109: . val - array of the dot products (does not allocate the array)
1111: Notes for Users of Complex Numbers:
1112: For complex vectors, VecMDot() computes
1113: $ val = (x,y) = y^H x,
1114: where y^H denotes the conjugate transpose of y.
1116: Use VecMTDot() for the indefinite form
1117: $ val = (x,y) = y^T x,
1118: where y^T denotes the transpose of y.
1120: Level: intermediate
1122: .seealso: `VecMTDot()`, `VecDot()`
1123: @*/
1124: PetscErrorCode VecMDot(Vec x, PetscInt nv, const Vec y[], PetscScalar val[])
1125: {
1129: if (!nv) return 0;
1131: for (PetscInt i = 0; i < nv; ++i) {
1135: VecCheckSameSize(x, 1, y[i], 3);
1136: VecLockReadPush(y[i]);
1137: }
1140: VecLockReadPush(x);
1141: PetscLogEventBegin(VEC_MDot, x, *y, 0, 0);
1142: PetscUseTypeMethod(x, mdot, nv, y, val);
1143: PetscLogEventEnd(VEC_MDot, x, *y, 0, 0);
1144: VecLockReadPop(x);
1145: for (PetscInt i = 0; i < nv; ++i) VecLockReadPop(y[i]);
1146: return 0;
1147: }
1149: /*@
1150: VecMAXPY - Computes y = y + sum alpha[i] x[i]
1152: Logically Collective on Vec
1154: Input Parameters:
1155: + nv - number of scalars and x-vectors
1156: . alpha - array of scalars
1157: . y - one vector
1158: - x - array of vectors
1160: Level: intermediate
1162: Notes:
1163: y cannot be any of the x vectors
1165: .seealso: `VecAYPX()`, `VecWAXPY()`, `VecAXPY()`, `VecAXPBYPCZ()`, `VecAXPBY()`
1166: @*/
1167: PetscErrorCode VecMAXPY(Vec y, PetscInt nv, const PetscScalar alpha[], Vec x[])
1168: {
1171: VecSetErrorIfLocked(y, 1);
1173: if (nv) {
1174: PetscInt zeros = 0;
1178: for (PetscInt i = 0; i < nv; ++i) {
1183: VecCheckSameSize(y, 1, x[i], 4);
1185: VecLockReadPush(x[i]);
1186: zeros += alpha[i] == (PetscScalar)0.0;
1187: }
1189: if (zeros < nv) {
1190: PetscLogEventBegin(VEC_MAXPY, y, *x, 0, 0);
1191: PetscUseTypeMethod(y, maxpy, nv, alpha, x);
1192: PetscLogEventEnd(VEC_MAXPY, y, *x, 0, 0);
1193: PetscObjectStateIncrease((PetscObject)y);
1194: }
1196: for (PetscInt i = 0; i < nv; ++i) VecLockReadPop(x[i]);
1197: }
1198: return 0;
1199: }
1201: /*@
1202: VecConcatenate - Creates a new vector that is a vertical concatenation of all the given array of vectors
1203: in the order they appear in the array. The concatenated vector resides on the same
1204: communicator and is the same type as the source vectors.
1206: Collective on X
1208: Input Parameters:
1209: + nx - number of vectors to be concatenated
1210: - X - array containing the vectors to be concatenated in the order of concatenation
1212: Output Parameters:
1213: + Y - concatenated vector
1214: - x_is - array of index sets corresponding to the concatenated components of Y (NULL if not needed)
1216: Notes:
1217: Concatenation is similar to the functionality of a VecNest object; they both represent combination of
1218: different vector spaces. However, concatenated vectors do not store any information about their
1219: sub-vectors and own their own data. Consequently, this function provides index sets to enable the
1220: manipulation of data in the concatenated vector that corresponds to the original components at creation.
1222: This is a useful tool for outer loop algorithms, particularly constrained optimizers, where the solver
1223: has to operate on combined vector spaces and cannot utilize VecNest objects due to incompatibility with
1224: bound projections.
1226: Level: advanced
1228: .seealso: `VECNEST`, `VECSCATTER`, `VecScatterCreate()`
1229: @*/
1230: PetscErrorCode VecConcatenate(PetscInt nx, const Vec X[], Vec *Y, IS *x_is[])
1231: {
1232: MPI_Comm comm;
1233: VecType vec_type;
1234: Vec Ytmp, Xtmp;
1235: IS *is_tmp;
1236: PetscInt i, shift = 0, Xnl, Xng, Xbegin;
1243: if ((*X)->ops->concatenate) {
1244: /* use the dedicated concatenation function if available */
1245: (*(*X)->ops->concatenate)(nx, X, Y, x_is);
1246: } else {
1247: /* loop over vectors and start creating IS */
1248: comm = PetscObjectComm((PetscObject)(*X));
1249: VecGetType(*X, &vec_type);
1250: PetscMalloc1(nx, &is_tmp);
1251: for (i = 0; i < nx; i++) {
1252: VecGetSize(X[i], &Xng);
1253: VecGetLocalSize(X[i], &Xnl);
1254: VecGetOwnershipRange(X[i], &Xbegin, NULL);
1255: ISCreateStride(comm, Xnl, shift + Xbegin, 1, &is_tmp[i]);
1256: shift += Xng;
1257: }
1258: /* create the concatenated vector */
1259: VecCreate(comm, &Ytmp);
1260: VecSetType(Ytmp, vec_type);
1261: VecSetSizes(Ytmp, PETSC_DECIDE, shift);
1262: VecSetUp(Ytmp);
1263: /* copy data from X array to Y and return */
1264: for (i = 0; i < nx; i++) {
1265: VecGetSubVector(Ytmp, is_tmp[i], &Xtmp);
1266: VecCopy(X[i], Xtmp);
1267: VecRestoreSubVector(Ytmp, is_tmp[i], &Xtmp);
1268: }
1269: *Y = Ytmp;
1270: if (x_is) {
1271: *x_is = is_tmp;
1272: } else {
1273: for (i = 0; i < nx; i++) ISDestroy(&is_tmp[i]);
1274: PetscFree(is_tmp);
1275: }
1276: }
1277: return 0;
1278: }
1280: /* A helper function for VecGetSubVector to check if we can implement it with no-copy (i.e. the subvector shares
1281: memory with the original vector), and the block size of the subvector.
1283: Input Parameters:
1284: + X - the original vector
1285: - is - the index set of the subvector
1287: Output Parameters:
1288: + contig - PETSC_TRUE if the index set refers to contiguous entries on this process, else PETSC_FALSE
1289: . start - start of contiguous block, as an offset from the start of the ownership range of the original vector
1290: - blocksize - the block size of the subvector
1292: */
1293: PetscErrorCode VecGetSubVectorContiguityAndBS_Private(Vec X, IS is, PetscBool *contig, PetscInt *start, PetscInt *blocksize)
1294: {
1295: PetscInt gstart, gend, lstart;
1296: PetscBool red[2] = {PETSC_TRUE /*contiguous*/, PETSC_TRUE /*validVBS*/};
1297: PetscInt n, N, ibs, vbs, bs = -1;
1299: ISGetLocalSize(is, &n);
1300: ISGetSize(is, &N);
1301: ISGetBlockSize(is, &ibs);
1302: VecGetBlockSize(X, &vbs);
1303: VecGetOwnershipRange(X, &gstart, &gend);
1304: ISContiguousLocal(is, gstart, gend, &lstart, &red[0]);
1305: /* block size is given by IS if ibs > 1; otherwise, check the vector */
1306: if (ibs > 1) {
1307: MPIU_Allreduce(MPI_IN_PLACE, red, 1, MPIU_BOOL, MPI_LAND, PetscObjectComm((PetscObject)is));
1308: bs = ibs;
1309: } else {
1310: if (n % vbs || vbs == 1) red[1] = PETSC_FALSE; /* this process invalidate the collectiveness of block size */
1311: MPIU_Allreduce(MPI_IN_PLACE, red, 2, MPIU_BOOL, MPI_LAND, PetscObjectComm((PetscObject)is));
1312: if (red[0] && red[1]) bs = vbs; /* all processes have a valid block size and the access will be contiguous */
1313: }
1315: *contig = red[0];
1316: *start = lstart;
1317: *blocksize = bs;
1318: return 0;
1319: }
1321: /* A helper function for VecGetSubVector, to be used when we have to build a standalone subvector through VecScatter
1323: Input Parameters:
1324: + X - the original vector
1325: . is - the index set of the subvector
1326: - bs - the block size of the subvector, gotten from VecGetSubVectorContiguityAndBS_Private()
1328: Output Parameters:
1329: . Z - the subvector, which will compose the VecScatter context on output
1330: */
1331: PetscErrorCode VecGetSubVectorThroughVecScatter_Private(Vec X, IS is, PetscInt bs, Vec *Z)
1332: {
1333: PetscInt n, N;
1334: VecScatter vscat;
1335: Vec Y;
1337: ISGetLocalSize(is, &n);
1338: ISGetSize(is, &N);
1339: VecCreate(PetscObjectComm((PetscObject)is), &Y);
1340: VecSetSizes(Y, n, N);
1341: VecSetBlockSize(Y, bs);
1342: VecSetType(Y, ((PetscObject)X)->type_name);
1343: VecScatterCreate(X, is, Y, NULL, &vscat);
1344: VecScatterBegin(vscat, X, Y, INSERT_VALUES, SCATTER_FORWARD);
1345: VecScatterEnd(vscat, X, Y, INSERT_VALUES, SCATTER_FORWARD);
1346: PetscObjectCompose((PetscObject)Y, "VecGetSubVector_Scatter", (PetscObject)vscat);
1347: VecScatterDestroy(&vscat);
1348: *Z = Y;
1349: return 0;
1350: }
1352: /*@
1353: VecGetSubVector - Gets a vector representing part of another vector
1355: Collective on X and IS
1357: Input Parameters:
1358: + X - vector from which to extract a subvector
1359: - is - index set representing portion of X to extract
1361: Output Parameter:
1362: . Y - subvector corresponding to is
1364: Level: advanced
1366: Notes:
1367: The subvector Y should be returned with VecRestoreSubVector().
1368: X and is must be defined on the same communicator
1370: This function may return a subvector without making a copy, therefore it is not safe to use the original vector while
1371: modifying the subvector. Other non-overlapping subvectors can still be obtained from X using this function.
1372: 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.
1374: .seealso: `MatCreateSubMatrix()`
1375: @*/
1376: PetscErrorCode VecGetSubVector(Vec X, IS is, Vec *Y)
1377: {
1378: Vec Z;
1384: if (X->ops->getsubvector) {
1385: PetscUseTypeMethod(X, getsubvector, is, &Z);
1386: } else { /* Default implementation currently does no caching */
1387: PetscBool contig;
1388: PetscInt n, N, start, bs;
1390: ISGetLocalSize(is, &n);
1391: ISGetSize(is, &N);
1392: VecGetSubVectorContiguityAndBS_Private(X, is, &contig, &start, &bs);
1393: if (contig) { /* We can do a no-copy implementation */
1394: const PetscScalar *x;
1395: PetscInt state = 0;
1396: PetscBool isstd, iscuda, iship;
1398: PetscObjectTypeCompareAny((PetscObject)X, &isstd, VECSEQ, VECMPI, VECSTANDARD, "");
1399: PetscObjectTypeCompareAny((PetscObject)X, &iscuda, VECSEQCUDA, VECMPICUDA, "");
1400: PetscObjectTypeCompareAny((PetscObject)X, &iship, VECSEQHIP, VECMPIHIP, "");
1401: if (iscuda) {
1402: #if defined(PETSC_HAVE_CUDA)
1403: const PetscScalar *x_d;
1404: PetscMPIInt size;
1405: PetscOffloadMask flg;
1407: VecCUDAGetArrays_Private(X, &x, &x_d, &flg);
1410: if (x) x += start;
1411: if (x_d) x_d += start;
1412: MPI_Comm_size(PetscObjectComm((PetscObject)X), &size);
1413: if (size == 1) {
1414: VecCreateSeqCUDAWithArrays(PetscObjectComm((PetscObject)X), bs, n, x, x_d, &Z);
1415: } else {
1416: VecCreateMPICUDAWithArrays(PetscObjectComm((PetscObject)X), bs, n, N, x, x_d, &Z);
1417: }
1418: Z->offloadmask = flg;
1419: #endif
1420: } else if (iship) {
1421: #if defined(PETSC_HAVE_HIP)
1422: const PetscScalar *x_d;
1423: PetscMPIInt size;
1424: PetscOffloadMask flg;
1426: VecHIPGetArrays_Private(X, &x, &x_d, &flg);
1429: if (x) x += start;
1430: if (x_d) x_d += start;
1431: MPI_Comm_size(PetscObjectComm((PetscObject)X), &size);
1432: if (size == 1) {
1433: VecCreateSeqHIPWithArrays(PetscObjectComm((PetscObject)X), bs, n, x, x_d, &Z);
1434: } else {
1435: VecCreateMPIHIPWithArrays(PetscObjectComm((PetscObject)X), bs, n, N, x, x_d, &Z);
1436: }
1437: Z->offloadmask = flg;
1438: #endif
1439: } else if (isstd) {
1440: PetscMPIInt size;
1442: MPI_Comm_size(PetscObjectComm((PetscObject)X), &size);
1443: VecGetArrayRead(X, &x);
1444: if (x) x += start;
1445: if (size == 1) {
1446: VecCreateSeqWithArray(PetscObjectComm((PetscObject)X), bs, n, x, &Z);
1447: } else {
1448: VecCreateMPIWithArray(PetscObjectComm((PetscObject)X), bs, n, N, x, &Z);
1449: }
1450: VecRestoreArrayRead(X, &x);
1451: } else { /* default implementation: use place array */
1452: VecGetArrayRead(X, &x);
1453: VecCreate(PetscObjectComm((PetscObject)X), &Z);
1454: VecSetType(Z, ((PetscObject)X)->type_name);
1455: VecSetSizes(Z, n, N);
1456: VecSetBlockSize(Z, bs);
1457: VecPlaceArray(Z, x ? x + start : NULL);
1458: VecRestoreArrayRead(X, &x);
1459: }
1461: /* this is relevant only in debug mode */
1462: VecLockGet(X, &state);
1463: if (state) VecLockReadPush(Z);
1464: Z->ops->placearray = NULL;
1465: Z->ops->replacearray = NULL;
1466: } else { /* Have to create a scatter and do a copy */
1467: VecGetSubVectorThroughVecScatter_Private(X, is, bs, &Z);
1468: }
1469: }
1470: /* Record the state when the subvector was gotten so we know whether its values need to be put back */
1471: if (VecGetSubVectorSavedStateId < 0) PetscObjectComposedDataRegister(&VecGetSubVectorSavedStateId);
1472: PetscObjectComposedDataSetInt((PetscObject)Z, VecGetSubVectorSavedStateId, 1);
1473: *Y = Z;
1474: return 0;
1475: }
1477: /*@
1478: VecRestoreSubVector - Restores a subvector extracted using VecGetSubVector()
1480: Collective on IS
1482: Input Parameters:
1483: + X - vector from which subvector was obtained
1484: . is - index set representing the subset of X
1485: - Y - subvector being restored
1487: Level: advanced
1489: .seealso: `VecGetSubVector()`
1490: @*/
1491: PetscErrorCode VecRestoreSubVector(Vec X, IS is, Vec *Y)
1492: {
1493: PETSC_UNUSED PetscObjectState dummystate = 0;
1494: PetscBool unchanged;
1502: if (X->ops->restoresubvector) PetscUseTypeMethod(X, restoresubvector, is, Y);
1503: else {
1504: PetscObjectComposedDataGetInt((PetscObject)*Y, VecGetSubVectorSavedStateId, dummystate, unchanged);
1505: if (!unchanged) { /* If Y's state has not changed since VecGetSubVector(), we only need to destroy Y */
1506: VecScatter scatter;
1507: PetscInt state;
1509: VecLockGet(X, &state);
1512: PetscObjectQuery((PetscObject)*Y, "VecGetSubVector_Scatter", (PetscObject *)&scatter);
1513: if (scatter) {
1514: VecScatterBegin(scatter, *Y, X, INSERT_VALUES, SCATTER_REVERSE);
1515: VecScatterEnd(scatter, *Y, X, INSERT_VALUES, SCATTER_REVERSE);
1516: } else {
1517: PetscBool iscuda, iship;
1518: PetscObjectTypeCompareAny((PetscObject)X, &iscuda, VECSEQCUDA, VECMPICUDA, "");
1519: PetscObjectTypeCompareAny((PetscObject)X, &iship, VECSEQHIP, VECMPIHIP, "");
1521: if (iscuda) {
1522: #if defined(PETSC_HAVE_CUDA)
1523: PetscOffloadMask ymask = (*Y)->offloadmask;
1525: /* The offloadmask of X dictates where to move memory
1526: If X GPU data is valid, then move Y data on GPU if needed
1527: Otherwise, move back to the CPU */
1528: switch (X->offloadmask) {
1529: case PETSC_OFFLOAD_BOTH:
1530: if (ymask == PETSC_OFFLOAD_CPU) {
1531: VecCUDAResetArray(*Y);
1532: } else if (ymask == PETSC_OFFLOAD_GPU) {
1533: X->offloadmask = PETSC_OFFLOAD_GPU;
1534: }
1535: break;
1536: case PETSC_OFFLOAD_GPU:
1537: if (ymask == PETSC_OFFLOAD_CPU) VecCUDAResetArray(*Y);
1538: break;
1539: case PETSC_OFFLOAD_CPU:
1540: if (ymask == PETSC_OFFLOAD_GPU) VecResetArray(*Y);
1541: break;
1542: case PETSC_OFFLOAD_UNALLOCATED:
1543: case PETSC_OFFLOAD_KOKKOS:
1544: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "This should not happen");
1545: }
1546: #endif
1547: } else if (iship) {
1548: #if defined(PETSC_HAVE_HIP)
1549: PetscOffloadMask ymask = (*Y)->offloadmask;
1551: /* The offloadmask of X dictates where to move memory
1552: If X GPU data is valid, then move Y data on GPU if needed
1553: Otherwise, move back to the CPU */
1554: switch (X->offloadmask) {
1555: case PETSC_OFFLOAD_BOTH:
1556: if (ymask == PETSC_OFFLOAD_CPU) {
1557: VecHIPResetArray(*Y);
1558: } else if (ymask == PETSC_OFFLOAD_GPU) {
1559: X->offloadmask = PETSC_OFFLOAD_GPU;
1560: }
1561: break;
1562: case PETSC_OFFLOAD_GPU:
1563: if (ymask == PETSC_OFFLOAD_CPU) VecHIPResetArray(*Y);
1564: break;
1565: case PETSC_OFFLOAD_CPU:
1566: if (ymask == PETSC_OFFLOAD_GPU) VecResetArray(*Y);
1567: break;
1568: case PETSC_OFFLOAD_UNALLOCATED:
1569: case PETSC_OFFLOAD_KOKKOS:
1570: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "This should not happen");
1571: }
1572: #endif
1573: } else {
1574: /* If OpenCL vecs updated the device memory, this triggers a copy on the CPU */
1575: VecResetArray(*Y);
1576: }
1577: PetscObjectStateIncrease((PetscObject)X);
1578: }
1579: }
1580: }
1581: VecDestroy(Y);
1582: return 0;
1583: }
1585: /*@
1586: VecCreateLocalVector - Creates a vector object suitable for use with VecGetLocalVector() and friends. You must call VecDestroy() when the
1587: vector is no longer needed.
1589: Not collective.
1591: Input parameter:
1592: . v - The vector for which the local vector is desired.
1594: Output parameter:
1595: . w - Upon exit this contains the local vector.
1597: Level: beginner
1599: .seealso: `VecGetLocalVectorRead()`, `VecRestoreLocalVectorRead()`, `VecGetLocalVector()`, `VecRestoreLocalVector()`
1600: @*/
1601: PetscErrorCode VecCreateLocalVector(Vec v, Vec *w)
1602: {
1603: PetscMPIInt size;
1607: MPI_Comm_size(PetscObjectComm((PetscObject)v), &size);
1608: if (size == 1) VecDuplicate(v, w);
1609: else if (v->ops->createlocalvector) PetscUseTypeMethod(v, createlocalvector, w);
1610: else {
1611: VecType type;
1612: PetscInt n;
1614: VecCreate(PETSC_COMM_SELF, w);
1615: VecGetLocalSize(v, &n);
1616: VecSetSizes(*w, n, n);
1617: VecGetBlockSize(v, &n);
1618: VecSetBlockSize(*w, n);
1619: VecGetType(v, &type);
1620: VecSetType(*w, type);
1621: }
1622: return 0;
1623: }
1625: /*@
1626: VecGetLocalVectorRead - Maps the local portion of a vector into a
1627: vector. You must call VecRestoreLocalVectorRead() when the local
1628: vector is no longer needed.
1630: Not collective.
1632: Input parameter:
1633: . v - The vector for which the local vector is desired.
1635: Output parameter:
1636: . w - Upon exit this contains the local vector.
1638: Level: beginner
1640: Notes:
1641: This function is similar to VecGetArrayRead() which maps the local
1642: portion into a raw pointer. VecGetLocalVectorRead() is usually
1643: almost as efficient as VecGetArrayRead() but in certain circumstances
1644: VecGetLocalVectorRead() can be much more efficient than
1645: VecGetArrayRead(). This is because the construction of a contiguous
1646: array representing the vector data required by VecGetArrayRead() can
1647: be an expensive operation for certain vector types. For example, for
1648: GPU vectors VecGetArrayRead() requires that the data between device
1649: and host is synchronized.
1651: Unlike VecGetLocalVector(), this routine is not collective and
1652: preserves cached information.
1654: .seealso: `VecCreateLocalVector()`, `VecRestoreLocalVectorRead()`, `VecGetLocalVector()`, `VecGetArrayRead()`, `VecGetArray()`
1655: @*/
1656: PetscErrorCode VecGetLocalVectorRead(Vec v, Vec w)
1657: {
1660: VecCheckSameLocalSize(v, 1, w, 2);
1661: if (v->ops->getlocalvectorread) {
1662: PetscUseTypeMethod(v, getlocalvectorread, w);
1663: } else {
1664: PetscScalar *a;
1666: VecGetArrayRead(v, (const PetscScalar **)&a);
1667: VecPlaceArray(w, a);
1668: }
1669: PetscObjectStateIncrease((PetscObject)w);
1670: VecLockReadPush(v);
1671: VecLockReadPush(w);
1672: return 0;
1673: }
1675: /*@
1676: VecRestoreLocalVectorRead - Unmaps the local portion of a vector
1677: previously mapped into a vector using VecGetLocalVectorRead().
1679: Not collective.
1681: Input parameter:
1682: + v - The local portion of this vector was previously mapped into w using VecGetLocalVectorRead().
1683: - w - The vector into which the local portion of v was mapped.
1685: Level: beginner
1687: .seealso: `VecCreateLocalVector()`, `VecGetLocalVectorRead()`, `VecGetLocalVector()`, `VecGetArrayRead()`, `VecGetArray()`
1688: @*/
1689: PetscErrorCode VecRestoreLocalVectorRead(Vec v, Vec w)
1690: {
1693: if (v->ops->restorelocalvectorread) {
1694: PetscUseTypeMethod(v, restorelocalvectorread, w);
1695: } else {
1696: const PetscScalar *a;
1698: VecGetArrayRead(w, &a);
1699: VecRestoreArrayRead(v, &a);
1700: VecResetArray(w);
1701: }
1702: VecLockReadPop(v);
1703: VecLockReadPop(w);
1704: PetscObjectStateIncrease((PetscObject)w);
1705: return 0;
1706: }
1708: /*@
1709: VecGetLocalVector - Maps the local portion of a vector into a
1710: vector.
1712: Collective on v, not collective on w.
1714: Input parameter:
1715: . v - The vector for which the local vector is desired.
1717: Output parameter:
1718: . w - Upon exit this contains the local vector.
1720: Level: beginner
1722: Notes:
1723: This function is similar to VecGetArray() which maps the local
1724: portion into a raw pointer. VecGetLocalVector() is usually about as
1725: efficient as VecGetArray() but in certain circumstances
1726: VecGetLocalVector() can be much more efficient than VecGetArray().
1727: This is because the construction of a contiguous array representing
1728: the vector data required by VecGetArray() can be an expensive
1729: operation for certain vector types. For example, for GPU vectors
1730: VecGetArray() requires that the data between device and host is
1731: synchronized.
1733: .seealso: `VecCreateLocalVector()`, `VecRestoreLocalVector()`, `VecGetLocalVectorRead()`, `VecGetArrayRead()`, `VecGetArray()`
1734: @*/
1735: PetscErrorCode VecGetLocalVector(Vec v, Vec w)
1736: {
1739: VecCheckSameLocalSize(v, 1, w, 2);
1740: if (v->ops->getlocalvector) {
1741: PetscUseTypeMethod(v, getlocalvector, w);
1742: } else {
1743: PetscScalar *a;
1745: VecGetArray(v, &a);
1746: VecPlaceArray(w, a);
1747: }
1748: PetscObjectStateIncrease((PetscObject)w);
1749: return 0;
1750: }
1752: /*@
1753: VecRestoreLocalVector - Unmaps the local portion of a vector
1754: previously mapped into a vector using VecGetLocalVector().
1756: Logically collective.
1758: Input parameter:
1759: + v - The local portion of this vector was previously mapped into w using VecGetLocalVector().
1760: - w - The vector into which the local portion of v was mapped.
1762: Level: beginner
1764: .seealso: `VecCreateLocalVector()`, `VecGetLocalVector()`, `VecGetLocalVectorRead()`, `VecRestoreLocalVectorRead()`, `LocalVectorRead()`, `VecGetArrayRead()`, `VecGetArray()`
1765: @*/
1766: PetscErrorCode VecRestoreLocalVector(Vec v, Vec w)
1767: {
1770: if (v->ops->restorelocalvector) {
1771: PetscUseTypeMethod(v, restorelocalvector, w);
1772: } else {
1773: PetscScalar *a;
1774: VecGetArray(w, &a);
1775: VecRestoreArray(v, &a);
1776: VecResetArray(w);
1777: }
1778: PetscObjectStateIncrease((PetscObject)w);
1779: PetscObjectStateIncrease((PetscObject)v);
1780: return 0;
1781: }
1783: /*@C
1784: VecGetArray - Returns a pointer to a contiguous array that contains this
1785: processor's portion of the vector data. For the standard PETSc
1786: vectors, VecGetArray() returns a pointer to the local data array and
1787: does not use any copies. If the underlying vector data is not stored
1788: in a contiguous array this routine will copy the data to a contiguous
1789: array and return a pointer to that. You MUST call VecRestoreArray()
1790: when you no longer need access to the array.
1792: Logically Collective on Vec
1794: Input Parameter:
1795: . x - the vector
1797: Output Parameter:
1798: . a - location to put pointer to the array
1800: Fortran Note:
1801: This routine is used differently from Fortran 77
1802: .vb
1803: Vec x
1804: PetscScalar x_array(1)
1805: PetscOffset i_x
1806: PetscErrorCode ierr
1807: call VecGetArray(x,x_array,i_x,ierr)
1809: Access first local entry in vector with
1810: value = x_array(i_x + 1)
1812: ...... other code
1813: call VecRestoreArray(x,x_array,i_x,ierr)
1814: .ve
1815: For Fortran 90 see `VecGetArrayF90()`
1817: See the Fortran chapter of the users manual and
1818: petsc/src/snes/tutorials/ex5f.F for details.
1820: Level: beginner
1822: .seealso: `VecRestoreArray()`, `VecGetArrayRead()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecGetArrayReadF90()`, `VecPlaceArray()`, `VecGetArray2d()`,
1823: `VecGetArrayPair()`, `VecRestoreArrayPair()`, `VecGetArrayWrite()`, `VecRestoreArrayWrite()`
1824: @*/
1825: PetscErrorCode VecGetArray(Vec x, PetscScalar **a)
1826: {
1828: VecSetErrorIfLocked(x, 1);
1829: if (x->ops->getarray) { /* The if-else order matters! VECNEST, VECCUDA etc should have ops->getarray while VECCUDA etc are petscnative */
1830: PetscUseTypeMethod(x, getarray, a);
1831: } else if (x->petscnative) { /* VECSTANDARD */
1832: *a = *((PetscScalar **)x->data);
1833: } else SETERRQ(PetscObjectComm((PetscObject)x), PETSC_ERR_SUP, "Cannot get array for vector type \"%s\"", ((PetscObject)x)->type_name);
1834: return 0;
1835: }
1837: /*@C
1838: VecRestoreArray - Restores a vector after VecGetArray() has been called.
1840: Logically Collective on Vec
1842: Input Parameters:
1843: + x - the vector
1844: - a - location of pointer to array obtained from VecGetArray()
1846: Level: beginner
1848: .seealso: `VecGetArray()`, `VecRestoreArrayRead()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecRestoreArrayReadF90()`, `VecPlaceArray()`, `VecRestoreArray2d()`,
1849: `VecGetArrayPair()`, `VecRestoreArrayPair()`
1850: @*/
1851: PetscErrorCode VecRestoreArray(Vec x, PetscScalar **a)
1852: {
1855: if (x->ops->restorearray) {
1856: PetscUseTypeMethod(x, restorearray, a);
1858: if (a) *a = NULL;
1859: PetscObjectStateIncrease((PetscObject)x);
1860: return 0;
1861: }
1862: /*@C
1863: VecGetArrayRead - Get read-only pointer to contiguous array containing this processor's portion of the vector data.
1865: Not Collective
1867: Input Parameter:
1868: . x - the vector
1870: Output Parameter:
1871: . a - the array
1873: Level: beginner
1875: Notes:
1876: The array must be returned using a matching call to VecRestoreArrayRead().
1878: Unlike VecGetArray(), this routine is not collective and preserves cached information like vector norms.
1880: Standard PETSc vectors use contiguous storage so that this routine does not perform a copy. Other vector
1881: implementations may require a copy, but must such implementations should cache the contiguous representation so that
1882: only one copy is performed when this routine is called multiple times in sequence.
1884: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayPair()`, `VecRestoreArrayPair()`
1885: @*/
1886: PetscErrorCode VecGetArrayRead(Vec x, const PetscScalar **a)
1887: {
1890: if (x->ops->getarrayread) {
1891: PetscUseTypeMethod(x, getarrayread, a);
1892: } else if (x->ops->getarray) {
1893: /* VECNEST, VECCUDA, VECKOKKOS etc */
1894: PetscUseTypeMethod(x, getarray, (PetscScalar **)a);
1895: } else if (x->petscnative) {
1896: /* VECSTANDARD */
1897: *a = *((PetscScalar **)x->data);
1898: } else SETERRQ(PetscObjectComm((PetscObject)x), PETSC_ERR_SUP, "Cannot get array read for vector type \"%s\"", ((PetscObject)x)->type_name);
1899: return 0;
1900: }
1902: /*@C
1903: VecRestoreArrayRead - Restore array obtained with VecGetArrayRead()
1905: Not Collective
1907: Input Parameters:
1908: + vec - the vector
1909: - array - the array
1911: Level: beginner
1913: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayPair()`, `VecRestoreArrayPair()`
1914: @*/
1915: PetscErrorCode VecRestoreArrayRead(Vec x, const PetscScalar **a)
1916: {
1919: if (x->petscnative) { /* VECSTANDARD, VECCUDA, VECKOKKOS etc */
1920: /* nothing */
1921: } else if (x->ops->restorearrayread) { /* VECNEST */
1922: PetscUseTypeMethod(x, restorearrayread, a);
1923: } else { /* No one? */
1924: PetscUseTypeMethod(x, restorearray, (PetscScalar **)a);
1925: }
1926: if (a) *a = NULL;
1927: return 0;
1928: }
1930: /*@C
1931: VecGetArrayWrite - Returns a pointer to a contiguous array that WILL contains this
1932: processor's portion of the vector data. The values in this array are NOT valid, the routine calling this
1933: routine is responsible for putting values into the array; any values it does not set will be invalid
1935: Logically Collective on Vec
1937: Input Parameter:
1938: . x - the vector
1940: Output Parameter:
1941: . a - location to put pointer to the array
1943: Level: intermediate
1945: This is for vectors associate with GPUs, the vector is not copied up before giving access. If you need correct
1946: values in the array use VecGetArray()
1948: .seealso: `VecRestoreArray()`, `VecGetArrayRead()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecGetArrayReadF90()`, `VecPlaceArray()`, `VecGetArray2d()`,
1949: `VecGetArrayPair()`, `VecRestoreArrayPair()`, `VecGetArray()`, `VecRestoreArrayWrite()`
1950: @*/
1951: PetscErrorCode VecGetArrayWrite(Vec x, PetscScalar **a)
1952: {
1955: VecSetErrorIfLocked(x, 1);
1956: if (x->ops->getarraywrite) {
1957: PetscUseTypeMethod(x, getarraywrite, a);
1958: } else {
1959: VecGetArray(x, a);
1960: }
1961: return 0;
1962: }
1964: /*@C
1965: VecRestoreArrayWrite - Restores a vector after VecGetArrayWrite() has been called.
1967: Logically Collective on Vec
1969: Input Parameters:
1970: + x - the vector
1971: - a - location of pointer to array obtained from VecGetArray()
1973: Level: beginner
1975: .seealso: `VecGetArray()`, `VecRestoreArrayRead()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecRestoreArrayReadF90()`, `VecPlaceArray()`, `VecRestoreArray2d()`,
1976: `VecGetArrayPair()`, `VecRestoreArrayPair()`, `VecGetArrayWrite()`
1977: @*/
1978: PetscErrorCode VecRestoreArrayWrite(Vec x, PetscScalar **a)
1979: {
1982: if (x->ops->restorearraywrite) {
1983: PetscUseTypeMethod(x, restorearraywrite, a);
1984: } else if (x->ops->restorearray) {
1985: PetscUseTypeMethod(x, restorearray, a);
1986: }
1987: if (a) *a = NULL;
1988: PetscObjectStateIncrease((PetscObject)x);
1989: return 0;
1990: }
1992: /*@C
1993: VecGetArrays - Returns a pointer to the arrays in a set of vectors
1994: that were created by a call to VecDuplicateVecs(). You MUST call
1995: VecRestoreArrays() when you no longer need access to the array.
1997: Logically Collective on Vec
1999: Input Parameters:
2000: + x - the vectors
2001: - n - the number of vectors
2003: Output Parameter:
2004: . a - location to put pointer to the array
2006: Fortran Note:
2007: This routine is not supported in Fortran.
2009: Level: intermediate
2011: .seealso: `VecGetArray()`, `VecRestoreArrays()`
2012: @*/
2013: PetscErrorCode VecGetArrays(const Vec x[], PetscInt n, PetscScalar **a[])
2014: {
2015: PetscInt i;
2016: PetscScalar **q;
2022: PetscMalloc1(n, &q);
2023: for (i = 0; i < n; ++i) VecGetArray(x[i], &q[i]);
2024: *a = q;
2025: return 0;
2026: }
2028: /*@C
2029: VecRestoreArrays - Restores a group of vectors after VecGetArrays()
2030: has been called.
2032: Logically Collective on Vec
2034: Input Parameters:
2035: + x - the vector
2036: . n - the number of vectors
2037: - a - location of pointer to arrays obtained from VecGetArrays()
2039: Notes:
2040: For regular PETSc vectors this routine does not involve any copies. For
2041: any special vectors that do not store local vector data in a contiguous
2042: array, this routine will copy the data back into the underlying
2043: vector data structure from the arrays obtained with VecGetArrays().
2045: Fortran Note:
2046: This routine is not supported in Fortran.
2048: Level: intermediate
2050: .seealso: `VecGetArrays()`, `VecRestoreArray()`
2051: @*/
2052: PetscErrorCode VecRestoreArrays(const Vec x[], PetscInt n, PetscScalar **a[])
2053: {
2054: PetscInt i;
2055: PetscScalar **q = *a;
2061: for (i = 0; i < n; ++i) VecRestoreArray(x[i], &q[i]);
2062: PetscFree(q);
2063: return 0;
2064: }
2066: /*@C
2067: VecGetArrayAndMemType - Like VecGetArray(), but if this is a standard device vector (e.g., VECCUDA), the returned pointer will be a device
2068: pointer to the device memory that contains this processor's portion of the vector data. Device data is guaranteed to have the latest value.
2069: Otherwise, when this is a host vector (e.g., VECMPI), this routine functions the same as VecGetArray() and returns a host pointer.
2071: For VECKOKKOS, if Kokkos is configured without device (e.g., use serial or openmp), per this function, the vector works like VECSEQ/VECMPI;
2072: otherwise, it works like VECCUDA or VECHIP etc.
2074: Logically Collective on Vec
2076: Input Parameter:
2077: . x - the vector
2079: Output Parameters:
2080: + a - location to put pointer to the array
2081: - mtype - memory type of the array
2083: Level: beginner
2085: .seealso: `VecRestoreArrayAndMemType()`, `VecGetArrayReadAndMemType()`, `VecGetArrayWriteAndMemType()`, `VecRestoreArray()`, `VecGetArrayRead()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecGetArrayReadF90()`,
2086: `VecPlaceArray()`, `VecGetArray2d()`, `VecGetArrayPair()`, `VecRestoreArrayPair()`, `VecGetArrayWrite()`, `VecRestoreArrayWrite()`
2087: @*/
2088: PetscErrorCode VecGetArrayAndMemType(Vec x, PetscScalar **a, PetscMemType *mtype)
2089: {
2094: VecSetErrorIfLocked(x, 1);
2095: if (x->ops->getarrayandmemtype) {
2096: /* VECCUDA, VECKOKKOS etc */
2097: PetscUseTypeMethod(x, getarrayandmemtype, a, mtype);
2098: } else {
2099: /* VECSTANDARD, VECNEST, VECVIENNACL */
2100: VecGetArray(x, a);
2101: if (mtype) *mtype = PETSC_MEMTYPE_HOST;
2102: }
2103: return 0;
2104: }
2106: /*@C
2107: VecRestoreArrayAndMemType - Restores a vector after VecGetArrayAndMemType() has been called.
2109: Logically Collective on Vec
2111: Input Parameters:
2112: + x - the vector
2113: - a - location of pointer to array obtained from VecGetArrayAndMemType()
2115: Level: beginner
2117: .seealso: `VecGetArrayAndMemType()`, `VecGetArray()`, `VecRestoreArrayRead()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecRestoreArrayReadF90()`,
2118: `VecPlaceArray()`, `VecRestoreArray2d()`, `VecGetArrayPair()`, `VecRestoreArrayPair()`
2119: @*/
2120: PetscErrorCode VecRestoreArrayAndMemType(Vec x, PetscScalar **a)
2121: {
2125: if (x->ops->restorearrayandmemtype) {
2126: /* VECCUDA, VECKOKKOS etc */
2127: PetscUseTypeMethod(x, restorearrayandmemtype, a);
2128: } else {
2129: /* VECNEST, VECVIENNACL */
2130: VecRestoreArray(x, a);
2131: } /* VECSTANDARD does nothing */
2132: if (a) *a = NULL;
2133: PetscObjectStateIncrease((PetscObject)x);
2134: return 0;
2135: }
2137: /*@C
2138: 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().
2140: Not Collective
2142: Input Parameter:
2143: . x - the vector
2145: Output Parameters:
2146: + a - the array
2147: - mtype - memory type of the array
2149: Level: beginner
2151: Notes:
2152: The array must be returned using a matching call to VecRestoreArrayReadAndMemType().
2154: .seealso: `VecRestoreArrayReadAndMemType()`, `VecGetArrayAndMemType()`, `VecGetArrayWriteAndMemType()`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayPair()`, `VecRestoreArrayPair()`, `VecGetArrayAndMemType()`
2155: @*/
2156: PetscErrorCode VecGetArrayReadAndMemType(Vec x, const PetscScalar **a, PetscMemType *mtype)
2157: {
2162: if (x->ops->getarrayreadandmemtype) {
2163: /* VECCUDA/VECHIP though they are also petscnative */
2164: PetscUseTypeMethod(x, getarrayreadandmemtype, a, mtype);
2165: } else if (x->ops->getarrayandmemtype) {
2166: /* VECKOKKOS */
2167: PetscUseTypeMethod(x, getarrayandmemtype, (PetscScalar **)a, mtype);
2168: } else {
2169: VecGetArrayRead(x, a);
2170: if (mtype) *mtype = PETSC_MEMTYPE_HOST;
2171: }
2172: return 0;
2173: }
2175: /*@C
2176: VecRestoreArrayReadAndMemType - Restore array obtained with VecGetArrayReadAndMemType()
2178: Not Collective
2180: Input Parameters:
2181: + vec - the vector
2182: - array - the array
2184: Level: beginner
2186: .seealso: `VecGetArrayReadAndMemType()`, `VecRestoreArrayAndMemType()`, `VecRestoreArrayWriteAndMemType()`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayPair()`, `VecRestoreArrayPair()`
2187: @*/
2188: PetscErrorCode VecRestoreArrayReadAndMemType(Vec x, const PetscScalar **a)
2189: {
2193: if (x->ops->restorearrayreadandmemtype) {
2194: /* VECCUDA/VECHIP */
2195: PetscUseTypeMethod(x, restorearrayreadandmemtype, a);
2196: } else if (!x->petscnative) {
2197: /* VECNEST */
2198: VecRestoreArrayRead(x, a);
2199: }
2200: if (a) *a = NULL;
2201: return 0;
2202: }
2204: /*@C
2205: VecGetArrayWriteAndMemType - Like VecGetArrayWrite(), but if this is a device vector it will aways return
2206: a device pointer to the device memory that contains this processor's portion of the vector data.
2208: Not Collective
2210: Input Parameter:
2211: . x - the vector
2213: Output Parameters:
2214: + a - the array
2215: - mtype - memory type of the array
2217: Level: beginner
2219: Notes:
2220: The array must be returned using a matching call to VecRestoreArrayWriteAndMemType(), where it will label the device memory as most recent.
2222: .seealso: `VecRestoreArrayWriteAndMemType()`, `VecGetArrayReadAndMemType()`, `VecGetArrayAndMemType()`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayPair()`, `VecRestoreArrayPair()`,
2223: @*/
2224: PetscErrorCode VecGetArrayWriteAndMemType(Vec x, PetscScalar **a, PetscMemType *mtype)
2225: {
2228: VecSetErrorIfLocked(x, 1);
2231: if (x->ops->getarraywriteandmemtype) {
2232: /* VECCUDA, VECHIP, VECKOKKOS etc, though they are also petscnative */
2233: PetscUseTypeMethod(x, getarraywriteandmemtype, a, mtype);
2234: } else if (x->ops->getarrayandmemtype) {
2235: VecGetArrayAndMemType(x, a, mtype);
2236: } else {
2237: /* VECNEST, VECVIENNACL */
2238: VecGetArrayWrite(x, a);
2239: if (mtype) *mtype = PETSC_MEMTYPE_HOST;
2240: }
2241: return 0;
2242: }
2244: /*@C
2245: VecRestoreArrayWriteAndMemType - Restore array obtained with VecGetArrayWriteAndMemType()
2247: Not Collective
2249: Input Parameters:
2250: + vec - the vector
2251: - array - the array
2253: Level: beginner
2255: .seealso: `VecGetArrayWriteAndMemType()`, `VecRestoreArrayAndMemType()`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayPair()`, `VecRestoreArrayPair()`
2256: @*/
2257: PetscErrorCode VecRestoreArrayWriteAndMemType(Vec x, PetscScalar **a)
2258: {
2261: VecSetErrorIfLocked(x, 1);
2263: if (x->ops->restorearraywriteandmemtype) {
2264: /* VECCUDA/VECHIP */
2265: PetscMemType PETSC_UNUSED mtype; // since this function doesn't accept a memtype?
2266: PetscUseTypeMethod(x, restorearraywriteandmemtype, a, &mtype);
2267: } else if (x->ops->restorearrayandmemtype) {
2268: VecRestoreArrayAndMemType(x, a);
2269: } else {
2270: VecRestoreArray(x, a);
2271: }
2272: if (a) *a = NULL;
2273: return 0;
2274: }
2276: /*@
2277: VecPlaceArray - Allows one to replace the array in a vector with an
2278: array provided by the user. This is useful to avoid copying an array
2279: into a vector.
2281: Not Collective
2283: Input Parameters:
2284: + vec - the vector
2285: - array - the array
2287: Notes:
2288: You can return to the original array with a call to `VecResetArray()`. `vec` does not take
2289: ownership of `array` in any way. The user must free `array` themselves but be careful not to
2290: do so before the vector has either been destroyed, had its original array restored with
2291: `VecResetArray()` or permanently replaced with `VecReplaceArray()`.
2293: Level: developer
2295: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecReplaceArray()`, `VecResetArray()`
2297: @*/
2298: PetscErrorCode VecPlaceArray(Vec vec, const PetscScalar array[])
2299: {
2303: PetscUseTypeMethod(vec, placearray, array);
2304: PetscObjectStateIncrease((PetscObject)vec);
2305: return 0;
2306: }
2308: /*@C
2309: VecReplaceArray - Allows one to replace the array in a vector with an
2310: array provided by the user. This is useful to avoid copying an array
2311: into a vector.
2313: Not Collective
2315: Input Parameters:
2316: + vec - the vector
2317: - array - the array
2319: Notes:
2320: This permanently replaces the array and frees the memory associated
2321: with the old array.
2323: The memory passed in MUST be obtained with PetscMalloc() and CANNOT be
2324: freed by the user. It will be freed when the vector is destroyed.
2326: Not supported from Fortran
2328: Level: developer
2330: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecPlaceArray()`, `VecResetArray()`
2332: @*/
2333: PetscErrorCode VecReplaceArray(Vec vec, const PetscScalar array[])
2334: {
2337: PetscUseTypeMethod(vec, replacearray, array);
2338: PetscObjectStateIncrease((PetscObject)vec);
2339: return 0;
2340: }
2342: /*@C
2343: VecCUDAGetArray - Provides access to the CUDA buffer inside a vector.
2345: This function has semantics similar to VecGetArray(): the pointer
2346: returned by this function points to a consistent view of the vector
2347: data. This may involve a copy operation of data from the host to the
2348: device if the data on the device is out of date. If the device
2349: memory hasn't been allocated previously it will be allocated as part
2350: of this function call. VecCUDAGetArray() assumes that
2351: the user will modify the vector data. This is similar to
2352: intent(inout) in fortran.
2354: The CUDA device pointer has to be released by calling
2355: VecCUDARestoreArray(). Upon restoring the vector data
2356: the data on the host will be marked as out of date. A subsequent
2357: access of the host data will thus incur a data transfer from the
2358: device to the host.
2360: Input Parameter:
2361: . v - the vector
2363: Output Parameter:
2364: . a - the CUDA device pointer
2366: Fortran note:
2367: This function is not currently available from Fortran.
2369: Level: intermediate
2371: .seealso: `VecCUDARestoreArray()`, `VecCUDAGetArrayRead()`, `VecCUDAGetArrayWrite()`, `VecGetArray()`, `VecGetArrayRead()`
2372: @*/
2373: PETSC_EXTERN PetscErrorCode VecCUDAGetArray(Vec v, PetscScalar **a)
2374: {
2376: #if defined(PETSC_HAVE_CUDA)
2377: {
2378: VecCUDACopyToGPU(v);
2379: *a = ((Vec_CUDA *)v->spptr)->GPUarray;
2380: }
2381: #endif
2382: return 0;
2383: }
2385: /*@C
2386: VecCUDARestoreArray - Restore a CUDA device pointer previously acquired with VecCUDAGetArray().
2388: This marks the host data as out of date. Subsequent access to the
2389: vector data on the host side with for instance VecGetArray() incurs a
2390: data transfer.
2392: Input Parameters:
2393: + v - the vector
2394: - a - the CUDA device pointer. This pointer is invalid after
2395: VecCUDARestoreArray() returns.
2397: Fortran note:
2398: This function is not currently available from Fortran.
2400: Level: intermediate
2402: .seealso: `VecCUDAGetArray()`, `VecCUDAGetArrayRead()`, `VecCUDAGetArrayWrite()`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayRead()`
2403: @*/
2404: PETSC_EXTERN PetscErrorCode VecCUDARestoreArray(Vec v, PetscScalar **a)
2405: {
2407: #if defined(PETSC_HAVE_CUDA)
2408: v->offloadmask = PETSC_OFFLOAD_GPU;
2409: #endif
2410: PetscObjectStateIncrease((PetscObject)v);
2411: return 0;
2412: }
2414: /*@C
2415: VecCUDAGetArrayRead - Provides read access to the CUDA buffer inside a vector.
2417: This function is analogous to VecGetArrayRead(): The pointer
2418: returned by this function points to a consistent view of the vector
2419: data. This may involve a copy operation of data from the host to the
2420: device if the data on the device is out of date. If the device
2421: memory hasn't been allocated previously it will be allocated as part
2422: of this function call. VecCUDAGetArrayRead() assumes that the
2423: user will not modify the vector data. This is analgogous to
2424: intent(in) in Fortran.
2426: The CUDA device pointer has to be released by calling
2427: VecCUDARestoreArrayRead(). If the data on the host side was
2428: previously up to date it will remain so, i.e. data on both the device
2429: and the host is up to date. Accessing data on the host side does not
2430: incur a device to host data transfer.
2432: Input Parameter:
2433: . v - the vector
2435: Output Parameter:
2436: . a - the CUDA pointer.
2438: Fortran note:
2439: This function is not currently available from Fortran.
2441: Level: intermediate
2443: .seealso: `VecCUDARestoreArrayRead()`, `VecCUDAGetArray()`, `VecCUDAGetArrayWrite()`, `VecGetArray()`, `VecGetArrayRead()`
2444: @*/
2445: PETSC_EXTERN PetscErrorCode VecCUDAGetArrayRead(Vec v, const PetscScalar **a)
2446: {
2447: VecCUDAGetArray(v, (PetscScalar **)a);
2448: return 0;
2449: }
2451: /*@C
2452: VecCUDARestoreArrayRead - Restore a CUDA device pointer previously acquired with VecCUDAGetArrayRead().
2454: If the data on the host side was previously up to date it will remain
2455: so, i.e. data on both the device and the host is up to date.
2456: Accessing data on the host side e.g. with VecGetArray() does not
2457: incur a device to host data transfer.
2459: Input Parameters:
2460: + v - the vector
2461: - a - the CUDA device pointer. This pointer is invalid after
2462: VecCUDARestoreArrayRead() returns.
2464: Fortran note:
2465: This function is not currently available from Fortran.
2467: Level: intermediate
2469: .seealso: `VecCUDAGetArrayRead()`, `VecCUDAGetArrayWrite()`, `VecCUDAGetArray()`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayRead()`
2470: @*/
2471: PETSC_EXTERN PetscErrorCode VecCUDARestoreArrayRead(Vec v, const PetscScalar **a)
2472: {
2474: *a = NULL;
2475: return 0;
2476: }
2478: /*@C
2479: VecCUDAGetArrayWrite - Provides write access to the CUDA buffer inside a vector.
2481: The data pointed to by the device pointer is uninitialized. The user
2482: may not read from this data. Furthermore, the entire array needs to
2483: be filled by the user to obtain well-defined behaviour. The device
2484: memory will be allocated by this function if it hasn't been allocated
2485: previously. This is analogous to intent(out) in Fortran.
2487: The device pointer needs to be released with
2488: VecCUDARestoreArrayWrite(). When the pointer is released the
2489: host data of the vector is marked as out of data. Subsequent access
2490: of the host data with e.g. VecGetArray() incurs a device to host data
2491: transfer.
2493: Input Parameter:
2494: . v - the vector
2496: Output Parameter:
2497: . a - the CUDA pointer
2499: Fortran note:
2500: This function is not currently available from Fortran.
2502: Level: advanced
2504: .seealso: `VecCUDARestoreArrayWrite()`, `VecCUDAGetArray()`, `VecCUDAGetArrayRead()`, `VecCUDAGetArrayWrite()`, `VecGetArray()`, `VecGetArrayRead()`
2505: @*/
2506: PETSC_EXTERN PetscErrorCode VecCUDAGetArrayWrite(Vec v, PetscScalar **a)
2507: {
2509: #if defined(PETSC_HAVE_CUDA)
2510: {
2511: VecCUDAAllocateCheck(v);
2512: *a = ((Vec_CUDA *)v->spptr)->GPUarray;
2513: }
2514: #endif
2515: return 0;
2516: }
2518: /*@C
2519: VecCUDARestoreArrayWrite - Restore a CUDA device pointer previously acquired with VecCUDAGetArrayWrite().
2521: Data on the host will be marked as out of date. Subsequent access of
2522: the data on the host side e.g. with VecGetArray() will incur a device
2523: to host data transfer.
2525: Input Parameters:
2526: + v - the vector
2527: - a - the CUDA device pointer. This pointer is invalid after
2528: VecCUDARestoreArrayWrite() returns.
2530: Fortran note:
2531: This function is not currently available from Fortran.
2533: Level: intermediate
2535: .seealso: `VecCUDAGetArrayWrite()`, `VecCUDAGetArray()`, `VecCUDAGetArrayRead()`, `VecCUDAGetArrayWrite()`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayRead()`
2536: @*/
2537: PETSC_EXTERN PetscErrorCode VecCUDARestoreArrayWrite(Vec v, PetscScalar **a)
2538: {
2540: #if defined(PETSC_HAVE_CUDA)
2541: v->offloadmask = PETSC_OFFLOAD_GPU;
2542: if (a) *a = NULL;
2543: #endif
2544: PetscObjectStateIncrease((PetscObject)v);
2545: return 0;
2546: }
2548: /*@C
2549: VecCUDAPlaceArray - Allows one to replace the GPU array in a vector with a
2550: GPU array provided by the user. This is useful to avoid copying an
2551: array into a vector.
2553: Not Collective
2555: Input Parameters:
2556: + vec - the vector
2557: - array - the GPU array
2559: Notes:
2560: You can return to the original GPU array with a call to VecCUDAResetArray()
2561: It is not possible to use VecCUDAPlaceArray() and VecPlaceArray() at the
2562: same time on the same vector.
2564: `vec` does not take ownership of `array` in any way. The user must free `array` themselves
2565: but be careful not to do so before the vector has either been destroyed, had its original
2566: array restored with `VecCUDAResetArray()` or permanently replaced with
2567: `VecCUDAReplaceArray()`.
2569: Level: developer
2571: .seealso: `VecPlaceArray()`, `VecGetArray()`, `VecRestoreArray()`, `VecReplaceArray()`, `VecResetArray()`, `VecCUDAResetArray()`, `VecCUDAReplaceArray()`
2573: @*/
2574: PetscErrorCode VecCUDAPlaceArray(Vec vin, const PetscScalar a[])
2575: {
2577: #if defined(PETSC_HAVE_CUDA)
2578: VecCUDACopyToGPU(vin);
2580: ((Vec_Seq *)vin->data)->unplacedarray = (PetscScalar *)((Vec_CUDA *)vin->spptr)->GPUarray; /* save previous GPU array so reset can bring it back */
2581: ((Vec_CUDA *)vin->spptr)->GPUarray = (PetscScalar *)a;
2582: vin->offloadmask = PETSC_OFFLOAD_GPU;
2583: #endif
2584: PetscObjectStateIncrease((PetscObject)vin);
2585: return 0;
2586: }
2588: /*@C
2589: VecCUDAReplaceArray - Allows one to replace the GPU array in a vector
2590: with a GPU array provided by the user. This is useful to avoid copying
2591: a GPU array into a vector.
2593: Not Collective
2595: Input Parameters:
2596: + vec - the vector
2597: - array - the GPU array
2599: Notes:
2600: This permanently replaces the GPU array and frees the memory associated
2601: with the old GPU array.
2603: The memory passed in CANNOT be freed by the user. It will be freed
2604: when the vector is destroyed.
2606: Not supported from Fortran
2608: Level: developer
2610: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecPlaceArray()`, `VecResetArray()`, `VecCUDAResetArray()`, `VecCUDAPlaceArray()`, `VecReplaceArray()`
2612: @*/
2613: PetscErrorCode VecCUDAReplaceArray(Vec vin, const PetscScalar a[])
2614: {
2615: #if defined(PETSC_HAVE_CUDA)
2616: #endif
2619: #if defined(PETSC_HAVE_CUDA)
2620: if (((Vec_CUDA *)vin->spptr)->GPUarray_allocated) cudaFree(((Vec_CUDA *)vin->spptr)->GPUarray_allocated);
2621: ((Vec_CUDA *)vin->spptr)->GPUarray_allocated = ((Vec_CUDA *)vin->spptr)->GPUarray = (PetscScalar *)a;
2622: vin->offloadmask = PETSC_OFFLOAD_GPU;
2623: #endif
2624: PetscObjectStateIncrease((PetscObject)vin);
2625: return 0;
2626: }
2628: /*@C
2629: VecCUDAResetArray - Resets a vector to use its default memory. Call this
2630: after the use of VecCUDAPlaceArray().
2632: Not Collective
2634: Input Parameters:
2635: . vec - the vector
2637: Level: developer
2639: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecReplaceArray()`, `VecPlaceArray()`, `VecResetArray()`, `VecCUDAPlaceArray()`, `VecCUDAReplaceArray()`
2641: @*/
2642: PetscErrorCode VecCUDAResetArray(Vec vin)
2643: {
2645: #if defined(PETSC_HAVE_CUDA)
2646: VecCUDACopyToGPU(vin);
2647: ((Vec_CUDA *)vin->spptr)->GPUarray = (PetscScalar *)((Vec_Seq *)vin->data)->unplacedarray;
2648: ((Vec_Seq *)vin->data)->unplacedarray = 0;
2649: vin->offloadmask = PETSC_OFFLOAD_GPU;
2650: #endif
2651: PetscObjectStateIncrease((PetscObject)vin);
2652: return 0;
2653: }
2655: /*@C
2656: VecHIPGetArray - Provides access to the HIP buffer inside a vector.
2658: This function has semantics similar to VecGetArray(): the pointer
2659: returned by this function points to a consistent view of the vector
2660: data. This may involve a copy operation of data from the host to the
2661: device if the data on the device is out of date. If the device
2662: memory hasn't been allocated previously it will be allocated as part
2663: of this function call. VecHIPGetArray() assumes that
2664: the user will modify the vector data. This is similar to
2665: intent(inout) in fortran.
2667: The HIP device pointer has to be released by calling
2668: VecHIPRestoreArray(). Upon restoring the vector data
2669: the data on the host will be marked as out of date. A subsequent
2670: access of the host data will thus incur a data transfer from the
2671: device to the host.
2673: Input Parameter:
2674: . v - the vector
2676: Output Parameter:
2677: . a - the HIP device pointer
2679: Fortran note:
2680: This function is not currently available from Fortran.
2682: Level: intermediate
2684: .seealso: `VecHIPRestoreArray()`, `VecHIPGetArrayRead()`, `VecHIPGetArrayWrite()`, `VecGetArray()`, `VecGetArrayRead()`
2685: @*/
2686: PETSC_EXTERN PetscErrorCode VecHIPGetArray(Vec v, PetscScalar **a)
2687: {
2689: #if defined(PETSC_HAVE_HIP)
2690: *a = 0;
2691: VecHIPCopyToGPU(v);
2692: *a = ((Vec_HIP *)v->spptr)->GPUarray;
2693: #endif
2694: return 0;
2695: }
2697: /*@C
2698: VecHIPRestoreArray - Restore a HIP device pointer previously acquired with VecHIPGetArray().
2700: This marks the host data as out of date. Subsequent access to the
2701: vector data on the host side with for instance VecGetArray() incurs a
2702: data transfer.
2704: Input Parameters:
2705: + v - the vector
2706: - a - the HIP device pointer. This pointer is invalid after
2707: VecHIPRestoreArray() returns.
2709: Fortran note:
2710: This function is not currently available from Fortran.
2712: Level: intermediate
2714: .seealso: `VecHIPGetArray()`, `VecHIPGetArrayRead()`, `VecHIPGetArrayWrite()`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayRead()`
2715: @*/
2716: PETSC_EXTERN PetscErrorCode VecHIPRestoreArray(Vec v, PetscScalar **a)
2717: {
2719: #if defined(PETSC_HAVE_HIP)
2720: v->offloadmask = PETSC_OFFLOAD_GPU;
2721: #endif
2723: PetscObjectStateIncrease((PetscObject)v);
2724: return 0;
2725: }
2727: /*@C
2728: VecHIPGetArrayRead - Provides read access to the HIP buffer inside a vector.
2730: This function is analogous to VecGetArrayRead(): The pointer
2731: returned by this function points to a consistent view of the vector
2732: data. This may involve a copy operation of data from the host to the
2733: device if the data on the device is out of date. If the device
2734: memory hasn't been allocated previously it will be allocated as part
2735: of this function call. VecHIPGetArrayRead() assumes that the
2736: user will not modify the vector data. This is analgogous to
2737: intent(in) in Fortran.
2739: The HIP device pointer has to be released by calling
2740: VecHIPRestoreArrayRead(). If the data on the host side was
2741: previously up to date it will remain so, i.e. data on both the device
2742: and the host is up to date. Accessing data on the host side does not
2743: incur a device to host data transfer.
2745: Input Parameter:
2746: . v - the vector
2748: Output Parameter:
2749: . a - the HIP pointer.
2751: Fortran note:
2752: This function is not currently available from Fortran.
2754: Level: intermediate
2756: .seealso: `VecHIPRestoreArrayRead()`, `VecHIPGetArray()`, `VecHIPGetArrayWrite()`, `VecGetArray()`, `VecGetArrayRead()`
2757: @*/
2758: PETSC_EXTERN PetscErrorCode VecHIPGetArrayRead(Vec v, const PetscScalar **a)
2759: {
2761: #if defined(PETSC_HAVE_HIP)
2762: *a = 0;
2763: VecHIPCopyToGPU(v);
2764: *a = ((Vec_HIP *)v->spptr)->GPUarray;
2765: #endif
2766: return 0;
2767: }
2769: /*@C
2770: VecHIPRestoreArrayRead - Restore a HIP device pointer previously acquired with VecHIPGetArrayRead().
2772: If the data on the host side was previously up to date it will remain
2773: so, i.e. data on both the device and the host is up to date.
2774: Accessing data on the host side e.g. with VecGetArray() does not
2775: incur a device to host data transfer.
2777: Input Parameters:
2778: + v - the vector
2779: - a - the HIP device pointer. This pointer is invalid after
2780: VecHIPRestoreArrayRead() returns.
2782: Fortran note:
2783: This function is not currently available from Fortran.
2785: Level: intermediate
2787: .seealso: `VecHIPGetArrayRead()`, `VecHIPGetArrayWrite()`, `VecHIPGetArray()`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayRead()`
2788: @*/
2789: PETSC_EXTERN PetscErrorCode VecHIPRestoreArrayRead(Vec v, const PetscScalar **a)
2790: {
2792: *a = NULL;
2793: return 0;
2794: }
2796: /*@C
2797: VecHIPGetArrayWrite - Provides write access to the HIP buffer inside a vector.
2799: The data pointed to by the device pointer is uninitialized. The user
2800: may not read from this data. Furthermore, the entire array needs to
2801: be filled by the user to obtain well-defined behaviour. The device
2802: memory will be allocated by this function if it hasn't been allocated
2803: previously. This is analogous to intent(out) in Fortran.
2805: The device pointer needs to be released with
2806: VecHIPRestoreArrayWrite(). When the pointer is released the
2807: host data of the vector is marked as out of data. Subsequent access
2808: of the host data with e.g. VecGetArray() incurs a device to host data
2809: transfer.
2811: Input Parameter:
2812: . v - the vector
2814: Output Parameter:
2815: . a - the HIP pointer
2817: Fortran note:
2818: This function is not currently available from Fortran.
2820: Level: advanced
2822: .seealso: `VecHIPRestoreArrayWrite()`, `VecHIPGetArray()`, `VecHIPGetArrayRead()`, `VecHIPGetArrayWrite()`, `VecGetArray()`, `VecGetArrayRead()`
2823: @*/
2824: PETSC_EXTERN PetscErrorCode VecHIPGetArrayWrite(Vec v, PetscScalar **a)
2825: {
2827: #if defined(PETSC_HAVE_HIP)
2828: *a = 0;
2829: VecHIPAllocateCheck(v);
2830: *a = ((Vec_HIP *)v->spptr)->GPUarray;
2831: #endif
2832: return 0;
2833: }
2835: /*@C
2836: VecHIPRestoreArrayWrite - Restore a HIP device pointer previously acquired with VecHIPGetArrayWrite().
2838: Data on the host will be marked as out of date. Subsequent access of
2839: the data on the host side e.g. with VecGetArray() will incur a device
2840: to host data transfer.
2842: Input Parameters:
2843: + v - the vector
2844: - a - the HIP device pointer. This pointer is invalid after
2845: VecHIPRestoreArrayWrite() returns.
2847: Fortran note:
2848: This function is not currently available from Fortran.
2850: Level: intermediate
2852: .seealso: `VecHIPGetArrayWrite()`, `VecHIPGetArray()`, `VecHIPGetArrayRead()`, `VecHIPGetArrayWrite()`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayRead()`
2853: @*/
2854: PETSC_EXTERN PetscErrorCode VecHIPRestoreArrayWrite(Vec v, PetscScalar **a)
2855: {
2857: #if defined(PETSC_HAVE_HIP)
2858: v->offloadmask = PETSC_OFFLOAD_GPU;
2859: #endif
2861: PetscObjectStateIncrease((PetscObject)v);
2862: return 0;
2863: }
2865: /*@C
2866: VecHIPPlaceArray - Allows one to replace the GPU array in a vector with a
2867: GPU array provided by the user. This is useful to avoid copying an
2868: array into a vector.
2870: Not Collective
2872: Input Parameters:
2873: + vec - the vector
2874: - array - the GPU array
2876: Notes:
2877: You can return to the original GPU array with a call to VecHIPResetArray()
2878: It is not possible to use VecHIPPlaceArray() and VecPlaceArray() at the
2879: same time on the same vector.
2881: `vec` does not take ownership of `array` in any way. The user must free `array` themselves
2882: but be careful not to do so before the vector has either been destroyed, had its original
2883: array restored with `VecHIPResetArray()` or permanently replaced with
2884: `VecHIPReplaceArray()`.
2886: Level: developer
2888: .seealso: `VecPlaceArray()`, `VecGetArray()`, `VecRestoreArray()`, `VecReplaceArray()`, `VecResetArray()`, `VecHIPResetArray()`, `VecHIPReplaceArray()`
2890: @*/
2891: PetscErrorCode VecHIPPlaceArray(Vec vin, const PetscScalar a[])
2892: {
2894: #if defined(PETSC_HAVE_HIP)
2895: VecHIPCopyToGPU(vin);
2897: ((Vec_Seq *)vin->data)->unplacedarray = (PetscScalar *)((Vec_HIP *)vin->spptr)->GPUarray; /* save previous GPU array so reset can bring it back */
2898: ((Vec_HIP *)vin->spptr)->GPUarray = (PetscScalar *)a;
2899: vin->offloadmask = PETSC_OFFLOAD_GPU;
2900: #endif
2901: PetscObjectStateIncrease((PetscObject)vin);
2902: return 0;
2903: }
2905: /*@C
2906: VecHIPReplaceArray - Allows one to replace the GPU array in a vector
2907: with a GPU array provided by the user. This is useful to avoid copying
2908: a GPU array into a vector.
2910: Not Collective
2912: Input Parameters:
2913: + vec - the vector
2914: - array - the GPU array
2916: Notes:
2917: This permanently replaces the GPU array and frees the memory associated
2918: with the old GPU array.
2920: The memory passed in CANNOT be freed by the user. It will be freed
2921: when the vector is destroyed.
2923: Not supported from Fortran
2925: Level: developer
2927: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecPlaceArray()`, `VecResetArray()`, `VecHIPResetArray()`, `VecHIPPlaceArray()`, `VecReplaceArray()`
2929: @*/
2930: PetscErrorCode VecHIPReplaceArray(Vec vin, const PetscScalar a[])
2931: {
2933: #if defined(PETSC_HAVE_HIP)
2934: hipFree(((Vec_HIP *)vin->spptr)->GPUarray);
2935: ((Vec_HIP *)vin->spptr)->GPUarray = (PetscScalar *)a;
2936: vin->offloadmask = PETSC_OFFLOAD_GPU;
2937: #endif
2938: PetscObjectStateIncrease((PetscObject)vin);
2939: return 0;
2940: }
2942: /*@C
2943: VecHIPResetArray - Resets a vector to use its default memory. Call this
2944: after the use of VecHIPPlaceArray().
2946: Not Collective
2948: Input Parameters:
2949: . vec - the vector
2951: Level: developer
2953: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecReplaceArray()`, `VecPlaceArray()`, `VecResetArray()`, `VecHIPPlaceArray()`, `VecHIPReplaceArray()`
2955: @*/
2956: PetscErrorCode VecHIPResetArray(Vec vin)
2957: {
2959: #if defined(PETSC_HAVE_HIP)
2960: VecHIPCopyToGPU(vin);
2961: ((Vec_HIP *)vin->spptr)->GPUarray = (PetscScalar *)((Vec_Seq *)vin->data)->unplacedarray;
2962: ((Vec_Seq *)vin->data)->unplacedarray = 0;
2963: vin->offloadmask = PETSC_OFFLOAD_GPU;
2964: #endif
2965: PetscObjectStateIncrease((PetscObject)vin);
2966: return 0;
2967: }
2969: /*MC
2970: VecDuplicateVecsF90 - Creates several vectors of the same type as an existing vector
2971: and makes them accessible via a Fortran90 pointer.
2973: Synopsis:
2974: VecDuplicateVecsF90(Vec x,PetscInt n,{Vec, pointer :: y(:)},integer ierr)
2976: Collective on Vec
2978: Input Parameters:
2979: + x - a vector to mimic
2980: - n - the number of vectors to obtain
2982: Output Parameters:
2983: + y - Fortran90 pointer to the array of vectors
2984: - ierr - error code
2986: Example of Usage:
2987: .vb
2988: #include <petsc/finclude/petscvec.h>
2989: use petscvec
2991: Vec x
2992: Vec, pointer :: y(:)
2993: ....
2994: call VecDuplicateVecsF90(x,2,y,ierr)
2995: call VecSet(y(2),alpha,ierr)
2996: call VecSet(y(2),alpha,ierr)
2997: ....
2998: call VecDestroyVecsF90(2,y,ierr)
2999: .ve
3001: Notes:
3002: Not yet supported for all F90 compilers
3004: Use VecDestroyVecsF90() to free the space.
3006: Level: beginner
3008: .seealso: `VecDestroyVecsF90()`, `VecDuplicateVecs()`
3010: M*/
3012: /*MC
3013: VecRestoreArrayF90 - Restores a vector to a usable state after a call to
3014: VecGetArrayF90().
3016: Synopsis:
3017: VecRestoreArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
3019: Logically Collective on Vec
3021: Input Parameters:
3022: + x - vector
3023: - xx_v - the Fortran90 pointer to the array
3025: Output Parameter:
3026: . ierr - error code
3028: Example of Usage:
3029: .vb
3030: #include <petsc/finclude/petscvec.h>
3031: use petscvec
3033: PetscScalar, pointer :: xx_v(:)
3034: ....
3035: call VecGetArrayF90(x,xx_v,ierr)
3036: xx_v(3) = a
3037: call VecRestoreArrayF90(x,xx_v,ierr)
3038: .ve
3040: Level: beginner
3042: .seealso: `VecGetArrayF90()`, `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrayReadF90()`
3044: M*/
3046: /*MC
3047: VecDestroyVecsF90 - Frees a block of vectors obtained with VecDuplicateVecsF90().
3049: Synopsis:
3050: VecDestroyVecsF90(PetscInt n,{Vec, pointer :: x(:)},PetscErrorCode ierr)
3052: Collective on Vec
3054: Input Parameters:
3055: + n - the number of vectors previously obtained
3056: - x - pointer to array of vector pointers
3058: Output Parameter:
3059: . ierr - error code
3061: Notes:
3062: Not yet supported for all F90 compilers
3064: Level: beginner
3066: .seealso: `VecDestroyVecs()`, `VecDuplicateVecsF90()`
3068: M*/
3070: /*MC
3071: VecGetArrayF90 - Accesses a vector array from Fortran90. For default PETSc
3072: vectors, VecGetArrayF90() returns a pointer to the local data array. Otherwise,
3073: this routine is implementation dependent. You MUST call VecRestoreArrayF90()
3074: when you no longer need access to the array.
3076: Synopsis:
3077: VecGetArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
3079: Logically Collective on Vec
3081: Input Parameter:
3082: . x - vector
3084: Output Parameters:
3085: + xx_v - the Fortran90 pointer to the array
3086: - ierr - error code
3088: Example of Usage:
3089: .vb
3090: #include <petsc/finclude/petscvec.h>
3091: use petscvec
3093: PetscScalar, pointer :: xx_v(:)
3094: ....
3095: call VecGetArrayF90(x,xx_v,ierr)
3096: xx_v(3) = a
3097: call VecRestoreArrayF90(x,xx_v,ierr)
3098: .ve
3100: If you ONLY intend to read entries from the array and not change any entries you should use VecGetArrayReadF90().
3102: Level: beginner
3104: .seealso: `VecRestoreArrayF90()`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayReadF90()`
3106: M*/
3108: /*MC
3109: VecGetArrayReadF90 - Accesses a read only array from Fortran90. For default PETSc
3110: vectors, VecGetArrayF90() returns a pointer to the local data array. Otherwise,
3111: this routine is implementation dependent. You MUST call VecRestoreArrayReadF90()
3112: when you no longer need access to the array.
3114: Synopsis:
3115: VecGetArrayReadF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
3117: Logically Collective on Vec
3119: Input Parameter:
3120: . x - vector
3122: Output Parameters:
3123: + xx_v - the Fortran90 pointer to the array
3124: - ierr - error code
3126: Example of Usage:
3127: .vb
3128: #include <petsc/finclude/petscvec.h>
3129: use petscvec
3131: PetscScalar, pointer :: xx_v(:)
3132: ....
3133: call VecGetArrayReadF90(x,xx_v,ierr)
3134: a = xx_v(3)
3135: call VecRestoreArrayReadF90(x,xx_v,ierr)
3136: .ve
3138: If you intend to write entries into the array you must use VecGetArrayF90().
3140: Level: beginner
3142: .seealso: `VecRestoreArrayReadF90()`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayRead()`, `VecRestoreArrayRead()`, `VecGetArrayF90()`
3144: M*/
3146: /*MC
3147: VecRestoreArrayReadF90 - Restores a readonly vector to a usable state after a call to
3148: VecGetArrayReadF90().
3150: Synopsis:
3151: VecRestoreArrayReadF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
3153: Logically Collective on Vec
3155: Input Parameters:
3156: + x - vector
3157: - xx_v - the Fortran90 pointer to the array
3159: Output Parameter:
3160: . ierr - error code
3162: Example of Usage:
3163: .vb
3164: #include <petsc/finclude/petscvec.h>
3165: use petscvec
3167: PetscScalar, pointer :: xx_v(:)
3168: ....
3169: call VecGetArrayReadF90(x,xx_v,ierr)
3170: a = xx_v(3)
3171: call VecRestoreArrayReadF90(x,xx_v,ierr)
3172: .ve
3174: Level: beginner
3176: .seealso: `VecGetArrayReadF90()`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayRead()`, `VecRestoreArrayRead()`, `VecRestoreArrayF90()`
3178: M*/
3180: /*@C
3181: VecGetArray2d - Returns a pointer to a 2d contiguous array that contains this
3182: processor's portion of the vector data. You MUST call VecRestoreArray2d()
3183: when you no longer need access to the array.
3185: Logically Collective
3187: Input Parameters:
3188: + x - the vector
3189: . m - first dimension of two dimensional array
3190: . n - second dimension of two dimensional array
3191: . mstart - first index you will use in first coordinate direction (often 0)
3192: - nstart - first index in the second coordinate direction (often 0)
3194: Output Parameter:
3195: . a - location to put pointer to the array
3197: Level: developer
3199: Notes:
3200: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
3201: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3202: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3203: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
3205: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3207: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
3208: `VecRestoreArray2d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
3209: `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3210: @*/
3211: PetscErrorCode VecGetArray2d(Vec x, PetscInt m, PetscInt n, PetscInt mstart, PetscInt nstart, PetscScalar **a[])
3212: {
3213: PetscInt i, N;
3214: PetscScalar *aa;
3219: VecGetLocalSize(x, &N);
3221: VecGetArray(x, &aa);
3223: PetscMalloc1(m, a);
3224: for (i = 0; i < m; i++) (*a)[i] = aa + i * n - nstart;
3225: *a -= mstart;
3226: return 0;
3227: }
3229: /*@C
3230: VecGetArray2dWrite - Returns a pointer to a 2d contiguous array that will contain this
3231: processor's portion of the vector data. You MUST call VecRestoreArray2dWrite()
3232: when you no longer need access to the array.
3234: Logically Collective
3236: Input Parameters:
3237: + x - the vector
3238: . m - first dimension of two dimensional array
3239: . n - second dimension of two dimensional array
3240: . mstart - first index you will use in first coordinate direction (often 0)
3241: - nstart - first index in the second coordinate direction (often 0)
3243: Output Parameter:
3244: . a - location to put pointer to the array
3246: Level: developer
3248: Notes:
3249: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
3250: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3251: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3252: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
3254: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3256: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
3257: `VecRestoreArray2d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
3258: `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3259: @*/
3260: PetscErrorCode VecGetArray2dWrite(Vec x, PetscInt m, PetscInt n, PetscInt mstart, PetscInt nstart, PetscScalar **a[])
3261: {
3262: PetscInt i, N;
3263: PetscScalar *aa;
3268: VecGetLocalSize(x, &N);
3270: VecGetArrayWrite(x, &aa);
3272: PetscMalloc1(m, a);
3273: for (i = 0; i < m; i++) (*a)[i] = aa + i * n - nstart;
3274: *a -= mstart;
3275: return 0;
3276: }
3278: /*@C
3279: VecRestoreArray2d - Restores a vector after VecGetArray2d() has been called.
3281: Logically Collective
3283: Input Parameters:
3284: + x - the vector
3285: . m - first dimension of two dimensional array
3286: . n - second dimension of the two dimensional array
3287: . mstart - first index you will use in first coordinate direction (often 0)
3288: . nstart - first index in the second coordinate direction (often 0)
3289: - a - location of pointer to array obtained from VecGetArray2d()
3291: Level: developer
3293: Notes:
3294: For regular PETSc vectors this routine does not involve any copies. For
3295: any special vectors that do not store local vector data in a contiguous
3296: array, this routine will copy the data back into the underlying
3297: vector data structure from the array obtained with VecGetArray().
3299: This routine actually zeros out the a pointer.
3301: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
3302: `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
3303: `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3304: @*/
3305: PetscErrorCode VecRestoreArray2d(Vec x, PetscInt m, PetscInt n, PetscInt mstart, PetscInt nstart, PetscScalar **a[])
3306: {
3307: void *dummy;
3312: dummy = (void *)(*a + mstart);
3313: PetscFree(dummy);
3314: VecRestoreArray(x, NULL);
3315: return 0;
3316: }
3318: /*@C
3319: VecRestoreArray2dWrite - Restores a vector after VecGetArray2dWrite() has been called.
3321: Logically Collective
3323: Input Parameters:
3324: + x - the vector
3325: . m - first dimension of two dimensional array
3326: . n - second dimension of the two dimensional array
3327: . mstart - first index you will use in first coordinate direction (often 0)
3328: . nstart - first index in the second coordinate direction (often 0)
3329: - a - location of pointer to array obtained from VecGetArray2d()
3331: Level: developer
3333: Notes:
3334: For regular PETSc vectors this routine does not involve any copies. For
3335: any special vectors that do not store local vector data in a contiguous
3336: array, this routine will copy the data back into the underlying
3337: vector data structure from the array obtained with VecGetArray().
3339: This routine actually zeros out the a pointer.
3341: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
3342: `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
3343: `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3344: @*/
3345: PetscErrorCode VecRestoreArray2dWrite(Vec x, PetscInt m, PetscInt n, PetscInt mstart, PetscInt nstart, PetscScalar **a[])
3346: {
3347: void *dummy;
3352: dummy = (void *)(*a + mstart);
3353: PetscFree(dummy);
3354: VecRestoreArrayWrite(x, NULL);
3355: return 0;
3356: }
3358: /*@C
3359: VecGetArray1d - Returns a pointer to a 1d contiguous array that contains this
3360: processor's portion of the vector data. You MUST call VecRestoreArray1d()
3361: when you no longer need access to the array.
3363: Logically Collective
3365: Input Parameters:
3366: + x - the vector
3367: . m - first dimension of two dimensional array
3368: - mstart - first index you will use in first coordinate direction (often 0)
3370: Output Parameter:
3371: . a - location to put pointer to the array
3373: Level: developer
3375: Notes:
3376: For a vector obtained from DMCreateLocalVector() mstart are likely
3377: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3378: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
3380: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3382: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
3383: `VecRestoreArray2d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
3384: `VecGetArray2d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3385: @*/
3386: PetscErrorCode VecGetArray1d(Vec x, PetscInt m, PetscInt mstart, PetscScalar *a[])
3387: {
3388: PetscInt N;
3393: VecGetLocalSize(x, &N);
3395: VecGetArray(x, a);
3396: *a -= mstart;
3397: return 0;
3398: }
3400: /*@C
3401: VecGetArray1dWrite - Returns a pointer to a 1d contiguous array that will contain this
3402: processor's portion of the vector data. You MUST call VecRestoreArray1dWrite()
3403: when you no longer need access to the array.
3405: Logically Collective
3407: Input Parameters:
3408: + x - the vector
3409: . m - first dimension of two dimensional array
3410: - mstart - first index you will use in first coordinate direction (often 0)
3412: Output Parameter:
3413: . a - location to put pointer to the array
3415: Level: developer
3417: Notes:
3418: For a vector obtained from DMCreateLocalVector() mstart are likely
3419: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3420: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
3422: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3424: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
3425: `VecRestoreArray2d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
3426: `VecGetArray2d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3427: @*/
3428: PetscErrorCode VecGetArray1dWrite(Vec x, PetscInt m, PetscInt mstart, PetscScalar *a[])
3429: {
3430: PetscInt N;
3435: VecGetLocalSize(x, &N);
3437: VecGetArrayWrite(x, a);
3438: *a -= mstart;
3439: return 0;
3440: }
3442: /*@C
3443: VecRestoreArray1d - Restores a vector after VecGetArray1d() has been called.
3445: Logically Collective
3447: Input Parameters:
3448: + x - the vector
3449: . m - first dimension of two dimensional array
3450: . mstart - first index you will use in first coordinate direction (often 0)
3451: - a - location of pointer to array obtained from VecGetArray21()
3453: Level: developer
3455: Notes:
3456: For regular PETSc vectors this routine does not involve any copies. For
3457: any special vectors that do not store local vector data in a contiguous
3458: array, this routine will copy the data back into the underlying
3459: vector data structure from the array obtained with VecGetArray1d().
3461: This routine actually zeros out the a pointer.
3463: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
3464: `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
3465: `VecGetArray1d()`, `VecRestoreArray2d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3466: @*/
3467: PetscErrorCode VecRestoreArray1d(Vec x, PetscInt m, PetscInt mstart, PetscScalar *a[])
3468: {
3471: VecRestoreArray(x, NULL);
3472: return 0;
3473: }
3475: /*@C
3476: VecRestoreArray1dWrite - Restores a vector after VecGetArray1dWrite() has been called.
3478: Logically Collective
3480: Input Parameters:
3481: + x - the vector
3482: . m - first dimension of two dimensional array
3483: . mstart - first index you will use in first coordinate direction (often 0)
3484: - a - location of pointer to array obtained from VecGetArray21()
3486: Level: developer
3488: Notes:
3489: For regular PETSc vectors this routine does not involve any copies. For
3490: any special vectors that do not store local vector data in a contiguous
3491: array, this routine will copy the data back into the underlying
3492: vector data structure from the array obtained with VecGetArray1d().
3494: This routine actually zeros out the a pointer.
3496: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
3497: `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
3498: `VecGetArray1d()`, `VecRestoreArray2d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3499: @*/
3500: PetscErrorCode VecRestoreArray1dWrite(Vec x, PetscInt m, PetscInt mstart, PetscScalar *a[])
3501: {
3504: VecRestoreArrayWrite(x, NULL);
3505: return 0;
3506: }
3508: /*@C
3509: VecGetArray3d - Returns a pointer to a 3d contiguous array that contains this
3510: processor's portion of the vector data. You MUST call VecRestoreArray3d()
3511: when you no longer need access to the array.
3513: Logically Collective
3515: Input Parameters:
3516: + x - the vector
3517: . m - first dimension of three dimensional array
3518: . n - second dimension of three dimensional array
3519: . p - third dimension of three dimensional array
3520: . mstart - first index you will use in first coordinate direction (often 0)
3521: . nstart - first index in the second coordinate direction (often 0)
3522: - pstart - first index in the third coordinate direction (often 0)
3524: Output Parameter:
3525: . a - location to put pointer to the array
3527: Level: developer
3529: Notes:
3530: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3531: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3532: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3533: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3535: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3537: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
3538: `VecRestoreArray2d()`, `DMDAVecGetarray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
3539: `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3540: @*/
3541: PetscErrorCode VecGetArray3d(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscScalar ***a[])
3542: {
3543: PetscInt i, N, j;
3544: PetscScalar *aa, **b;
3549: VecGetLocalSize(x, &N);
3551: VecGetArray(x, &aa);
3553: PetscMalloc(m * sizeof(PetscScalar **) + m * n * sizeof(PetscScalar *), a);
3554: b = (PetscScalar **)((*a) + m);
3555: for (i = 0; i < m; i++) (*a)[i] = b + i * n - nstart;
3556: for (i = 0; i < m; i++)
3557: for (j = 0; j < n; j++) b[i * n + j] = aa + i * n * p + j * p - pstart;
3558: *a -= mstart;
3559: return 0;
3560: }
3562: /*@C
3563: VecGetArray3dWrite - Returns a pointer to a 3d contiguous array that will contain this
3564: processor's portion of the vector data. You MUST call VecRestoreArray3dWrite()
3565: when you no longer need access to the array.
3567: Logically Collective
3569: Input Parameters:
3570: + x - the vector
3571: . m - first dimension of three dimensional array
3572: . n - second dimension of three dimensional array
3573: . p - third dimension of three dimensional array
3574: . mstart - first index you will use in first coordinate direction (often 0)
3575: . nstart - first index in the second coordinate direction (often 0)
3576: - pstart - first index in the third coordinate direction (often 0)
3578: Output Parameter:
3579: . a - location to put pointer to the array
3581: Level: developer
3583: Notes:
3584: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3585: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3586: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3587: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3589: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3591: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
3592: `VecRestoreArray2d()`, `DMDAVecGetarray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
3593: `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3594: @*/
3595: PetscErrorCode VecGetArray3dWrite(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscScalar ***a[])
3596: {
3597: PetscInt i, N, j;
3598: PetscScalar *aa, **b;
3603: VecGetLocalSize(x, &N);
3605: VecGetArrayWrite(x, &aa);
3607: PetscMalloc(m * sizeof(PetscScalar **) + m * n * sizeof(PetscScalar *), a);
3608: b = (PetscScalar **)((*a) + m);
3609: for (i = 0; i < m; i++) (*a)[i] = b + i * n - nstart;
3610: for (i = 0; i < m; i++)
3611: for (j = 0; j < n; j++) b[i * n + j] = aa + i * n * p + j * p - pstart;
3613: *a -= mstart;
3614: return 0;
3615: }
3617: /*@C
3618: VecRestoreArray3d - Restores a vector after VecGetArray3d() has been called.
3620: Logically Collective
3622: Input Parameters:
3623: + x - the vector
3624: . m - first dimension of three dimensional array
3625: . n - second dimension of the three dimensional array
3626: . p - third dimension of the three dimensional array
3627: . mstart - first index you will use in first coordinate direction (often 0)
3628: . nstart - first index in the second coordinate direction (often 0)
3629: . pstart - first index in the third coordinate direction (often 0)
3630: - a - location of pointer to array obtained from VecGetArray3d()
3632: Level: developer
3634: Notes:
3635: For regular PETSc vectors this routine does not involve any copies. For
3636: any special vectors that do not store local vector data in a contiguous
3637: array, this routine will copy the data back into the underlying
3638: vector data structure from the array obtained with VecGetArray().
3640: This routine actually zeros out the a pointer.
3642: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
3643: `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
3644: `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`, `VecGet`
3645: @*/
3646: PetscErrorCode VecRestoreArray3d(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscScalar ***a[])
3647: {
3648: void *dummy;
3653: dummy = (void *)(*a + mstart);
3654: PetscFree(dummy);
3655: VecRestoreArray(x, NULL);
3656: return 0;
3657: }
3659: /*@C
3660: VecRestoreArray3dWrite - Restores a vector after VecGetArray3dWrite() has been called.
3662: Logically Collective
3664: Input Parameters:
3665: + x - the vector
3666: . m - first dimension of three dimensional array
3667: . n - second dimension of the three dimensional array
3668: . p - third dimension of the three dimensional array
3669: . mstart - first index you will use in first coordinate direction (often 0)
3670: . nstart - first index in the second coordinate direction (often 0)
3671: . pstart - first index in the third coordinate direction (often 0)
3672: - a - location of pointer to array obtained from VecGetArray3d()
3674: Level: developer
3676: Notes:
3677: For regular PETSc vectors this routine does not involve any copies. For
3678: any special vectors that do not store local vector data in a contiguous
3679: array, this routine will copy the data back into the underlying
3680: vector data structure from the array obtained with VecGetArray().
3682: This routine actually zeros out the a pointer.
3684: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
3685: `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
3686: `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`, `VecGet`
3687: @*/
3688: PetscErrorCode VecRestoreArray3dWrite(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscScalar ***a[])
3689: {
3690: void *dummy;
3695: dummy = (void *)(*a + mstart);
3696: PetscFree(dummy);
3697: VecRestoreArrayWrite(x, NULL);
3698: return 0;
3699: }
3701: /*@C
3702: VecGetArray4d - Returns a pointer to a 4d contiguous array that contains this
3703: processor's portion of the vector data. You MUST call VecRestoreArray4d()
3704: when you no longer need access to the array.
3706: Logically Collective
3708: Input Parameters:
3709: + x - the vector
3710: . m - first dimension of four dimensional array
3711: . n - second dimension of four dimensional array
3712: . p - third dimension of four dimensional array
3713: . q - fourth dimension of four dimensional array
3714: . mstart - first index you will use in first coordinate direction (often 0)
3715: . nstart - first index in the second coordinate direction (often 0)
3716: . pstart - first index in the third coordinate direction (often 0)
3717: - qstart - first index in the fourth coordinate direction (often 0)
3719: Output Parameter:
3720: . a - location to put pointer to the array
3722: Level: beginner
3724: Notes:
3725: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3726: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3727: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3728: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3730: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3732: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
3733: `VecRestoreArray2d()`, `DMDAVecGetarray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
3734: `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3735: @*/
3736: PetscErrorCode VecGetArray4d(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt q, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscInt qstart, PetscScalar ****a[])
3737: {
3738: PetscInt i, N, j, k;
3739: PetscScalar *aa, ***b, **c;
3744: VecGetLocalSize(x, &N);
3746: VecGetArray(x, &aa);
3748: PetscMalloc(m * sizeof(PetscScalar ***) + m * n * sizeof(PetscScalar **) + m * n * p * sizeof(PetscScalar *), a);
3749: b = (PetscScalar ***)((*a) + m);
3750: c = (PetscScalar **)(b + m * n);
3751: for (i = 0; i < m; i++) (*a)[i] = b + i * n - nstart;
3752: for (i = 0; i < m; i++)
3753: for (j = 0; j < n; j++) b[i * n + j] = c + i * n * p + j * p - pstart;
3754: for (i = 0; i < m; i++)
3755: for (j = 0; j < n; j++)
3756: for (k = 0; k < p; k++) c[i * n * p + j * p + k] = aa + i * n * p * q + j * p * q + k * q - qstart;
3757: *a -= mstart;
3758: return 0;
3759: }
3761: /*@C
3762: VecGetArray4dWrite - Returns a pointer to a 4d contiguous array that will contain this
3763: processor's portion of the vector data. You MUST call VecRestoreArray4dWrite()
3764: when you no longer need access to the array.
3766: Logically Collective
3768: Input Parameters:
3769: + x - the vector
3770: . m - first dimension of four dimensional array
3771: . n - second dimension of four dimensional array
3772: . p - third dimension of four dimensional array
3773: . q - fourth dimension of four dimensional array
3774: . mstart - first index you will use in first coordinate direction (often 0)
3775: . nstart - first index in the second coordinate direction (often 0)
3776: . pstart - first index in the third coordinate direction (often 0)
3777: - qstart - first index in the fourth coordinate direction (often 0)
3779: Output Parameter:
3780: . a - location to put pointer to the array
3782: Level: beginner
3784: Notes:
3785: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3786: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3787: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3788: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3790: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3792: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
3793: `VecRestoreArray2d()`, `DMDAVecGetarray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
3794: `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3795: @*/
3796: PetscErrorCode VecGetArray4dWrite(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt q, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscInt qstart, PetscScalar ****a[])
3797: {
3798: PetscInt i, N, j, k;
3799: PetscScalar *aa, ***b, **c;
3804: VecGetLocalSize(x, &N);
3806: VecGetArrayWrite(x, &aa);
3808: PetscMalloc(m * sizeof(PetscScalar ***) + m * n * sizeof(PetscScalar **) + m * n * p * sizeof(PetscScalar *), a);
3809: b = (PetscScalar ***)((*a) + m);
3810: c = (PetscScalar **)(b + m * n);
3811: for (i = 0; i < m; i++) (*a)[i] = b + i * n - nstart;
3812: for (i = 0; i < m; i++)
3813: for (j = 0; j < n; j++) b[i * n + j] = c + i * n * p + j * p - pstart;
3814: for (i = 0; i < m; i++)
3815: for (j = 0; j < n; j++)
3816: for (k = 0; k < p; k++) c[i * n * p + j * p + k] = aa + i * n * p * q + j * p * q + k * q - qstart;
3817: *a -= mstart;
3818: return 0;
3819: }
3821: /*@C
3822: VecRestoreArray4d - Restores a vector after VecGetArray3d() has been called.
3824: Logically Collective
3826: Input Parameters:
3827: + x - the vector
3828: . m - first dimension of four dimensional array
3829: . n - second dimension of the four dimensional array
3830: . p - third dimension of the four dimensional array
3831: . q - fourth dimension of the four dimensional array
3832: . mstart - first index you will use in first coordinate direction (often 0)
3833: . nstart - first index in the second coordinate direction (often 0)
3834: . pstart - first index in the third coordinate direction (often 0)
3835: . qstart - first index in the fourth coordinate direction (often 0)
3836: - a - location of pointer to array obtained from VecGetArray4d()
3838: Level: beginner
3840: Notes:
3841: For regular PETSc vectors this routine does not involve any copies. For
3842: any special vectors that do not store local vector data in a contiguous
3843: array, this routine will copy the data back into the underlying
3844: vector data structure from the array obtained with VecGetArray().
3846: This routine actually zeros out the a pointer.
3848: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
3849: `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
3850: `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`, `VecGet`
3851: @*/
3852: PetscErrorCode VecRestoreArray4d(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt q, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscInt qstart, PetscScalar ****a[])
3853: {
3854: void *dummy;
3859: dummy = (void *)(*a + mstart);
3860: PetscFree(dummy);
3861: VecRestoreArray(x, NULL);
3862: return 0;
3863: }
3865: /*@C
3866: VecRestoreArray4dWrite - Restores a vector after VecGetArray3dWrite() has been called.
3868: Logically Collective
3870: Input Parameters:
3871: + x - the vector
3872: . m - first dimension of four dimensional array
3873: . n - second dimension of the four dimensional array
3874: . p - third dimension of the four dimensional array
3875: . q - fourth dimension of the four dimensional array
3876: . mstart - first index you will use in first coordinate direction (often 0)
3877: . nstart - first index in the second coordinate direction (often 0)
3878: . pstart - first index in the third coordinate direction (often 0)
3879: . qstart - first index in the fourth coordinate direction (often 0)
3880: - a - location of pointer to array obtained from VecGetArray4d()
3882: Level: beginner
3884: Notes:
3885: For regular PETSc vectors this routine does not involve any copies. For
3886: any special vectors that do not store local vector data in a contiguous
3887: array, this routine will copy the data back into the underlying
3888: vector data structure from the array obtained with VecGetArray().
3890: This routine actually zeros out the a pointer.
3892: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
3893: `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
3894: `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`, `VecGet`
3895: @*/
3896: PetscErrorCode VecRestoreArray4dWrite(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt q, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscInt qstart, PetscScalar ****a[])
3897: {
3898: void *dummy;
3903: dummy = (void *)(*a + mstart);
3904: PetscFree(dummy);
3905: VecRestoreArrayWrite(x, NULL);
3906: return 0;
3907: }
3909: /*@C
3910: VecGetArray2dRead - Returns a pointer to a 2d contiguous array that contains this
3911: processor's portion of the vector data. You MUST call VecRestoreArray2dRead()
3912: when you no longer need access to the array.
3914: Logically Collective
3916: Input Parameters:
3917: + x - the vector
3918: . m - first dimension of two dimensional array
3919: . n - second dimension of two dimensional array
3920: . mstart - first index you will use in first coordinate direction (often 0)
3921: - nstart - first index in the second 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 and nstart are likely
3930: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3931: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3932: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
3934: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3936: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
3937: `VecRestoreArray2d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
3938: `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3939: @*/
3940: PetscErrorCode VecGetArray2dRead(Vec x, PetscInt m, PetscInt n, PetscInt mstart, PetscInt nstart, PetscScalar **a[])
3941: {
3942: PetscInt i, N;
3943: const PetscScalar *aa;
3948: VecGetLocalSize(x, &N);
3950: VecGetArrayRead(x, &aa);
3952: PetscMalloc1(m, a);
3953: for (i = 0; i < m; i++) (*a)[i] = (PetscScalar *)aa + i * n - nstart;
3954: *a -= mstart;
3955: return 0;
3956: }
3958: /*@C
3959: VecRestoreArray2dRead - Restores a vector after VecGetArray2dRead() has been called.
3961: Logically Collective
3963: Input Parameters:
3964: + x - the vector
3965: . m - first dimension of two dimensional array
3966: . n - second dimension of the two dimensional array
3967: . mstart - first index you will use in first coordinate direction (often 0)
3968: . nstart - first index in the second coordinate direction (often 0)
3969: - a - location of pointer to array obtained from VecGetArray2d()
3971: Level: developer
3973: Notes:
3974: For regular PETSc vectors this routine does not involve any copies. For
3975: any special vectors that do not store local vector data in a contiguous
3976: array, this routine will copy the data back into the underlying
3977: vector data structure from the array obtained with VecGetArray().
3979: This routine actually zeros out the a pointer.
3981: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
3982: `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
3983: `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3984: @*/
3985: PetscErrorCode VecRestoreArray2dRead(Vec x, PetscInt m, PetscInt n, PetscInt mstart, PetscInt nstart, PetscScalar **a[])
3986: {
3987: void *dummy;
3992: dummy = (void *)(*a + mstart);
3993: PetscFree(dummy);
3994: VecRestoreArrayRead(x, NULL);
3995: return 0;
3996: }
3998: /*@C
3999: VecGetArray1dRead - Returns a pointer to a 1d contiguous array that contains this
4000: processor's portion of the vector data. You MUST call VecRestoreArray1dRead()
4001: when you no longer need access to the array.
4003: Logically Collective
4005: Input Parameters:
4006: + x - the vector
4007: . m - first dimension of two dimensional array
4008: - mstart - first index you will use in first coordinate direction (often 0)
4010: Output Parameter:
4011: . a - location to put pointer to the array
4013: Level: developer
4015: Notes:
4016: For a vector obtained from DMCreateLocalVector() mstart are likely
4017: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
4018: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
4020: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
4022: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
4023: `VecRestoreArray2d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
4024: `VecGetArray2d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
4025: @*/
4026: PetscErrorCode VecGetArray1dRead(Vec x, PetscInt m, PetscInt mstart, PetscScalar *a[])
4027: {
4028: PetscInt N;
4033: VecGetLocalSize(x, &N);
4035: VecGetArrayRead(x, (const PetscScalar **)a);
4036: *a -= mstart;
4037: return 0;
4038: }
4040: /*@C
4041: VecRestoreArray1dRead - Restores a vector after VecGetArray1dRead() has been called.
4043: Logically Collective
4045: Input Parameters:
4046: + x - the vector
4047: . m - first dimension of two dimensional array
4048: . mstart - first index you will use in first coordinate direction (often 0)
4049: - a - location of pointer to array obtained from VecGetArray21()
4051: Level: developer
4053: Notes:
4054: For regular PETSc vectors this routine does not involve any copies. For
4055: any special vectors that do not store local vector data in a contiguous
4056: array, this routine will copy the data back into the underlying
4057: vector data structure from the array obtained with VecGetArray1dRead().
4059: This routine actually zeros out the a pointer.
4061: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
4062: `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
4063: `VecGetArray1d()`, `VecRestoreArray2d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
4064: @*/
4065: PetscErrorCode VecRestoreArray1dRead(Vec x, PetscInt m, PetscInt mstart, PetscScalar *a[])
4066: {
4069: VecRestoreArrayRead(x, NULL);
4070: return 0;
4071: }
4073: /*@C
4074: VecGetArray3dRead - Returns a pointer to a 3d contiguous array that contains this
4075: processor's portion of the vector data. You MUST call VecRestoreArray3dRead()
4076: when you no longer need access to the array.
4078: Logically Collective
4080: Input Parameters:
4081: + x - the vector
4082: . m - first dimension of three dimensional array
4083: . n - second dimension of three dimensional array
4084: . p - third dimension of three dimensional array
4085: . mstart - first index you will use in first coordinate direction (often 0)
4086: . nstart - first index in the second coordinate direction (often 0)
4087: - pstart - first index in the third coordinate direction (often 0)
4089: Output Parameter:
4090: . a - location to put pointer to the array
4092: Level: developer
4094: Notes:
4095: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
4096: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
4097: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
4098: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3dRead().
4100: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
4102: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
4103: `VecRestoreArray2d()`, `DMDAVecGetarray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
4104: `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
4105: @*/
4106: PetscErrorCode VecGetArray3dRead(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscScalar ***a[])
4107: {
4108: PetscInt i, N, j;
4109: const PetscScalar *aa;
4110: PetscScalar **b;
4115: VecGetLocalSize(x, &N);
4117: VecGetArrayRead(x, &aa);
4119: PetscMalloc(m * sizeof(PetscScalar **) + m * n * sizeof(PetscScalar *), a);
4120: b = (PetscScalar **)((*a) + m);
4121: for (i = 0; i < m; i++) (*a)[i] = b + i * n - nstart;
4122: for (i = 0; i < m; i++)
4123: for (j = 0; j < n; j++) b[i * n + j] = (PetscScalar *)aa + i * n * p + j * p - pstart;
4124: *a -= mstart;
4125: return 0;
4126: }
4128: /*@C
4129: VecRestoreArray3dRead - Restores a vector after VecGetArray3dRead() has been called.
4131: Logically Collective
4133: Input Parameters:
4134: + x - the vector
4135: . m - first dimension of three dimensional array
4136: . n - second dimension of the three dimensional array
4137: . p - third dimension of the three dimensional array
4138: . mstart - first index you will use in first coordinate direction (often 0)
4139: . nstart - first index in the second coordinate direction (often 0)
4140: . pstart - first index in the third coordinate direction (often 0)
4141: - a - location of pointer to array obtained from VecGetArray3dRead()
4143: Level: developer
4145: Notes:
4146: For regular PETSc vectors this routine does not involve any copies. For
4147: any special vectors that do not store local vector data in a contiguous
4148: array, this routine will copy the data back into the underlying
4149: vector data structure from the array obtained with VecGetArray().
4151: This routine actually zeros out the a pointer.
4153: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
4154: `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
4155: `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`, `VecGet`
4156: @*/
4157: PetscErrorCode VecRestoreArray3dRead(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscScalar ***a[])
4158: {
4159: void *dummy;
4164: dummy = (void *)(*a + mstart);
4165: PetscFree(dummy);
4166: VecRestoreArrayRead(x, NULL);
4167: return 0;
4168: }
4170: /*@C
4171: VecGetArray4dRead - Returns a pointer to a 4d contiguous array that contains this
4172: processor's portion of the vector data. You MUST call VecRestoreArray4dRead()
4173: when you no longer need access to the array.
4175: Logically Collective
4177: Input Parameters:
4178: + x - the vector
4179: . m - first dimension of four dimensional array
4180: . n - second dimension of four dimensional array
4181: . p - third dimension of four dimensional array
4182: . q - fourth dimension of four dimensional array
4183: . mstart - first index you will use in first coordinate direction (often 0)
4184: . nstart - first index in the second coordinate direction (often 0)
4185: . pstart - first index in the third coordinate direction (often 0)
4186: - qstart - first index in the fourth coordinate direction (often 0)
4188: Output Parameter:
4189: . a - location to put pointer to the array
4191: Level: beginner
4193: Notes:
4194: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
4195: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
4196: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
4197: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
4199: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
4201: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
4202: `VecRestoreArray2d()`, `DMDAVecGetarray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
4203: `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
4204: @*/
4205: PetscErrorCode VecGetArray4dRead(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt q, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscInt qstart, PetscScalar ****a[])
4206: {
4207: PetscInt i, N, j, k;
4208: const PetscScalar *aa;
4209: PetscScalar ***b, **c;
4214: VecGetLocalSize(x, &N);
4216: VecGetArrayRead(x, &aa);
4218: PetscMalloc(m * sizeof(PetscScalar ***) + m * n * sizeof(PetscScalar **) + m * n * p * sizeof(PetscScalar *), a);
4219: b = (PetscScalar ***)((*a) + m);
4220: c = (PetscScalar **)(b + m * n);
4221: for (i = 0; i < m; i++) (*a)[i] = b + i * n - nstart;
4222: for (i = 0; i < m; i++)
4223: for (j = 0; j < n; j++) b[i * n + j] = c + i * n * p + j * p - pstart;
4224: for (i = 0; i < m; i++)
4225: for (j = 0; j < n; j++)
4226: for (k = 0; k < p; k++) c[i * n * p + j * p + k] = (PetscScalar *)aa + i * n * p * q + j * p * q + k * q - qstart;
4227: *a -= mstart;
4228: return 0;
4229: }
4231: /*@C
4232: VecRestoreArray4dRead - Restores a vector after VecGetArray3d() has been called.
4234: Logically Collective
4236: Input Parameters:
4237: + x - the vector
4238: . m - first dimension of four dimensional array
4239: . n - second dimension of the four dimensional array
4240: . p - third dimension of the four dimensional array
4241: . q - fourth dimension of the four dimensional array
4242: . mstart - first index you will use in first coordinate direction (often 0)
4243: . nstart - first index in the second coordinate direction (often 0)
4244: . pstart - first index in the third coordinate direction (often 0)
4245: . qstart - first index in the fourth coordinate direction (often 0)
4246: - a - location of pointer to array obtained from VecGetArray4dRead()
4248: Level: beginner
4250: Notes:
4251: For regular PETSc vectors this routine does not involve any copies. For
4252: any special vectors that do not store local vector data in a contiguous
4253: array, this routine will copy the data back into the underlying
4254: vector data structure from the array obtained with VecGetArray().
4256: This routine actually zeros out the a pointer.
4258: .seealso: `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
4259: `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
4260: `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`, `VecGet`
4261: @*/
4262: PetscErrorCode VecRestoreArray4dRead(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt q, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscInt qstart, PetscScalar ****a[])
4263: {
4264: void *dummy;
4269: dummy = (void *)(*a + mstart);
4270: PetscFree(dummy);
4271: VecRestoreArrayRead(x, NULL);
4272: return 0;
4273: }
4275: #if defined(PETSC_USE_DEBUG)
4277: /*@
4278: VecLockGet - Gets the current lock status of a vector
4280: Logically Collective on Vec
4282: Input Parameter:
4283: . x - the vector
4285: Output Parameter:
4286: . state - greater than zero indicates the vector is locked for read; less then zero indicates the vector is
4287: locked for write; equal to zero means the vector is unlocked, that is, it is free to read or write.
4289: Level: beginner
4291: .seealso: `VecRestoreArray()`, `VecGetArrayRead()`, `VecLockReadPush()`, `VecLockReadPop()`
4292: @*/
4293: PetscErrorCode VecLockGet(Vec x, PetscInt *state)
4294: {
4296: *state = x->lock;
4297: return 0;
4298: }
4300: PetscErrorCode VecLockGetLocation(Vec x, const char *file[], const char *func[], int *line)
4301: {
4306: {
4307: const int index = x->lockstack.currentsize - 1;
4310: *file = x->lockstack.file[index];
4311: *func = x->lockstack.function[index];
4312: *line = x->lockstack.line[index];
4313: }
4314: return 0;
4315: }
4317: /*@
4318: VecLockReadPush - Pushes a read-only lock on a vector to prevent it from writing
4320: Logically Collective on Vec
4322: Input Parameter:
4323: . x - the vector
4325: Notes:
4326: If this is set then calls to VecGetArray() or VecSetValues() or any other routines that change the vectors values will fail.
4328: The call can be nested, i.e., called multiple times on the same vector, but each VecLockReadPush(x) has to have one matching
4329: VecLockReadPop(x), which removes the latest read-only lock.
4331: Level: beginner
4333: .seealso: `VecRestoreArray()`, `VecGetArrayRead()`, `VecLockReadPop()`, `VecLockGet()`
4334: @*/
4335: PetscErrorCode VecLockReadPush(Vec x)
4336: {
4337: const char *file, *func;
4338: int index, line;
4342: if ((index = petscstack.currentsize - 2) == -1) {
4343: // vec was locked "outside" of petsc, either in user-land or main. the error message will
4344: // now show this function as the culprit, but it will include the stacktrace
4345: file = "unknown user-file";
4346: func = "unknown_user_function";
4347: line = 0;
4348: } else {
4350: file = petscstack.file[index];
4351: func = petscstack.function[index];
4352: line = petscstack.line[index];
4353: }
4354: PetscStackPush_Private(x->lockstack, file, func, line, petscstack.petscroutine[index], PETSC_FALSE);
4355: return 0;
4356: }
4358: /*@
4359: VecLockReadPop - Pops a read-only lock from a vector
4361: Logically Collective on Vec
4363: Input Parameter:
4364: . x - the vector
4366: Level: beginner
4368: .seealso: `VecRestoreArray()`, `VecGetArrayRead()`, `VecLockReadPush()`, `VecLockGet()`
4369: @*/
4370: PetscErrorCode VecLockReadPop(Vec x)
4371: {
4374: {
4375: const char *previous = x->lockstack.function[x->lockstack.currentsize - 1];
4377: PetscStackPop_Private(x->lockstack, previous);
4378: }
4379: return 0;
4380: }
4382: /*@C
4383: VecLockWriteSet - Lock or unlock a vector for exclusive read/write access
4385: Logically Collective on Vec
4387: Input Parameters:
4388: + x - the vector
4389: - flg - PETSC_TRUE to lock the vector for exclusive read/write access; PETSC_FALSE to unlock it.
4391: Notes:
4392: The function is usefull in split-phase computations, which usually have a begin phase and an end phase.
4393: One can call VecLockWriteSet(x,PETSC_TRUE) in the begin phase to lock a vector for exclusive
4394: access, and call VecLockWriteSet(x,PETSC_FALSE) in the end phase to unlock the vector from exclusive
4395: access. In this way, one is ensured no other operations can access the vector in between. The code may like
4397: VecGetArray(x,&xdata); // begin phase
4398: VecLockWriteSet(v,PETSC_TRUE);
4400: Other operations, which can not access x anymore (they can access xdata, of course)
4402: VecRestoreArray(x,&vdata); // end phase
4403: VecLockWriteSet(v,PETSC_FALSE);
4405: The call can not be nested on the same vector, in other words, one can not call VecLockWriteSet(x,PETSC_TRUE)
4406: again before calling VecLockWriteSet(v,PETSC_FALSE).
4408: Level: beginner
4410: .seealso: `VecRestoreArray()`, `VecGetArrayRead()`, `VecLockReadPush()`, `VecLockReadPop()`, `VecLockGet()`
4411: @*/
4412: PetscErrorCode VecLockWriteSet(Vec x, PetscBool flg)
4413: {
4415: if (flg) {
4418: x->lock = -1;
4419: } else {
4421: x->lock = 0;
4422: }
4423: return 0;
4424: }
4426: /*@
4427: VecLockPush - Pushes a read-only lock on a vector to prevent it from writing
4429: Level: deprecated
4431: .seealso: `VecLockReadPush()`
4432: @*/
4433: PetscErrorCode VecLockPush(Vec x)
4434: {
4435: VecLockReadPush(x);
4436: return 0;
4437: }
4439: /*@
4440: VecLockPop - Pops a read-only lock from a vector
4442: Level: deprecated
4444: .seealso: `VecLockReadPop()`
4445: @*/
4446: PetscErrorCode VecLockPop(Vec x)
4447: {
4448: VecLockReadPop(x);
4449: return 0;
4450: }
4452: #endif