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

 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: [](chapter_vectors), `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

 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:    Level: intermediate

101:    Notes for Users of Complex Numbers:
102:    For complex vectors, `VecDot()` computes
103: $     val = (x,y) = y^H x,
104:    where y^H denotes the conjugate transpose of y. Note that this corresponds to the usual "mathematicians" complex
105:    inner product where the SECOND argument gets the complex conjugate. Since the `BLASdot()` complex conjugates the first
106:    first argument we call the `BLASdot()` with the arguments reversed.

108:    Use `VecTDot()` for the indefinite form
109: $     val = (x,y) = y^T x,
110:    where y^T denotes the transpose of y.

112: .seealso: [](chapter_vectors), `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

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: [](chapter_vectors), `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

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:     Level: intermediate

196:    Notes:
197:       For complex numbers `NORM_1` will return the traditional 1 norm of the 2 norm of the complex numbers; that is the 1
198:       norm of the absolute values of the complex entries. In PETSc 3.6 and earlier releases it returned the 1 norm of
199:       the 1 norm of the complex entries (what is returned by the BLAS routine asum()). Both are valid norms but most
200:       people expect the former.

202:       This routine stashes the computed norm value, repeated calls before the vector entries are changed are then rapid since the
203:       precomputed value is immediately available. Certain vector operations such as VecSet() store the norms so the value is
204:       immediately available and does not need to be explicitly computed. `VecScale()` updates any stashed norm values, thus calls after `VecScale()`
205:       do not need to explicitly recompute the norm.

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: [](chapter_vectors), `Vec`, `VecDot()`, `VecTDot()`, `VecDotBegin()`, `VecDotEnd()`, `VecNormAvailable()`,
214:           `VecNormBegin()`, `VecNormEnd()`, `NormType()`
215: @*/
216: PetscErrorCode VecNorm(Vec x, NormType type, PetscReal *val)
217: {

222:   /* Cached data? */
223:   if (type != NORM_1_AND_2) {
224:     PetscBool flg;

226:     PetscObjectComposedDataGetReal((PetscObject)x, NormIds[type], *val, flg);
227:     if (flg) return 0;
228:   }

230:   VecLockReadPush(x);
231:   PetscLogEventBegin(VEC_Norm, x, 0, 0, 0);
232:   PetscUseTypeMethod(x, norm, type, val);
233:   PetscLogEventEnd(VEC_Norm, x, 0, 0, 0);
234:   VecLockReadPop(x);

236:   if (type != NORM_1_AND_2) PetscObjectComposedDataSetReal((PetscObject)x, NormIds[type], *val);
237:   return 0;
238: }

240: /*@
241:    VecNormAvailable  - Returns the vector norm if it is already known.

243:    Not Collective

245:    Input Parameters:
246: +  x - the vector
247: -  type - one of `NORM_1` (sum_i |x[i]|), `NORM_2` sqrt(sum_i (x[i])^2), `NORM_INFINITY` max_i |x[i]|.  Also available
248:           `NORM_1_AND_2`, which computes both norms and stores them
249:           in a two element array.

251:    Output Parameters:
252: +  available - `PETSC_TRUE` if the val returned is valid
253: -  val - the norm

255:    Level: intermediate

257:    Performance Issues:
258: .vb
259:     per-processor memory bandwidth
260:     interprocessor latency
261:     work load imbalance that causes certain processes to arrive much earlier than others
262: .ve

264:    Developer Note:
265:    PETSC_HAVE_SLOW_BLAS_NORM2 will cause a C (loop unrolled) version of the norm to be used, rather
266:    than the BLAS. This should probably only be used when one is using the FORTRAN BLAS routines
267:    (as opposed to vendor provided) because the FORTRAN BLAS NRM2() routine is very slow.

269: .seealso: [](chapter_vectors), `Vec`, `VecDot()`, `VecTDot()`, `VecNorm()`, `VecDotBegin()`, `VecDotEnd()`,
270:           `VecNormBegin()`, `VecNormEnd()`
271: @*/
272: PetscErrorCode VecNormAvailable(Vec x, NormType type, PetscBool *available, PetscReal *val)
273: {

279:   if (type == NORM_1_AND_2) {
280:     *available = PETSC_FALSE;
281:   } else {
282:     PetscObjectComposedDataGetReal((PetscObject)x, NormIds[type], *val, *available);
283:   }
284:   return 0;
285: }

287: /*@
288:    VecNormalize - Normalizes a vector by its 2-norm.

290:    Collective

292:    Input Parameter:
293: .  x - the vector

295:    Output Parameter:
296: .  val - the vector norm before normalization. May be `NULL` if the value is not needed.

298:    Level: intermediate

300: .seealso: [](chapter_vectors), `Vec`, `VecNorm()`
301: @*/
302: PetscErrorCode VecNormalize(Vec x, PetscReal *val)
303: {
304:   PetscReal norm;

308:   VecSetErrorIfLocked(x, 1);
310:   PetscLogEventBegin(VEC_Normalize, x, 0, 0, 0);
311:   VecNorm(x, NORM_2, &norm);
312:   if (norm == 0.0) {
313:     PetscInfo(x, "Vector of zero norm can not be normalized; Returning only the zero norm\n");
314:   } else if (norm != 1.0) {
315:     VecScale(x, 1.0 / norm);
316:   }
317:   PetscLogEventEnd(VEC_Normalize, x, 0, 0, 0);
318:   if (val) *val = norm;
319:   return 0;
320: }

322: /*@C
323:    VecMax - Determines the vector component with maximum real part and its location.

325:    Collective

327:    Input Parameter:
328: .  x - the vector

330:    Output Parameters:
331: +  p - the location of val (pass `NULL` if you don't want this)
332: -  val - the maximum component

334:    Level: intermediate

336:  Notes:
337:    Returns the value `PETSC_MIN_REAL` and negative p if the vector is of length 0.

339:    Returns the smallest index with the maximum value

341: .seealso: [](chapter_vectors), `Vec`, `VecNorm()`, `VecMin()`
342: @*/
343: PetscErrorCode VecMax(Vec x, PetscInt *p, PetscReal *val)
344: {
349:   VecLockReadPush(x);
350:   PetscLogEventBegin(VEC_Max, x, 0, 0, 0);
351:   PetscUseTypeMethod(x, max, p, val);
352:   PetscLogEventEnd(VEC_Max, x, 0, 0, 0);
353:   VecLockReadPop(x);
354:   return 0;
355: }

357: /*@C
358:    VecMin - Determines the vector component with minimum real part and its location.

360:    Collective

362:    Input Parameter:
363: .  x - the vector

365:    Output Parameters:
366: +  p - the location of val (pass `NULL` if you don't want this location)
367: -  val - the minimum component

369:    Level: intermediate

371:    Notes:
372:    Returns the value `PETSC_MAX_REAL` and negative p if the vector is of length 0.

374:    This returns the smallest index with the minimum value

376: .seealso: [](chapter_vectors), `Vec`, `VecMax()`
377: @*/
378: PetscErrorCode VecMin(Vec x, PetscInt *p, PetscReal *val)
379: {
384:   VecLockReadPush(x);
385:   PetscLogEventBegin(VEC_Min, x, 0, 0, 0);
386:   PetscUseTypeMethod(x, min, p, val);
387:   PetscLogEventEnd(VEC_Min, x, 0, 0, 0);
388:   VecLockReadPop(x);
389:   return 0;
390: }

392: /*@
393:    VecTDot - Computes an indefinite vector dot product. That is, this
394:    routine does NOT use the complex conjugate.

396:    Collective

398:    Input Parameters:
399: .  x, y - the vectors

401:    Output Parameter:
402: .  val - the dot product

404:    Level: intermediate

406:    Notes for Users of Complex Numbers:
407:    For complex vectors, VecTDot() computes the indefinite form
408: $     val = (x,y) = y^T x,
409:    where y^T denotes the transpose of y.

411:    Use VecDot() for the inner product
412: $     val = (x,y) = y^H x,
413:    where y^H denotes the conjugate transpose of y.

415: .seealso: [](chapter_vectors), `Vec`, `VecDot()`, `VecMTDot()`
416: @*/
417: PetscErrorCode VecTDot(Vec x, Vec y, PetscScalar *val)
418: {
425:   VecCheckSameSize(x, 1, y, 2);

427:   VecLockReadPush(x);
428:   VecLockReadPush(y);
429:   PetscLogEventBegin(VEC_TDot, x, y, 0, 0);
430:   PetscUseTypeMethod(x, tdot, y, val);
431:   PetscLogEventEnd(VEC_TDot, x, y, 0, 0);
432:   VecLockReadPop(x);
433:   VecLockReadPop(y);
434:   return 0;
435: }

437: /*@
438:    VecScale - Scales a vector.

440:    Not collective

442:    Input Parameters:
443: +  x - the vector
444: -  alpha - the scalar

446:    Level: intermediate

448:  Note:
449:    For a vector with n components, `VecScale()` computes
450: $      x[i] = alpha * x[i], for i=1,...,n.

452: .seealso: [](chapter_vectors), `Vec`, `VecSet()`
453: @*/
454: PetscErrorCode VecScale(Vec x, PetscScalar alpha)
455: {
456:   PetscReal norms[4];
457:   PetscBool flgs[4];

462:   VecSetErrorIfLocked(x, 1);
463:   if (alpha == (PetscScalar)1.0) return 0;

465:   /* get current stashed norms */
466:   for (PetscInt i = 0; i < 4; i++) PetscObjectComposedDataGetReal((PetscObject)x, NormIds[i], norms[i], flgs[i]);

468:   PetscLogEventBegin(VEC_Scale, x, 0, 0, 0);
469:   PetscUseTypeMethod(x, scale, alpha);
470:   PetscLogEventEnd(VEC_Scale, x, 0, 0, 0);

472:   PetscObjectStateIncrease((PetscObject)x);
473:   /* put the scaled stashed norms back into the Vec */
474:   for (PetscInt i = 0; i < 4; i++) {
475:     if (flgs[i]) PetscObjectComposedDataSetReal((PetscObject)x, NormIds[i], PetscAbsScalar(alpha) * norms[i]);
476:   }
477:   return 0;
478: }

480: /*@
481:    VecSet - Sets all components of a vector to a single scalar value.

483:    Logically Collective

485:    Input Parameters:
486: +  x  - the vector
487: -  alpha - the scalar

489:    Output Parameter:
490: .  x  - the vector

492:    Level: beginner

494:    Notes:
495:    For a vector of dimension n, `VecSet()` computes
496: $     x[i] = alpha, for i=1,...,n,
497:    so that all vector entries then equal the identical
498:    scalar value, alpha.  Use the more general routine
499:    `VecSetValues()` to set different vector entries.

501:    You CANNOT call this after you have called `VecSetValues()` but before you call
502:    `VecAssemblyBegin()`

504: .seealso: [](chapter_vectors), `Vec`, `VecSetValues()`, `VecSetValuesBlocked()`, `VecSetRandom()`
505: @*/
506: PetscErrorCode VecSet(Vec x, PetscScalar alpha)
507: {
512:   VecSetErrorIfLocked(x, 1);

514:   PetscLogEventBegin(VEC_Set, x, 0, 0, 0);
515:   PetscUseTypeMethod(x, set, alpha);
516:   PetscLogEventEnd(VEC_Set, x, 0, 0, 0);
517:   PetscObjectStateIncrease((PetscObject)x);

519:   /*  norms can be simply set (if |alpha|*N not too large) */

521:   {
522:     PetscReal      val = PetscAbsScalar(alpha);
523:     const PetscInt N   = x->map->N;

525:     if (N == 0) {
526:       PetscObjectComposedDataSetReal((PetscObject)x, NormIds[NORM_1], 0.0l);
527:       PetscObjectComposedDataSetReal((PetscObject)x, NormIds[NORM_INFINITY], 0.0);
528:       PetscObjectComposedDataSetReal((PetscObject)x, NormIds[NORM_2], 0.0);
529:       PetscObjectComposedDataSetReal((PetscObject)x, NormIds[NORM_FROBENIUS], 0.0);
530:     } else if (val > PETSC_MAX_REAL / N) {
531:       PetscObjectComposedDataSetReal((PetscObject)x, NormIds[NORM_INFINITY], val);
532:     } else {
533:       PetscObjectComposedDataSetReal((PetscObject)x, NormIds[NORM_1], N * val);
534:       PetscObjectComposedDataSetReal((PetscObject)x, NormIds[NORM_INFINITY], val);
535:       val = PetscSqrtReal((PetscReal)N) * val;
536:       PetscObjectComposedDataSetReal((PetscObject)x, NormIds[NORM_2], val);
537:       PetscObjectComposedDataSetReal((PetscObject)x, NormIds[NORM_FROBENIUS], val);
538:     }
539:   }
540:   return 0;
541: }

543: /*@
544:    VecAXPY - Computes `y = alpha x + y`.

546:    Logically Collective

548:    Input Parameters:
549: +  alpha - the scalar
550: -  x, y  - the vectors

552:    Output Parameter:
553: .  y - output vector

555:    Level: intermediate

557:    Notes:
558:     `x` and `y` MUST be different vectors
559:     This routine is optimized for alpha of 0.0, otherwise it calls the BLAS routine
560: .vb
561:     VecAXPY(y,alpha,x)                   y = alpha x           +      y
562:     VecAYPX(y,beta,x)                    y =       x           + beta y
563:     VecAXPBY(y,alpha,beta,x)             y = alpha x           + beta y
564:     VecWAXPY(w,alpha,x,y)                w = alpha x           +      y
565:     VecAXPBYPCZ(w,alpha,beta,gamma,x,y)  z = alpha x           + beta y + gamma z
566:     VecMAXPY(y,nv,alpha[],x[])           y = sum alpha[i] x[i] +      y
567: .ve

569: .seealso: [](chapter_vectors), `Vec`, `VecAYPX()`, `VecMAXPY()`, `VecWAXPY()`, `VecAXPBYPCZ()`, `VecAXPBY()`
570: @*/
571: PetscErrorCode VecAXPY(Vec y, PetscScalar alpha, Vec x)
572: {
578:   VecCheckSameSize(x, 3, y, 1);
581:   if (alpha == (PetscScalar)0.0) return 0;

583:   VecSetErrorIfLocked(y, 1);
584:   VecLockReadPush(x);
585:   PetscLogEventBegin(VEC_AXPY, x, y, 0, 0);
586:   PetscUseTypeMethod(y, axpy, alpha, x);
587:   PetscLogEventEnd(VEC_AXPY, x, y, 0, 0);
588:   VecLockReadPop(x);
589:   PetscObjectStateIncrease((PetscObject)y);
590:   return 0;
591: }

593: /*@
594:    VecAYPX - Computes `y = x + beta y`.

596:    Logically Collective

598:    Input Parameters:
599: +  beta - the scalar
600: .  x - the unscaled vector
601: -  y - the vector to be scaled

603:    Output Parameter:
604: .  y - output vector

606:    Level: intermediate

608:    Note:
609:    `x` and `y` MUST be different vectors

611:    Developer Note:
612:     The implementation is optimized for beta of -1.0, 0.0, and 1.0

614: .seealso: [](chapter_vectors), `Vec`, `VecMAXPY()`, `VecWAXPY()`, `VecAXPY()`, `VecAXPBYPCZ()`, `VecAXPBY()`
615: @*/
616: PetscErrorCode VecAYPX(Vec y, PetscScalar beta, Vec x)
617: {
623:   VecCheckSameSize(x, 1, y, 3);
626:   VecSetErrorIfLocked(y, 1);
627:   VecLockReadPush(x);
628:   if (beta == (PetscScalar)0.0) {
629:     VecCopy(x, y);
630:   } else {
631:     PetscLogEventBegin(VEC_AYPX, x, y, 0, 0);
632:     PetscUseTypeMethod(y, aypx, beta, x);
633:     PetscLogEventEnd(VEC_AYPX, x, y, 0, 0);
634:     PetscObjectStateIncrease((PetscObject)y);
635:   }
636:   VecLockReadPop(x);
637:   return 0;
638: }

640: /*@
641:    VecAXPBY - Computes `y = alpha x + beta y`.

643:    Logically Collective

645:    Input Parameters:
646: +  alpha,beta - the scalars
647: .  x - the first scaled vector
648: -  y - the second scaled vector

650:    Output Parameter:
651: .  y - output vector

653:    Level: intermediate

655:    Note:
656:    `x` and `y` MUST be different vectors

658:    Developer Note:
659:    The implementation is optimized for alpha and/or beta values of 0.0 and 1.0

661: .seealso: [](chapter_vectors), `Vec`, `VecAYPX()`, `VecMAXPY()`, `VecWAXPY()`, `VecAXPY()`, `VecAXPBYPCZ()`
662: @*/
663: PetscErrorCode VecAXPBY(Vec y, PetscScalar alpha, PetscScalar beta, Vec x)
664: {
670:   VecCheckSameSize(y, 1, x, 4);
674:   if (alpha == (PetscScalar)0.0 && beta == (PetscScalar)1.0) return 0;

676:   VecSetErrorIfLocked(y, 1);
677:   VecLockReadPush(x);
678:   PetscLogEventBegin(VEC_AXPY, y, x, 0, 0);
679:   PetscUseTypeMethod(y, axpby, alpha, beta, x);
680:   PetscLogEventEnd(VEC_AXPY, y, x, 0, 0);
681:   PetscObjectStateIncrease((PetscObject)y);
682:   VecLockReadPop(x);
683:   return 0;
684: }

686: /*@
687:    VecAXPBYPCZ - Computes `z = alpha x + beta y + gamma z`

689:    Logically Collective

691:    Input Parameters:
692: +  alpha,beta, gamma - the scalars
693: -  x, y, z  - the vectors

695:    Output Parameter:
696: .  z - output vector

698:    Level: intermediate

700:    Note:
701:    `x`, `y` and `z` must be different vectors

703:    Developer Note:
704:     The implementation is optimized for alpha of 1.0 and gamma of 1.0 or 0.0

706: .seealso: [](chapter_vectors), `Vec`, `VecAYPX()`, `VecMAXPY()`, `VecWAXPY()`, `VecAXPY()`, `VecAXPBY()`
707: @*/
708: PetscErrorCode VecAXPBYPCZ(Vec z, PetscScalar alpha, PetscScalar beta, PetscScalar gamma, Vec x, Vec y)
709: {
718:   VecCheckSameSize(x, 1, y, 5);
719:   VecCheckSameSize(x, 1, z, 6);
725:   if (alpha == (PetscScalar)0.0 && beta == (PetscScalar)0.0 && gamma == (PetscScalar)1.0) return 0;

727:   VecSetErrorIfLocked(z, 1);
728:   VecLockReadPush(x);
729:   VecLockReadPush(y);
730:   PetscLogEventBegin(VEC_AXPBYPCZ, x, y, z, 0);
731:   PetscUseTypeMethod(z, axpbypcz, alpha, beta, gamma, x, y);
732:   PetscLogEventEnd(VEC_AXPBYPCZ, x, y, z, 0);
733:   PetscObjectStateIncrease((PetscObject)z);
734:   VecLockReadPop(x);
735:   VecLockReadPop(y);
736:   return 0;
737: }

739: /*@
740:    VecWAXPY - Computes `w = alpha x + y`.

742:    Logically Collective

744:    Input Parameters:
745: +  alpha - the scalar
746: -  x, y  - the vectors

748:    Output Parameter:
749: .  w - the result

751:    Level: intermediate

753:    Note:
754:     `w` cannot be either `x` or `y`, but `x` and `y` can be the same

756:    Developer Note:
757:     The implementation is optimized for alpha of -1.0, 0.0, and 1.0

759: .seealso: [](chapter_vectors), `Vec`, `VecAXPY()`, `VecAYPX()`, `VecAXPBY()`, `VecMAXPY()`, `VecAXPBYPCZ()`
760: @*/
761: PetscErrorCode VecWAXPY(Vec w, PetscScalar alpha, Vec x, Vec y)
762: {
771:   VecCheckSameSize(x, 3, y, 4);
772:   VecCheckSameSize(x, 3, w, 1);
776:   VecSetErrorIfLocked(w, 1);

778:   VecLockReadPush(x);
779:   VecLockReadPush(y);
780:   if (alpha == (PetscScalar)0.0) {
781:     VecCopy(y, w);
782:   } else {
783:     PetscLogEventBegin(VEC_WAXPY, x, y, w, 0);
784:     PetscUseTypeMethod(w, waxpy, alpha, x, y);
785:     PetscLogEventEnd(VEC_WAXPY, x, y, w, 0);
786:     PetscObjectStateIncrease((PetscObject)w);
787:   }
788:   VecLockReadPop(x);
789:   VecLockReadPop(y);
790:   return 0;
791: }

793: /*@C
794:    VecSetValues - Inserts or adds values into certain locations of a vector.

796:    Not Collective

798:    Input Parameters:
799: +  x - vector to insert in
800: .  ni - number of elements to add
801: .  ix - indices where to add
802: .  y - array of values
803: -  iora - either `INSERT_VALUES` or `ADD_VALUES`, where
804:    `ADD_VALUES` adds values to any existing entries, and
805:    `INSERT_VALUES` replaces existing entries with new values

807:    Level: beginner

809:    Notes:
810:    `VecSetValues()` sets x[ix[i]] = y[i], for i=0,...,ni-1.

812:    Calls to `VecSetValues()` with the `INSERT_VALUES` and `ADD_VALUES`
813:    options cannot be mixed without intervening calls to the assembly
814:    routines.

816:    These values may be cached, so `VecAssemblyBegin()` and `VecAssemblyEnd()`
817:    MUST be called after all calls to `VecSetValues()` have been completed.

819:    VecSetValues() uses 0-based indices in Fortran as well as in C.

821:    If you call `VecSetOption`(x, `VEC_IGNORE_NEGATIVE_INDICES`,`PETSC_TRUE`),
822:    negative indices may be passed in ix. These rows are
823:    simply ignored. This allows easily inserting element load matrices
824:    with homogeneous Dirchlet boundary conditions that you don't want represented
825:    in the vector.

827: .seealso: [](chapter_vectors), `Vec`, `VecAssemblyBegin()`, `VecAssemblyEnd()`, `VecSetValuesLocal()`,
828:           `VecSetValue()`, `VecSetValuesBlocked()`, `InsertMode`, `INSERT_VALUES`, `ADD_VALUES`, `VecGetValues()`
829: @*/
830: PetscErrorCode VecSetValues(Vec x, PetscInt ni, const PetscInt ix[], const PetscScalar y[], InsertMode iora)
831: {
834:   if (!ni) return 0;

839:   PetscLogEventBegin(VEC_SetValues, x, 0, 0, 0);
840:   PetscUseTypeMethod(x, setvalues, ni, ix, y, iora);
841:   PetscLogEventEnd(VEC_SetValues, x, 0, 0, 0);
842:   PetscObjectStateIncrease((PetscObject)x);
843:   return 0;
844: }

846: /*@C
847:    VecGetValues - Gets values from certain locations of a vector. Currently
848:           can only get values on the same processor

850:     Not Collective

852:    Input Parameters:
853: +  x - vector to get values from
854: .  ni - number of elements to get
855: -  ix - indices where to get them from (in global 1d numbering)

857:    Output Parameter:
858: .   y - array of values

860:    Level: beginner

862:    Notes:
863:    The user provides the allocated array y; it is NOT allocated in this routine

865:    `VecGetValues()` gets y[i] = x[ix[i]], for i=0,...,ni-1.

867:    `VecAssemblyBegin()` and `VecAssemblyEnd()`  MUST be called before calling this if `VecSetValues()` or related routine has been called

869:    VecGetValues() uses 0-based indices in Fortran as well as in C.

871:    If you call `VecSetOption`(x, `VEC_IGNORE_NEGATIVE_INDICES`,`PETSC_TRUE`),
872:    negative indices may be passed in ix. These rows are
873:    simply ignored.

875: .seealso: [](chapter_vectors), `Vec`, `VecAssemblyBegin()`, `VecAssemblyEnd()`, `VecSetValues()`
876: @*/
877: PetscErrorCode VecGetValues(Vec x, PetscInt ni, const PetscInt ix[], PetscScalar y[])
878: {
880:   if (!ni) return 0;
884:   PetscUseTypeMethod(x, getvalues, ni, ix, y);
885:   return 0;
886: }

888: /*@C
889:    VecSetValuesBlocked - Inserts or adds blocks of values into certain locations of a vector.

891:    Not Collective

893:    Input Parameters:
894: +  x - vector to insert in
895: .  ni - number of blocks to add
896: .  ix - indices where to add in block count, rather than element count
897: .  y - array of values
898: -  iora - either `INSERT_VALUES` or `ADD_VALUES`, where
899:    `ADD_VALUES` adds values to any existing entries, and
900:    `INSERT_VALUES` replaces existing entries with new values

902:    Level: intermediate

904:    Notes:
905:    `VecSetValuesBlocked()` sets x[bs*ix[i]+j] = y[bs*i+j],
906:    for j=0,...,bs-1, for i=0,...,ni-1. where bs was set with VecSetBlockSize().

908:    Calls to `VecSetValuesBlocked()` with the `INSERT_VALUES` and `ADD_VALUES`
909:    options cannot be mixed without intervening calls to the assembly
910:    routines.

912:    These values may be cached, so `VecAssemblyBegin()` and `VecAssemblyEnd()`
913:    MUST be called after all calls to `VecSetValuesBlocked()` have been completed.

915:    `VecSetValuesBlocked()` uses 0-based indices in Fortran as well as in C.

917:    Negative indices may be passed in ix, these rows are
918:    simply ignored. This allows easily inserting element load matrices
919:    with homogeneous Dirchlet boundary conditions that you don't want represented
920:    in the vector.

922: .seealso: [](chapter_vectors), `Vec`, `VecAssemblyBegin()`, `VecAssemblyEnd()`, `VecSetValuesBlockedLocal()`,
923:           `VecSetValues()`
924: @*/
925: PetscErrorCode VecSetValuesBlocked(Vec x, PetscInt ni, const PetscInt ix[], const PetscScalar y[], InsertMode iora)
926: {
929:   if (!ni) return 0;

934:   PetscLogEventBegin(VEC_SetValues, x, 0, 0, 0);
935:   PetscUseTypeMethod(x, setvaluesblocked, ni, ix, y, iora);
936:   PetscLogEventEnd(VEC_SetValues, x, 0, 0, 0);
937:   PetscObjectStateIncrease((PetscObject)x);
938:   return 0;
939: }

941: /*@C
942:    VecSetValuesLocal - Inserts or adds values into certain locations of a vector,
943:    using a local ordering of the nodes.

945:    Not Collective

947:    Input Parameters:
948: +  x - vector to insert in
949: .  ni - number of elements to add
950: .  ix - indices where to add
951: .  y - array of values
952: -  iora - either `INSERT_VALUES` or `ADD_VALUES`, where
953:    `ADD_VALUES` adds values to any existing entries, and
954:    `INSERT_VALUES` replaces existing entries with new values

956:    Level: intermediate

958:    Notes:
959:    `VecSetValuesLocal()` sets x[ix[i]] = y[i], for i=0,...,ni-1.

961:    Calls to `VecSetValues()` with the `INSERT_VALUES` and `ADD_VALUES`
962:    options cannot be mixed without intervening calls to the assembly
963:    routines.

965:    These values may be cached, so `VecAssemblyBegin()` and `VecAssemblyEnd()`
966:    MUST be called after all calls to `VecSetValuesLocal()` have been completed.

968:    `VecSetValuesLocal()` uses 0-based indices in Fortran as well as in C.

970: .seealso: [](chapter_vectors), `Vec`, `VecAssemblyBegin()`, `VecAssemblyEnd()`, `VecSetValues()`, `VecSetLocalToGlobalMapping()`,
971:           `VecSetValuesBlockedLocal()`
972: @*/
973: PetscErrorCode VecSetValuesLocal(Vec x, PetscInt ni, const PetscInt ix[], const PetscScalar y[], InsertMode iora)
974: {
975:   PetscInt lixp[128], *lix = lixp;

979:   if (!ni) return 0;

984:   PetscLogEventBegin(VEC_SetValues, x, 0, 0, 0);
985:   if (!x->ops->setvalueslocal) {
986:     if (x->map->mapping) {
987:       if (ni > 128) PetscMalloc1(ni, &lix);
988:       ISLocalToGlobalMappingApply(x->map->mapping, ni, (PetscInt *)ix, lix);
989:       PetscUseTypeMethod(x, setvalues, ni, lix, y, iora);
990:       if (ni > 128) PetscFree(lix);
991:     } else PetscUseTypeMethod(x, setvalues, ni, ix, y, iora);
992:   } else PetscUseTypeMethod(x, setvalueslocal, ni, ix, y, iora);
993:   PetscLogEventEnd(VEC_SetValues, x, 0, 0, 0);
994:   PetscObjectStateIncrease((PetscObject)x);
995:   return 0;
996: }

998: /*@
999:    VecSetValuesBlockedLocal - Inserts or adds values into certain locations of a vector,
1000:    using a local ordering of the nodes.

1002:    Not Collective

1004:    Input Parameters:
1005: +  x - vector to insert in
1006: .  ni - number of blocks to add
1007: .  ix - indices where to add in block count, not element count
1008: .  y - array of values
1009: -  iora - either `INSERT_VALUES` or `ADD_VALUES`, where
1010:    `ADD_VALUES` adds values to any existing entries, and
1011:    `INSERT_VALUES` replaces existing entries with new values

1013:    Level: intermediate

1015:    Notes:
1016:    `VecSetValuesBlockedLocal()` sets x[bs*ix[i]+j] = y[bs*i+j],
1017:    for j=0,..bs-1, for i=0,...,ni-1, where bs has been set with `VecSetBlockSize()`.

1019:    Calls to `VecSetValuesBlockedLocal()` with the `INSERT_VALUES` and `ADD_VALUES`
1020:    options cannot be mixed without intervening calls to the assembly
1021:    routines.

1023:    These values may be cached, so `VecAssemblyBegin()` and `VecAssemblyEnd()`
1024:    MUST be called after all calls to `VecSetValuesBlockedLocal()` have been completed.

1026:    `VecSetValuesBlockedLocal()` uses 0-based indices in Fortran as well as in C.

1028: .seealso: [](chapter_vectors), `Vec`, `VecAssemblyBegin()`, `VecAssemblyEnd()`, `VecSetValues()`, `VecSetValuesBlocked()`,
1029:           `VecSetLocalToGlobalMapping()`
1030: @*/
1031: PetscErrorCode VecSetValuesBlockedLocal(Vec x, PetscInt ni, const PetscInt ix[], const PetscScalar y[], InsertMode iora)
1032: {
1033:   PetscInt lixp[128], *lix = lixp;

1037:   if (!ni) return 0;
1041:   PetscLogEventBegin(VEC_SetValues, x, 0, 0, 0);
1042:   if (x->map->mapping) {
1043:     if (ni > 128) PetscMalloc1(ni, &lix);
1044:     ISLocalToGlobalMappingApplyBlock(x->map->mapping, ni, (PetscInt *)ix, lix);
1045:     PetscUseTypeMethod(x, setvaluesblocked, ni, lix, y, iora);
1046:     if (ni > 128) PetscFree(lix);
1047:   } else {
1048:     PetscUseTypeMethod(x, setvaluesblocked, ni, ix, y, iora);
1049:   }
1050:   PetscLogEventEnd(VEC_SetValues, x, 0, 0, 0);
1051:   PetscObjectStateIncrease((PetscObject)x);
1052:   return 0;
1053: }

1055: /*@
1056:    VecMTDot - Computes indefinite vector multiple dot products.
1057:    That is, it does NOT use the complex conjugate.

1059:    Collective

1061:    Input Parameters:
1062: +  x - one vector
1063: .  nv - number of vectors
1064: -  y - array of vectors.  Note that vectors are pointers

1066:    Output Parameter:
1067: .  val - array of the dot products

1069:    Level: intermediate

1071:    Notes for Users of Complex Numbers:
1072:    For complex vectors, `VecMTDot()` computes the indefinite form
1073: $      val = (x,y) = y^T x,
1074:    where y^T denotes the transpose of y.

1076:    Use `VecMDot()` for the inner product
1077: $      val = (x,y) = y^H x,
1078:    where y^H denotes the conjugate transpose of y.

1080: .seealso: [](chapter_vectors), `Vec`, `VecMDot()`, `VecTDot()`
1081: @*/
1082: PetscErrorCode VecMTDot(Vec x, PetscInt nv, const Vec y[], PetscScalar val[])
1083: {
1087:   if (!nv) return 0;
1089:   for (PetscInt i = 0; i < nv; ++i) {
1093:     VecCheckSameSize(x, 1, y[i], 3);
1094:     VecLockReadPush(y[i]);
1095:   }

1098:   VecLockReadPush(x);
1099:   PetscLogEventBegin(VEC_MTDot, x, *y, 0, 0);
1100:   PetscUseTypeMethod(x, mtdot, nv, y, val);
1101:   PetscLogEventEnd(VEC_MTDot, x, *y, 0, 0);
1102:   VecLockReadPop(x);
1103:   for (PetscInt i = 0; i < nv; ++i) VecLockReadPop(y[i]);
1104:   return 0;
1105: }

1107: /*@
1108:    VecMDot - Computes vector multiple dot products.

1110:    Collective

1112:    Input Parameters:
1113: +  x - one vector
1114: .  nv - number of vectors
1115: -  y - array of vectors.

1117:    Output Parameter:
1118: .  val - array of the dot products (does not allocate the array)

1120:    Level: intermediate

1122:    Notes for Users of Complex Numbers:
1123:    For complex vectors, `VecMDot()` computes
1124: $     val = (x,y) = y^H x,
1125:    where y^H denotes the conjugate transpose of y.

1127:    Use `VecMTDot()` for the indefinite form
1128: $     val = (x,y) = y^T x,
1129:    where y^T denotes the transpose of y.

1131: .seealso: [](chapter_vectors), `Vec`, `VecMTDot()`, `VecDot()`
1132: @*/
1133: PetscErrorCode VecMDot(Vec x, PetscInt nv, const Vec y[], PetscScalar val[])
1134: {
1138:   if (!nv) return 0;
1140:   for (PetscInt i = 0; i < nv; ++i) {
1144:     VecCheckSameSize(x, 1, y[i], 3);
1145:     VecLockReadPush(y[i]);
1146:   }

1149:   VecLockReadPush(x);
1150:   PetscLogEventBegin(VEC_MDot, x, *y, 0, 0);
1151:   PetscUseTypeMethod(x, mdot, nv, y, val);
1152:   PetscLogEventEnd(VEC_MDot, x, *y, 0, 0);
1153:   VecLockReadPop(x);
1154:   for (PetscInt i = 0; i < nv; ++i) VecLockReadPop(y[i]);
1155:   return 0;
1156: }

1158: /*@
1159:    VecMAXPY - Computes `y = y + sum alpha[i] x[i]`

1161:    Logically Collective

1163:    Input Parameters:
1164: +  nv - number of scalars and x-vectors
1165: .  alpha - array of scalars
1166: .  y - one vector
1167: -  x - array of vectors

1169:    Level: intermediate

1171:    Note:
1172:     `y` cannot be any of the `x` vectors

1174: .seealso: [](chapter_vectors), `Vec`, `VecAYPX()`, `VecWAXPY()`, `VecAXPY()`, `VecAXPBYPCZ()`, `VecAXPBY()`
1175: @*/
1176: PetscErrorCode VecMAXPY(Vec y, PetscInt nv, const PetscScalar alpha[], Vec x[])
1177: {
1180:   VecSetErrorIfLocked(y, 1);
1182:   if (nv) {
1183:     PetscInt zeros = 0;

1187:     for (PetscInt i = 0; i < nv; ++i) {
1192:       VecCheckSameSize(y, 1, x[i], 4);
1194:       VecLockReadPush(x[i]);
1195:       zeros += alpha[i] == (PetscScalar)0.0;
1196:     }

1198:     if (zeros < nv) {
1199:       PetscLogEventBegin(VEC_MAXPY, y, *x, 0, 0);
1200:       PetscUseTypeMethod(y, maxpy, nv, alpha, x);
1201:       PetscLogEventEnd(VEC_MAXPY, y, *x, 0, 0);
1202:       PetscObjectStateIncrease((PetscObject)y);
1203:     }

1205:     for (PetscInt i = 0; i < nv; ++i) VecLockReadPop(x[i]);
1206:   }
1207:   return 0;
1208: }

1210: /*@
1211:    VecConcatenate - Creates a new vector that is a vertical concatenation of all the given array of vectors
1212:                     in the order they appear in the array. The concatenated vector resides on the same
1213:                     communicator and is the same type as the source vectors.

1215:    Collective

1217:    Input Parameters:
1218: +  nx   - number of vectors to be concatenated
1219: -  X    - array containing the vectors to be concatenated in the order of concatenation

1221:    Output Parameters:
1222: +  Y    - concatenated vector
1223: -  x_is - array of index sets corresponding to the concatenated components of Y (`NULL` if not needed)

1225:    Level: advanced

1227:    Notes:
1228:    Concatenation is similar to the functionality of a `VECNEST` object; they both represent combination of
1229:    different vector spaces. However, concatenated vectors do not store any information about their
1230:    sub-vectors and own their own data. Consequently, this function provides index sets to enable the
1231:    manipulation of data in the concatenated vector that corresponds to the original components at creation.

1233:    This is a useful tool for outer loop algorithms, particularly constrained optimizers, where the solver
1234:    has to operate on combined vector spaces and cannot utilize `VECNEST` objects due to incompatibility with
1235:    bound projections.

1237: .seealso: [](chapter_vectors), `Vec`, `VECNEST`, `VECSCATTER`, `VecScatterCreate()`
1238: @*/
1239: PetscErrorCode VecConcatenate(PetscInt nx, const Vec X[], Vec *Y, IS *x_is[])
1240: {
1241:   MPI_Comm comm;
1242:   VecType  vec_type;
1243:   Vec      Ytmp, Xtmp;
1244:   IS      *is_tmp;
1245:   PetscInt i, shift = 0, Xnl, Xng, Xbegin;


1252:   if ((*X)->ops->concatenate) {
1253:     /* use the dedicated concatenation function if available */
1254:     (*(*X)->ops->concatenate)(nx, X, Y, x_is);
1255:   } else {
1256:     /* loop over vectors and start creating IS */
1257:     comm = PetscObjectComm((PetscObject)(*X));
1258:     VecGetType(*X, &vec_type);
1259:     PetscMalloc1(nx, &is_tmp);
1260:     for (i = 0; i < nx; i++) {
1261:       VecGetSize(X[i], &Xng);
1262:       VecGetLocalSize(X[i], &Xnl);
1263:       VecGetOwnershipRange(X[i], &Xbegin, NULL);
1264:       ISCreateStride(comm, Xnl, shift + Xbegin, 1, &is_tmp[i]);
1265:       shift += Xng;
1266:     }
1267:     /* create the concatenated vector */
1268:     VecCreate(comm, &Ytmp);
1269:     VecSetType(Ytmp, vec_type);
1270:     VecSetSizes(Ytmp, PETSC_DECIDE, shift);
1271:     VecSetUp(Ytmp);
1272:     /* copy data from X array to Y and return */
1273:     for (i = 0; i < nx; i++) {
1274:       VecGetSubVector(Ytmp, is_tmp[i], &Xtmp);
1275:       VecCopy(X[i], Xtmp);
1276:       VecRestoreSubVector(Ytmp, is_tmp[i], &Xtmp);
1277:     }
1278:     *Y = Ytmp;
1279:     if (x_is) {
1280:       *x_is = is_tmp;
1281:     } else {
1282:       for (i = 0; i < nx; i++) ISDestroy(&is_tmp[i]);
1283:       PetscFree(is_tmp);
1284:     }
1285:   }
1286:   return 0;
1287: }

1289: /* A helper function for VecGetSubVector to check if we can implement it with no-copy (i.e. the subvector shares
1290:    memory with the original vector), and the block size of the subvector.

1292:     Input Parameters:
1293: +   X - the original vector
1294: -   is - the index set of the subvector

1296:     Output Parameters:
1297: +   contig - PETSC_TRUE if the index set refers to contiguous entries on this process, else PETSC_FALSE
1298: .   start  - start of contiguous block, as an offset from the start of the ownership range of the original vector
1299: -   blocksize - the block size of the subvector

1301: */
1302: PetscErrorCode VecGetSubVectorContiguityAndBS_Private(Vec X, IS is, PetscBool *contig, PetscInt *start, PetscInt *blocksize)
1303: {
1304:   PetscInt  gstart, gend, lstart;
1305:   PetscBool red[2] = {PETSC_TRUE /*contiguous*/, PETSC_TRUE /*validVBS*/};
1306:   PetscInt  n, N, ibs, vbs, bs = -1;

1308:   ISGetLocalSize(is, &n);
1309:   ISGetSize(is, &N);
1310:   ISGetBlockSize(is, &ibs);
1311:   VecGetBlockSize(X, &vbs);
1312:   VecGetOwnershipRange(X, &gstart, &gend);
1313:   ISContiguousLocal(is, gstart, gend, &lstart, &red[0]);
1314:   /* block size is given by IS if ibs > 1; otherwise, check the vector */
1315:   if (ibs > 1) {
1316:     MPIU_Allreduce(MPI_IN_PLACE, red, 1, MPIU_BOOL, MPI_LAND, PetscObjectComm((PetscObject)is));
1317:     bs = ibs;
1318:   } else {
1319:     if (n % vbs || vbs == 1) red[1] = PETSC_FALSE; /* this process invalidate the collectiveness of block size */
1320:     MPIU_Allreduce(MPI_IN_PLACE, red, 2, MPIU_BOOL, MPI_LAND, PetscObjectComm((PetscObject)is));
1321:     if (red[0] && red[1]) bs = vbs; /* all processes have a valid block size and the access will be contiguous */
1322:   }

1324:   *contig    = red[0];
1325:   *start     = lstart;
1326:   *blocksize = bs;
1327:   return 0;
1328: }

1330: /* A helper function for VecGetSubVector, to be used when we have to build a standalone subvector through VecScatter

1332:     Input Parameters:
1333: +   X - the original vector
1334: .   is - the index set of the subvector
1335: -   bs - the block size of the subvector, gotten from VecGetSubVectorContiguityAndBS_Private()

1337:     Output Parameters:
1338: .   Z  - the subvector, which will compose the VecScatter context on output
1339: */
1340: PetscErrorCode VecGetSubVectorThroughVecScatter_Private(Vec X, IS is, PetscInt bs, Vec *Z)
1341: {
1342:   PetscInt   n, N;
1343:   VecScatter vscat;
1344:   Vec        Y;

1346:   ISGetLocalSize(is, &n);
1347:   ISGetSize(is, &N);
1348:   VecCreate(PetscObjectComm((PetscObject)is), &Y);
1349:   VecSetSizes(Y, n, N);
1350:   VecSetBlockSize(Y, bs);
1351:   VecSetType(Y, ((PetscObject)X)->type_name);
1352:   VecScatterCreate(X, is, Y, NULL, &vscat);
1353:   VecScatterBegin(vscat, X, Y, INSERT_VALUES, SCATTER_FORWARD);
1354:   VecScatterEnd(vscat, X, Y, INSERT_VALUES, SCATTER_FORWARD);
1355:   PetscObjectCompose((PetscObject)Y, "VecGetSubVector_Scatter", (PetscObject)vscat);
1356:   VecScatterDestroy(&vscat);
1357:   *Z = Y;
1358:   return 0;
1359: }

1361: /*@
1362:    VecGetSubVector - Gets a vector representing part of another vector

1364:    Collective

1366:    Input Parameters:
1367: +  X - vector from which to extract a subvector
1368: -  is - index set representing portion of X to extract

1370:    Output Parameter:
1371: .  Y - subvector corresponding to is

1373:    Level: advanced

1375:    Notes:
1376:    The subvector `Y` should be returned with `VecRestoreSubVector()`.
1377:    `X` and must be defined on the same communicator

1379:    This function may return a subvector without making a copy, therefore it is not safe to use the original vector while
1380:    modifying the subvector.  Other non-overlapping subvectors can still be obtained from X using this function.
1381:    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.

1383: .seealso: [](chapter_vectors), `Vec`, `IS`, `VECNEST`, `MatCreateSubMatrix()`
1384: @*/
1385: PetscErrorCode VecGetSubVector(Vec X, IS is, Vec *Y)
1386: {
1387:   Vec Z;

1393:   if (X->ops->getsubvector) {
1394:     PetscUseTypeMethod(X, getsubvector, is, &Z);
1395:   } else { /* Default implementation currently does no caching */
1396:     PetscBool contig;
1397:     PetscInt  n, N, start, bs;

1399:     ISGetLocalSize(is, &n);
1400:     ISGetSize(is, &N);
1401:     VecGetSubVectorContiguityAndBS_Private(X, is, &contig, &start, &bs);
1402:     if (contig) { /* We can do a no-copy implementation */
1403:       const PetscScalar *x;
1404:       PetscInt           state = 0;
1405:       PetscBool          isstd, iscuda, iship;

1407:       PetscObjectTypeCompareAny((PetscObject)X, &isstd, VECSEQ, VECMPI, VECSTANDARD, "");
1408:       PetscObjectTypeCompareAny((PetscObject)X, &iscuda, VECSEQCUDA, VECMPICUDA, "");
1409:       PetscObjectTypeCompareAny((PetscObject)X, &iship, VECSEQHIP, VECMPIHIP, "");
1410:       if (iscuda) {
1411: #if defined(PETSC_HAVE_CUDA)
1412:         const PetscScalar *x_d;
1413:         PetscMPIInt        size;
1414:         PetscOffloadMask   flg;

1416:         VecCUDAGetArrays_Private(X, &x, &x_d, &flg);
1419:         if (x) x += start;
1420:         if (x_d) x_d += start;
1421:         MPI_Comm_size(PetscObjectComm((PetscObject)X), &size);
1422:         if (size == 1) {
1423:           VecCreateSeqCUDAWithArrays(PetscObjectComm((PetscObject)X), bs, n, x, x_d, &Z);
1424:         } else {
1425:           VecCreateMPICUDAWithArrays(PetscObjectComm((PetscObject)X), bs, n, N, x, x_d, &Z);
1426:         }
1427:         Z->offloadmask = flg;
1428: #endif
1429:       } else if (iship) {
1430: #if defined(PETSC_HAVE_HIP)
1431:         const PetscScalar *x_d;
1432:         PetscMPIInt        size;
1433:         PetscOffloadMask   flg;

1435:         VecHIPGetArrays_Private(X, &x, &x_d, &flg);
1438:         if (x) x += start;
1439:         if (x_d) x_d += start;
1440:         MPI_Comm_size(PetscObjectComm((PetscObject)X), &size);
1441:         if (size == 1) {
1442:           VecCreateSeqHIPWithArrays(PetscObjectComm((PetscObject)X), bs, n, x, x_d, &Z);
1443:         } else {
1444:           VecCreateMPIHIPWithArrays(PetscObjectComm((PetscObject)X), bs, n, N, x, x_d, &Z);
1445:         }
1446:         Z->offloadmask = flg;
1447: #endif
1448:       } else if (isstd) {
1449:         PetscMPIInt size;

1451:         MPI_Comm_size(PetscObjectComm((PetscObject)X), &size);
1452:         VecGetArrayRead(X, &x);
1453:         if (x) x += start;
1454:         if (size == 1) {
1455:           VecCreateSeqWithArray(PetscObjectComm((PetscObject)X), bs, n, x, &Z);
1456:         } else {
1457:           VecCreateMPIWithArray(PetscObjectComm((PetscObject)X), bs, n, N, x, &Z);
1458:         }
1459:         VecRestoreArrayRead(X, &x);
1460:       } else { /* default implementation: use place array */
1461:         VecGetArrayRead(X, &x);
1462:         VecCreate(PetscObjectComm((PetscObject)X), &Z);
1463:         VecSetType(Z, ((PetscObject)X)->type_name);
1464:         VecSetSizes(Z, n, N);
1465:         VecSetBlockSize(Z, bs);
1466:         VecPlaceArray(Z, x ? x + start : NULL);
1467:         VecRestoreArrayRead(X, &x);
1468:       }

1470:       /* this is relevant only in debug mode */
1471:       VecLockGet(X, &state);
1472:       if (state) VecLockReadPush(Z);
1473:       Z->ops->placearray   = NULL;
1474:       Z->ops->replacearray = NULL;
1475:     } else { /* Have to create a scatter and do a copy */
1476:       VecGetSubVectorThroughVecScatter_Private(X, is, bs, &Z);
1477:     }
1478:   }
1479:   /* Record the state when the subvector was gotten so we know whether its values need to be put back */
1480:   if (VecGetSubVectorSavedStateId < 0) PetscObjectComposedDataRegister(&VecGetSubVectorSavedStateId);
1481:   PetscObjectComposedDataSetInt((PetscObject)Z, VecGetSubVectorSavedStateId, 1);
1482:   *Y = Z;
1483:   return 0;
1484: }

1486: /*@
1487:    VecRestoreSubVector - Restores a subvector extracted using `VecGetSubVector()`

1489:    Collective

1491:    Input Parameters:
1492: + X - vector from which subvector was obtained
1493: . is - index set representing the subset of `X`
1494: - Y - subvector being restored

1496:    Level: advanced

1498: .seealso: [](chapter_vectors), `Vec`, `IS`, `VecGetSubVector()`
1499: @*/
1500: PetscErrorCode VecRestoreSubVector(Vec X, IS is, Vec *Y)
1501: {
1502:   PETSC_UNUSED PetscObjectState dummystate = 0;
1503:   PetscBool                     unchanged;


1511:   if (X->ops->restoresubvector) PetscUseTypeMethod(X, restoresubvector, is, Y);
1512:   else {
1513:     PetscObjectComposedDataGetInt((PetscObject)*Y, VecGetSubVectorSavedStateId, dummystate, unchanged);
1514:     if (!unchanged) { /* If Y's state has not changed since VecGetSubVector(), we only need to destroy Y */
1515:       VecScatter scatter;
1516:       PetscInt   state;

1518:       VecLockGet(X, &state);

1521:       PetscObjectQuery((PetscObject)*Y, "VecGetSubVector_Scatter", (PetscObject *)&scatter);
1522:       if (scatter) {
1523:         VecScatterBegin(scatter, *Y, X, INSERT_VALUES, SCATTER_REVERSE);
1524:         VecScatterEnd(scatter, *Y, X, INSERT_VALUES, SCATTER_REVERSE);
1525:       } else {
1526:         PetscBool iscuda, iship;
1527:         PetscObjectTypeCompareAny((PetscObject)X, &iscuda, VECSEQCUDA, VECMPICUDA, "");
1528:         PetscObjectTypeCompareAny((PetscObject)X, &iship, VECSEQHIP, VECMPIHIP, "");

1530:         if (iscuda) {
1531: #if defined(PETSC_HAVE_CUDA)
1532:           PetscOffloadMask ymask = (*Y)->offloadmask;

1534:           /* The offloadmask of X dictates where to move memory
1535:               If X GPU data is valid, then move Y data on GPU if needed
1536:               Otherwise, move back to the CPU */
1537:           switch (X->offloadmask) {
1538:           case PETSC_OFFLOAD_BOTH:
1539:             if (ymask == PETSC_OFFLOAD_CPU) {
1540:               VecCUDAResetArray(*Y);
1541:             } else if (ymask == PETSC_OFFLOAD_GPU) {
1542:               X->offloadmask = PETSC_OFFLOAD_GPU;
1543:             }
1544:             break;
1545:           case PETSC_OFFLOAD_GPU:
1546:             if (ymask == PETSC_OFFLOAD_CPU) VecCUDAResetArray(*Y);
1547:             break;
1548:           case PETSC_OFFLOAD_CPU:
1549:             if (ymask == PETSC_OFFLOAD_GPU) VecResetArray(*Y);
1550:             break;
1551:           case PETSC_OFFLOAD_UNALLOCATED:
1552:           case PETSC_OFFLOAD_KOKKOS:
1553:             SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "This should not happen");
1554:           }
1555: #endif
1556:         } else if (iship) {
1557: #if defined(PETSC_HAVE_HIP)
1558:           PetscOffloadMask ymask = (*Y)->offloadmask;

1560:           /* The offloadmask of X dictates where to move memory
1561:               If X GPU data is valid, then move Y data on GPU if needed
1562:               Otherwise, move back to the CPU */
1563:           switch (X->offloadmask) {
1564:           case PETSC_OFFLOAD_BOTH:
1565:             if (ymask == PETSC_OFFLOAD_CPU) {
1566:               VecHIPResetArray(*Y);
1567:             } else if (ymask == PETSC_OFFLOAD_GPU) {
1568:               X->offloadmask = PETSC_OFFLOAD_GPU;
1569:             }
1570:             break;
1571:           case PETSC_OFFLOAD_GPU:
1572:             if (ymask == PETSC_OFFLOAD_CPU) VecHIPResetArray(*Y);
1573:             break;
1574:           case PETSC_OFFLOAD_CPU:
1575:             if (ymask == PETSC_OFFLOAD_GPU) VecResetArray(*Y);
1576:             break;
1577:           case PETSC_OFFLOAD_UNALLOCATED:
1578:           case PETSC_OFFLOAD_KOKKOS:
1579:             SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "This should not happen");
1580:           }
1581: #endif
1582:         } else {
1583:           /* If OpenCL vecs updated the device memory, this triggers a copy on the CPU */
1584:           VecResetArray(*Y);
1585:         }
1586:         PetscObjectStateIncrease((PetscObject)X);
1587:       }
1588:     }
1589:   }
1590:   VecDestroy(Y);
1591:   return 0;
1592: }

1594: /*@
1595:    VecCreateLocalVector - Creates a vector object suitable for use with `VecGetLocalVector()` and friends. You must call `VecDestroy()` when the
1596:    vector is no longer needed.

1598:    Not collective.

1600:    Input parameter:
1601: .  v - The vector for which the local vector is desired.

1603:    Output parameter:
1604: .  w - Upon exit this contains the local vector.

1606:    Level: beginner

1608: .seealso: [](chapter_vectors), `Vec`, `VecGetLocalVectorRead()`, `VecRestoreLocalVectorRead()`, `VecGetLocalVector()`, `VecRestoreLocalVector()`
1609: @*/
1610: PetscErrorCode VecCreateLocalVector(Vec v, Vec *w)
1611: {
1612:   PetscMPIInt size;

1616:   MPI_Comm_size(PetscObjectComm((PetscObject)v), &size);
1617:   if (size == 1) VecDuplicate(v, w);
1618:   else if (v->ops->createlocalvector) PetscUseTypeMethod(v, createlocalvector, w);
1619:   else {
1620:     VecType  type;
1621:     PetscInt n;

1623:     VecCreate(PETSC_COMM_SELF, w);
1624:     VecGetLocalSize(v, &n);
1625:     VecSetSizes(*w, n, n);
1626:     VecGetBlockSize(v, &n);
1627:     VecSetBlockSize(*w, n);
1628:     VecGetType(v, &type);
1629:     VecSetType(*w, type);
1630:   }
1631:   return 0;
1632: }

1634: /*@
1635:    VecGetLocalVectorRead - Maps the local portion of a vector into a
1636:    vector.  You must call `VecRestoreLocalVectorRead()` when the local
1637:    vector is no longer needed.

1639:    Not collective.

1641:    Input parameter:
1642: .  v - The vector for which the local vector is desired.

1644:    Output parameter:
1645: .  w - Upon exit this contains the local vector.

1647:    Level: beginner

1649:    Notes:
1650:    This function is similar to `VecGetArrayRead()` which maps the local
1651:    portion into a raw pointer.  `VecGetLocalVectorRead()` is usually
1652:    almost as efficient as `VecGetArrayRead()` but in certain circumstances
1653:    `VecGetLocalVectorRead()` can be much more efficient than
1654:    `VecGetArrayRead()`.  This is because the construction of a contiguous
1655:    array representing the vector data required by `VecGetArrayRead()` can
1656:    be an expensive operation for certain vector types.  For example, for
1657:    GPU vectors `VecGetArrayRead()` requires that the data between device
1658:    and host is synchronized.

1660:    Unlike `VecGetLocalVector()`, this routine is not collective and
1661:    preserves cached information.

1663: .seealso: [](chapter_vectors), `Vec`, `VecCreateLocalVector()`, `VecRestoreLocalVectorRead()`, `VecGetLocalVector()`, `VecGetArrayRead()`, `VecGetArray()`
1664: @*/
1665: PetscErrorCode VecGetLocalVectorRead(Vec v, Vec w)
1666: {
1669:   VecCheckSameLocalSize(v, 1, w, 2);
1670:   if (v->ops->getlocalvectorread) {
1671:     PetscUseTypeMethod(v, getlocalvectorread, w);
1672:   } else {
1673:     PetscScalar *a;

1675:     VecGetArrayRead(v, (const PetscScalar **)&a);
1676:     VecPlaceArray(w, a);
1677:   }
1678:   PetscObjectStateIncrease((PetscObject)w);
1679:   VecLockReadPush(v);
1680:   VecLockReadPush(w);
1681:   return 0;
1682: }

1684: /*@
1685:    VecRestoreLocalVectorRead - Unmaps the local portion of a vector
1686:    previously mapped into a vector using `VecGetLocalVectorRead()`.

1688:    Not collective.

1690:    Input parameter:
1691: +  v - The local portion of this vector was previously mapped into w using `VecGetLocalVectorRead()`.
1692: -  w - The vector into which the local portion of v was mapped.

1694:    Level: beginner

1696: .seealso: [](chapter_vectors), `Vec`, `VecCreateLocalVector()`, `VecGetLocalVectorRead()`, `VecGetLocalVector()`, `VecGetArrayRead()`, `VecGetArray()`
1697: @*/
1698: PetscErrorCode VecRestoreLocalVectorRead(Vec v, Vec w)
1699: {
1702:   if (v->ops->restorelocalvectorread) {
1703:     PetscUseTypeMethod(v, restorelocalvectorread, w);
1704:   } else {
1705:     const PetscScalar *a;

1707:     VecGetArrayRead(w, &a);
1708:     VecRestoreArrayRead(v, &a);
1709:     VecResetArray(w);
1710:   }
1711:   VecLockReadPop(v);
1712:   VecLockReadPop(w);
1713:   PetscObjectStateIncrease((PetscObject)w);
1714:   return 0;
1715: }

1717: /*@
1718:    VecGetLocalVector - Maps the local portion of a vector into a
1719:    vector.

1721:    Collective on v, not collective on w.

1723:    Input parameter:
1724: .  v - The vector for which the local vector is desired.

1726:    Output parameter:
1727: .  w - Upon exit this contains the local vector.

1729:    Level: beginner

1731:    Notes:
1732:    This function is similar to `VecGetArray()` which maps the local
1733:    portion into a raw pointer.  `VecGetLocalVector()` is usually about as
1734:    efficient as `VecGetArray()` but in certain circumstances
1735:    `VecGetLocalVector()` can be much more efficient than `VecGetArray()`.
1736:    This is because the construction of a contiguous array representing
1737:    the vector data required by `VecGetArray()` can be an expensive
1738:    operation for certain vector types.  For example, for GPU vectors
1739:    `VecGetArray()` requires that the data between device and host is
1740:    synchronized.

1742: .seealso: [](chapter_vectors), `Vec`, `VecCreateLocalVector()`, `VecRestoreLocalVector()`, `VecGetLocalVectorRead()`, `VecGetArrayRead()`, `VecGetArray()`
1743: @*/
1744: PetscErrorCode VecGetLocalVector(Vec v, Vec w)
1745: {
1748:   VecCheckSameLocalSize(v, 1, w, 2);
1749:   if (v->ops->getlocalvector) {
1750:     PetscUseTypeMethod(v, getlocalvector, w);
1751:   } else {
1752:     PetscScalar *a;

1754:     VecGetArray(v, &a);
1755:     VecPlaceArray(w, a);
1756:   }
1757:   PetscObjectStateIncrease((PetscObject)w);
1758:   return 0;
1759: }

1761: /*@
1762:    VecRestoreLocalVector - Unmaps the local portion of a vector
1763:    previously mapped into a vector using `VecGetLocalVector()`.

1765:    Logically collective.

1767:    Input parameter:
1768: +  v - The local portion of this vector was previously mapped into `w` using `VecGetLocalVector()`.
1769: -  w - The vector into which the local portion of v was mapped.

1771:    Level: beginner

1773: .seealso: [](chapter_vectors), `Vec`, `VecCreateLocalVector()`, `VecGetLocalVector()`, `VecGetLocalVectorRead()`, `VecRestoreLocalVectorRead()`, `LocalVectorRead()`, `VecGetArrayRead()`, `VecGetArray()`
1774: @*/
1775: PetscErrorCode VecRestoreLocalVector(Vec v, Vec w)
1776: {
1779:   if (v->ops->restorelocalvector) {
1780:     PetscUseTypeMethod(v, restorelocalvector, w);
1781:   } else {
1782:     PetscScalar *a;
1783:     VecGetArray(w, &a);
1784:     VecRestoreArray(v, &a);
1785:     VecResetArray(w);
1786:   }
1787:   PetscObjectStateIncrease((PetscObject)w);
1788:   PetscObjectStateIncrease((PetscObject)v);
1789:   return 0;
1790: }

1792: /*@C
1793:    VecGetArray - Returns a pointer to a contiguous array that contains this
1794:    processor's portion of the vector data. For the standard PETSc
1795:    vectors, VecGetArray() returns a pointer to the local data array and
1796:    does not use any copies. If the underlying vector data is not stored
1797:    in a contiguous array this routine will copy the data to a contiguous
1798:    array and return a pointer to that. You MUST call `VecRestoreArray()`
1799:    when you no longer need access to the array.

1801:    Logically Collective

1803:    Input Parameter:
1804: .  x - the vector

1806:    Output Parameter:
1807: .  a - location to put pointer to the array

1809:    Level: beginner

1811:    Fortran Note:
1812:    Use `VecGetArrayF90()`, this routine is deprecated

1814:    This routine is used differently from Fortran
1815: .vb
1816:     Vec         x
1817:     PetscScalar x_array(1)
1818:     PetscOffset i_x
1819:     PetscErrorCode ierr
1820:     call VecGetArray(x,x_array,i_x,ierr)

1822:     Access first local entry in vector with
1823:     value = x_array(i_x + 1)

1825:     ...... other code
1826:     call VecRestoreArray(x,x_array,i_x,ierr)
1827: .ve
1828: .seealso: [](chapter_vectors), `Vec`, `VecRestoreArray()`, `VecGetArrayRead()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecGetArrayReadF90()`, `VecPlaceArray()`, `VecGetArray2d()`,
1829:           `VecGetArrayPair()`, `VecRestoreArrayPair()`, `VecGetArrayWrite()`, `VecRestoreArrayWrite()`
1830: @*/
1831: PetscErrorCode VecGetArray(Vec x, PetscScalar **a)
1832: {
1834:   VecSetErrorIfLocked(x, 1);
1835:   if (x->ops->getarray) { /* The if-else order matters! VECNEST, VECCUDA etc should have ops->getarray while VECCUDA etc are petscnative */
1836:     PetscUseTypeMethod(x, getarray, a);
1837:   } else if (x->petscnative) { /* VECSTANDARD */
1838:     *a = *((PetscScalar **)x->data);
1839:   } else SETERRQ(PetscObjectComm((PetscObject)x), PETSC_ERR_SUP, "Cannot get array for vector type \"%s\"", ((PetscObject)x)->type_name);
1840:   return 0;
1841: }

1843: /*@C
1844:    VecRestoreArray - Restores a vector after `VecGetArray()` has been called.

1846:    Logically Collective

1848:    Input Parameters:
1849: +  x - the vector
1850: -  a - location of pointer to array obtained from `VecGetArray()`

1852:    Level: beginner

1854: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArrayRead()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecRestoreArrayReadF90()`, `VecPlaceArray()`, `VecRestoreArray2d()`,
1855:           `VecGetArrayPair()`, `VecRestoreArrayPair()`
1856: @*/
1857: PetscErrorCode VecRestoreArray(Vec x, PetscScalar **a)
1858: {
1861:   if (x->ops->restorearray) {
1862:     PetscUseTypeMethod(x, restorearray, a);
1864:   if (a) *a = NULL;
1865:   PetscObjectStateIncrease((PetscObject)x);
1866:   return 0;
1867: }
1868: /*@C
1869:    VecGetArrayRead - Get read-only pointer to contiguous array containing this processor's portion of the vector data.

1871:    Not Collective

1873:    Input Parameter:
1874: .  x - the vector

1876:    Output Parameter:
1877: .  a - the array

1879:    Level: beginner

1881:    Notes:
1882:    The array must be returned using a matching call to `VecRestoreArrayRead()`.

1884:    Unlike `VecGetArray()`, this routine is not collective and preserves cached information like vector norms.

1886:    Standard PETSc vectors use contiguous storage so that this routine does not perform a copy.  Other vector
1887:    implementations may require a copy, but must such implementations should cache the contiguous representation so that
1888:    only one copy is performed when this routine is called multiple times in sequence.

1890: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayPair()`, `VecRestoreArrayPair()`
1891: @*/
1892: PetscErrorCode VecGetArrayRead(Vec x, const PetscScalar **a)
1893: {
1896:   if (x->ops->getarrayread) {
1897:     PetscUseTypeMethod(x, getarrayread, a);
1898:   } else if (x->ops->getarray) {
1899:     /* VECNEST, VECCUDA, VECKOKKOS etc */
1900:     PetscUseTypeMethod(x, getarray, (PetscScalar **)a);
1901:   } else if (x->petscnative) {
1902:     /* VECSTANDARD */
1903:     *a = *((PetscScalar **)x->data);
1904:   } else SETERRQ(PetscObjectComm((PetscObject)x), PETSC_ERR_SUP, "Cannot get array read for vector type \"%s\"", ((PetscObject)x)->type_name);
1905:   return 0;
1906: }

1908: /*@C
1909:    VecRestoreArrayRead - Restore array obtained with `VecGetArrayRead()`

1911:    Not Collective

1913:    Input Parameters:
1914: +  vec - the vector
1915: -  array - the array

1917:    Level: beginner

1919: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayPair()`, `VecRestoreArrayPair()`
1920: @*/
1921: PetscErrorCode VecRestoreArrayRead(Vec x, const PetscScalar **a)
1922: {
1925:   if (x->petscnative) { /* VECSTANDARD, VECCUDA, VECKOKKOS etc */
1926:     /* nothing */
1927:   } else if (x->ops->restorearrayread) { /* VECNEST */
1928:     PetscUseTypeMethod(x, restorearrayread, a);
1929:   } else { /* No one? */
1930:     PetscUseTypeMethod(x, restorearray, (PetscScalar **)a);
1931:   }
1932:   if (a) *a = NULL;
1933:   return 0;
1934: }

1936: /*@C
1937:    VecGetArrayWrite - Returns a pointer to a contiguous array that WILL contains this
1938:    processor's portion of the vector data. The values in this array are NOT valid, the routine calling this
1939:    routine is responsible for putting values into the array; any values it does not set will be invalid

1941:    Logically Collective

1943:    Input Parameter:
1944: .  x - the vector

1946:    Output Parameter:
1947: .  a - location to put pointer to the array

1949:    Level: intermediate

1951:    This is for vectors associate with GPUs, the vector is not copied up before giving access. If you need correct
1952:    values in the array use `VecGetArray()`

1954: .seealso: [](chapter_vectors), `Vec`, `VecRestoreArray()`, `VecGetArrayRead()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecGetArrayReadF90()`, `VecPlaceArray()`, `VecGetArray2d()`,
1955:           `VecGetArrayPair()`, `VecRestoreArrayPair()`, `VecGetArray()`, `VecRestoreArrayWrite()`
1956: @*/
1957: PetscErrorCode VecGetArrayWrite(Vec x, PetscScalar **a)
1958: {
1961:   VecSetErrorIfLocked(x, 1);
1962:   if (x->ops->getarraywrite) {
1963:     PetscUseTypeMethod(x, getarraywrite, a);
1964:   } else {
1965:     VecGetArray(x, a);
1966:   }
1967:   return 0;
1968: }

1970: /*@C
1971:    VecRestoreArrayWrite - Restores a vector after `VecGetArrayWrite()` has been called.

1973:    Logically Collective

1975:    Input Parameters:
1976: +  x - the vector
1977: -  a - location of pointer to array obtained from `VecGetArray()`

1979:    Level: beginner

1981: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArrayRead()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecRestoreArrayReadF90()`, `VecPlaceArray()`, `VecRestoreArray2d()`,
1982:           `VecGetArrayPair()`, `VecRestoreArrayPair()`, `VecGetArrayWrite()`
1983: @*/
1984: PetscErrorCode VecRestoreArrayWrite(Vec x, PetscScalar **a)
1985: {
1988:   if (x->ops->restorearraywrite) {
1989:     PetscUseTypeMethod(x, restorearraywrite, a);
1990:   } else if (x->ops->restorearray) {
1991:     PetscUseTypeMethod(x, restorearray, a);
1992:   }
1993:   if (a) *a = NULL;
1994:   PetscObjectStateIncrease((PetscObject)x);
1995:   return 0;
1996: }

1998: /*@C
1999:    VecGetArrays - Returns a pointer to the arrays in a set of vectors
2000:    that were created by a call to `VecDuplicateVecs()`.  You MUST call
2001:    `VecRestoreArrays()` when you no longer need access to the array.

2003:    Logically Collective; No Fortran Support

2005:    Input Parameters:
2006: +  x - the vectors
2007: -  n - the number of vectors

2009:    Output Parameter:
2010: .  a - location to put pointer to the array

2012:    Level: intermediate

2014: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArrays()`
2015: @*/
2016: PetscErrorCode VecGetArrays(const Vec x[], PetscInt n, PetscScalar **a[])
2017: {
2018:   PetscInt      i;
2019:   PetscScalar **q;

2025:   PetscMalloc1(n, &q);
2026:   for (i = 0; i < n; ++i) VecGetArray(x[i], &q[i]);
2027:   *a = q;
2028:   return 0;
2029: }

2031: /*@C
2032:    VecRestoreArrays - Restores a group of vectors after `VecGetArrays()`
2033:    has been called.

2035:    Logically Collective; No Fortran Support

2037:    Input Parameters:
2038: +  x - the vector
2039: .  n - the number of vectors
2040: -  a - location of pointer to arrays obtained from `VecGetArrays()`

2042:    Notes:
2043:    For regular PETSc vectors this routine does not involve any copies. For
2044:    any special vectors that do not store local vector data in a contiguous
2045:    array, this routine will copy the data back into the underlying
2046:    vector data structure from the arrays obtained with `VecGetArrays()`.

2048:    Level: intermediate

2050: .seealso: [](chapter_vectors), `Vec`, `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

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: [](chapter_vectors), `Vec`, `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

2111:    Input Parameters:
2112: +  x - the vector
2113: -  a - location of pointer to array obtained from `VecGetArrayAndMemType()`

2115:    Level: beginner

2117: .seealso: [](chapter_vectors), `Vec`, `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 guaranteed 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: [](chapter_vectors), `Vec`, `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: [](chapter_vectors), `Vec`, `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:    Note:
2220:    The array must be returned using a matching call to `VecRestoreArrayWriteAndMemType()`, where it will label the device memory as most recent.

2222: .seealso: [](chapter_vectors), `Vec`, `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: [](chapter_vectors), `Vec`, `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: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecReplaceArray()`, `VecResetArray()`
2296: @*/
2297: PetscErrorCode VecPlaceArray(Vec vec, const PetscScalar array[])
2298: {
2302:   PetscUseTypeMethod(vec, placearray, array);
2303:   PetscObjectStateIncrease((PetscObject)vec);
2304:   return 0;
2305: }

2307: /*@C
2308:    VecReplaceArray - Allows one to replace the array in a vector with an
2309:    array provided by the user. This is useful to avoid copying an array
2310:    into a vector.

2312:    Not Collective; No Fortran Support

2314:    Input Parameters:
2315: +  vec - the vector
2316: -  array - the array

2318:    Level: developer

2320:    Notes:
2321:    This permanently replaces the array and frees the memory associated
2322:    with the old array.

2324:    The memory passed in MUST be obtained with `PetscMalloc()` and CANNOT be
2325:    freed by the user. It will be freed when the vector is destroyed.

2327: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecPlaceArray()`, `VecResetArray()`
2328: @*/
2329: PetscErrorCode VecReplaceArray(Vec vec, const PetscScalar array[])
2330: {
2333:   PetscUseTypeMethod(vec, replacearray, array);
2334:   PetscObjectStateIncrease((PetscObject)vec);
2335:   return 0;
2336: }

2338: /*@C
2339:    VecCUDAGetArray - Provides access to the CUDA buffer inside a vector.

2341:    Not Collective; No Fortran Support

2343:    Input Parameter:
2344: .  v - the vector

2346:    Output Parameter:
2347: .  a - the CUDA device pointer

2349:    Level: intermediate

2351:    Notes:
2352:    This function has semantics similar to `VecGetArray()`:  the pointer
2353:    returned by this function points to a consistent view of the vector
2354:    data.  This may involve a copy operation of data from the host to the
2355:    device if the data on the device is out of date.  If the device
2356:    memory hasn't been allocated previously it will be allocated as part
2357:    of this function call.  `VecCUDAGetArray()` assumes that
2358:    the user will modify the vector data.  This is similar to
2359:    intent(inout) in fortran.

2361:    The CUDA device pointer has to be released by calling
2362:    `VecCUDARestoreArray()`.  Upon restoring the vector data
2363:    the data on the host will be marked as out of date.  A subsequent
2364:    access of the host data will thus incur a data transfer from the
2365:    device to the host.

2367: .seealso: [](chapter_vectors), `Vec`, `VecCUDARestoreArray()`, `VecCUDAGetArrayRead()`, `VecCUDAGetArrayWrite()`, `VecGetArray()`, `VecGetArrayRead()`
2368: @*/
2369: PETSC_EXTERN PetscErrorCode VecCUDAGetArray(Vec v, PetscScalar **a)
2370: {
2372: #if defined(PETSC_HAVE_CUDA)
2373:   {
2374:     VecCUDACopyToGPU(v);
2375:     *a = ((Vec_CUDA *)v->spptr)->GPUarray;
2376:   }
2377: #endif
2378:   return 0;
2379: }

2381: /*@C
2382:    VecCUDARestoreArray - Restore a CUDA device pointer previously acquired with VecCUDAGetArray().

2384:    Not Collective; No Fortran Support

2386:    Input Parameters:
2387: +  v - the vector
2388: -  a - the CUDA device pointer.  This pointer is invalid after
2389:        `VecCUDARestoreArray()` returns.

2391:    Level: intermediate

2393:    Note:
2394:     This marks the host data as out of date.  Subsequent access to the
2395:    vector data on the host side with for instance `VecGetArray()` incurs a
2396:    data transfer.

2398: .seealso: [](chapter_vectors), `Vec`, `VecCUDAGetArray()`, `VecCUDAGetArrayRead()`, `VecCUDAGetArrayWrite()`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayRead()`
2399: @*/
2400: PETSC_EXTERN PetscErrorCode VecCUDARestoreArray(Vec v, PetscScalar **a)
2401: {
2403: #if defined(PETSC_HAVE_CUDA)
2404:   v->offloadmask = PETSC_OFFLOAD_GPU;
2405: #endif
2406:   PetscObjectStateIncrease((PetscObject)v);
2407:   return 0;
2408: }

2410: /*@C
2411:    VecCUDAGetArrayRead - Provides read access to the CUDA buffer inside a vector.

2413:    Not Collective; No Fortran Support

2415:    Input Parameter:
2416: .  v - the vector

2418:    Output Parameter:
2419: .  a - the CUDA pointer.

2421:    Level: intermediate

2423:    Notes:
2424:    This function is analogous to `VecGetArrayRead()`:  The pointer
2425:    returned by this function points to a consistent view of the vector
2426:    data.  This may involve a copy operation of data from the host to the
2427:    device if the data on the device is out of date.  If the device
2428:    memory hasn't been allocated previously it will be allocated as part
2429:    of this function call.  `VecCUDAGetArrayRead()` assumes that the
2430:    user will not modify the vector data.  This is analgogous to
2431:    intent(in) in Fortran.

2433:    The CUDA device pointer has to be released by calling
2434:    `VecCUDARestoreArrayRead()`.  If the data on the host side was
2435:    previously up to date it will remain so, i.e. data on both the device
2436:    and the host is up to date.  Accessing data on the host side does not
2437:    incur a device to host data transfer.

2439: .seealso: [](chapter_vectors), `Vec`, `VecCUDARestoreArrayRead()`, `VecCUDAGetArray()`, `VecCUDAGetArrayWrite()`, `VecGetArray()`, `VecGetArrayRead()`
2440: @*/
2441: PETSC_EXTERN PetscErrorCode VecCUDAGetArrayRead(Vec v, const PetscScalar **a)
2442: {
2443:   VecCUDAGetArray(v, (PetscScalar **)a);
2444:   return 0;
2445: }

2447: /*@C
2448:    VecCUDARestoreArrayRead - Restore a CUDA device pointer previously acquired with `VecCUDAGetArrayRead()`.

2450:    Not Collective; No Fortran Support

2452:    Input Parameters:
2453: +  v - the vector
2454: -  a - the CUDA device pointer.  This pointer is invalid after `VecCUDARestoreArrayRead()` returns.

2456:    Level: intermediate

2458:    Note:
2459:    If the data on the host side was previously up to date it will remain
2460:    so, i.e. data on both the device and the host is up to date.
2461:    Accessing data on the host side e.g. with `VecGetArray()` does not
2462:    incur a device to host data transfer.

2464: .seealso: [](chapter_vectors), `Vec`, `VecCUDAGetArrayRead()`, `VecCUDAGetArrayWrite()`, `VecCUDAGetArray()`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayRead()`
2465: @*/
2466: PETSC_EXTERN PetscErrorCode VecCUDARestoreArrayRead(Vec v, const PetscScalar **a)
2467: {
2469:   *a = NULL;
2470:   return 0;
2471: }

2473: /*@C
2474:    VecCUDAGetArrayWrite - Provides write access to the CUDA buffer inside a vector.

2476:    Not Collective; No Fortran Support

2478:    Input Parameter:
2479: .  v - the vector

2481:    Output Parameter:
2482: .  a - the CUDA pointer

2484:    Level: advanced

2486:    Notes:
2487:    The data pointed to by the device pointer is uninitialized.  The user
2488:    may not read from this data.  Furthermore, the entire array needs to
2489:    be filled by the user to obtain well-defined behaviour.  The device
2490:    memory will be allocated by this function if it hasn't been allocated
2491:    previously.  This is analogous to intent(out) in Fortran.

2493:    The device pointer needs to be released with
2494:    `VecCUDARestoreArrayWrite()`.  When the pointer is released the
2495:    host data of the vector is marked as out of data.  Subsequent access
2496:    of the host data with e.g. `VecGetArray()` incurs a device to host data
2497:    transfer.

2499: .seealso: [](chapter_vectors), `Vec`, `VecCUDARestoreArrayWrite()`, `VecCUDAGetArray()`, `VecCUDAGetArrayRead()`, `VecCUDAGetArrayWrite()`, `VecGetArray()`, `VecGetArrayRead()`
2500: @*/
2501: PETSC_EXTERN PetscErrorCode VecCUDAGetArrayWrite(Vec v, PetscScalar **a)
2502: {
2504: #if defined(PETSC_HAVE_CUDA)
2505:   {
2506:     VecCUDAAllocateCheck(v);
2507:     *a = ((Vec_CUDA *)v->spptr)->GPUarray;
2508:   }
2509: #endif
2510:   return 0;
2511: }

2513: /*@C
2514:    VecCUDARestoreArrayWrite - Restore a CUDA device pointer previously acquired with `VecCUDAGetArrayWrite()`.

2516:    Not Collective; No Fortran Support

2518:    Input Parameters:
2519: +  v - the vector
2520: -  a - the CUDA device pointer.  This pointer is invalid after
2521:        `VecCUDARestoreArrayWrite()` returns.

2523:    Level: intermediate

2525:    Note:
2526:     Data on the host will be marked as out of date.  Subsequent access of
2527:    the data on the host side e.g. with `VecGetArray()` will incur a device
2528:    to host data transfer.

2530: .seealso: [](chapter_vectors), `Vec`, `VecCUDAGetArrayWrite()`, `VecCUDAGetArray()`, `VecCUDAGetArrayRead()`, `VecCUDAGetArrayWrite()`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayRead()`
2531: @*/
2532: PETSC_EXTERN PetscErrorCode VecCUDARestoreArrayWrite(Vec v, PetscScalar **a)
2533: {
2535: #if defined(PETSC_HAVE_CUDA)
2536:   v->offloadmask = PETSC_OFFLOAD_GPU;
2537:   if (a) *a = NULL;
2538: #endif
2539:   PetscObjectStateIncrease((PetscObject)v);
2540:   return 0;
2541: }

2543: /*@C
2544:    VecCUDAPlaceArray - Allows one to replace the GPU array in a vector with a
2545:    GPU array provided by the user. This is useful to avoid copying an
2546:    array into a vector.

2548:    Not Collective; No Fortran Support

2550:    Input Parameters:
2551: +  vec - the vector
2552: -  array - the GPU array

2554:    Level: developer

2556:    Notes:
2557:    You can return to the original GPU array with a call to `VecCUDAResetArray()`
2558:    It is not possible to use `VecCUDAPlaceArray()` and `VecPlaceArray()` at the
2559:    same time on the same vector.

2561:    `vec` does not take ownership of `array` in any way. The user must free `array` themselves
2562:    but be careful not to do so before the vector has either been destroyed, had its original
2563:    array restored with `VecCUDAResetArray()` or permanently replaced with
2564:    `VecCUDAReplaceArray()`.

2566: .seealso: [](chapter_vectors), `Vec`, `VecPlaceArray()`, `VecGetArray()`, `VecRestoreArray()`, `VecReplaceArray()`, `VecResetArray()`, `VecCUDAResetArray()`, `VecCUDAReplaceArray()`
2567: @*/
2568: PetscErrorCode VecCUDAPlaceArray(Vec vin, const PetscScalar a[])
2569: {
2571: #if defined(PETSC_HAVE_CUDA)
2572:   VecCUDACopyToGPU(vin);
2574:   ((Vec_Seq *)vin->data)->unplacedarray = (PetscScalar *)((Vec_CUDA *)vin->spptr)->GPUarray; /* save previous GPU array so reset can bring it back */
2575:   ((Vec_CUDA *)vin->spptr)->GPUarray    = (PetscScalar *)a;
2576:   vin->offloadmask                      = PETSC_OFFLOAD_GPU;
2577: #endif
2578:   PetscObjectStateIncrease((PetscObject)vin);
2579:   return 0;
2580: }

2582: /*@C
2583:    VecCUDAReplaceArray - Allows one to replace the GPU array in a vector
2584:    with a GPU array provided by the user. This is useful to avoid copying
2585:    a GPU array into a vector.

2587:    Not Collective; No Fortran Support

2589:    Input Parameters:
2590: +  vec - the vector
2591: -  array - the GPU array

2593:    Level: developer

2595:    Notes:
2596:    This permanently replaces the GPU array and frees the memory associated
2597:    with the old GPU array.

2599:    The memory passed in CANNOT be freed by the user. It will be freed
2600:    when the vector is destroyed.

2602: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecPlaceArray()`, `VecResetArray()`, `VecCUDAResetArray()`, `VecCUDAPlaceArray()`, `VecReplaceArray()`
2603: @*/
2604: PetscErrorCode VecCUDAReplaceArray(Vec vin, const PetscScalar a[])
2605: {
2606: #if defined(PETSC_HAVE_CUDA)
2607: #endif

2610: #if defined(PETSC_HAVE_CUDA)
2611:   if (((Vec_CUDA *)vin->spptr)->GPUarray_allocated) cudaFree(((Vec_CUDA *)vin->spptr)->GPUarray_allocated);
2612:   ((Vec_CUDA *)vin->spptr)->GPUarray_allocated = ((Vec_CUDA *)vin->spptr)->GPUarray = (PetscScalar *)a;
2613:   vin->offloadmask                                                                  = PETSC_OFFLOAD_GPU;
2614: #endif
2615:   PetscObjectStateIncrease((PetscObject)vin);
2616:   return 0;
2617: }

2619: /*@C
2620:    VecCUDAResetArray - Resets a vector to use its default memory. Call this
2621:    after the use of `VecCUDAPlaceArray()`.

2623:    Not Collective; No Fortran Support

2625:    Input Parameters:
2626: .  vec - the vector

2628:    Level: developer

2630: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecReplaceArray()`, `VecPlaceArray()`, `VecResetArray()`, `VecCUDAPlaceArray()`, `VecCUDAReplaceArray()`
2631: @*/
2632: PetscErrorCode VecCUDAResetArray(Vec vin)
2633: {
2635: #if defined(PETSC_HAVE_CUDA)
2636:   VecCUDACopyToGPU(vin);
2637:   ((Vec_CUDA *)vin->spptr)->GPUarray    = (PetscScalar *)((Vec_Seq *)vin->data)->unplacedarray;
2638:   ((Vec_Seq *)vin->data)->unplacedarray = 0;
2639:   vin->offloadmask                      = PETSC_OFFLOAD_GPU;
2640: #endif
2641:   PetscObjectStateIncrease((PetscObject)vin);
2642:   return 0;
2643: }

2645: /*@C
2646:    VecHIPGetArray - Provides access to the HIP buffer inside a vector.

2648:    Not Collective; No Fortran Support

2650:    Input Parameter:
2651: .  v - the vector

2653:    Output Parameter:
2654: .  a - the HIP device pointer

2656:    Fortran note:
2657:    This function is not currently available from Fortran.

2659:    Level: intermediate

2661:    Notes:
2662:    This function has semantics similar to `VecGetArray()`:  the pointer
2663:    returned by this function points to a consistent view of the vector
2664:    data.  This may involve a copy operation of data from the host to the
2665:    device if the data on the device is out of date.  If the device
2666:    memory hasn't been allocated previously it will be allocated as part
2667:    of this function call.  `VecHIPGetArray()` assumes that
2668:    the user will modify the vector data.  This is similar to
2669:    intent(inout) in fortran.

2671:    The HIP device pointer has to be released by calling
2672:    `VecHIPRestoreArray()`.  Upon restoring the vector data
2673:    the data on the host will be marked as out of date.  A subsequent
2674:    access of the host data will thus incur a data transfer from the
2675:    device to the host.

2677: .seealso: [](chapter_vectors), `Vec`, `VecHIPRestoreArray()`, `VecHIPGetArrayRead()`, `VecHIPGetArrayWrite()`, `VecGetArray()`, `VecGetArrayRead()`
2678: @*/
2679: PETSC_EXTERN PetscErrorCode VecHIPGetArray(Vec v, PetscScalar **a)
2680: {
2682: #if defined(PETSC_HAVE_HIP)
2683:   *a = 0;
2684:   VecHIPCopyToGPU(v);
2685:   *a = ((Vec_HIP *)v->spptr)->GPUarray;
2686: #endif
2687:   return 0;
2688: }

2690: /*@C
2691:    VecHIPRestoreArray - Restore a HIP device pointer previously acquired with `VecHIPGetArray()`.

2693:    Not Collective; No Fortran Support

2695:    Input Parameters:
2696: +  v - the vector
2697: -  a - the HIP device pointer.  This pointer is invalid after
2698:        `VecHIPRestoreArray()` returns.

2700:    Level: intermediate

2702:    Note:
2703:    This marks the host data as out of date.  Subsequent access to the
2704:    vector data on the host side with for instance `VecGetArray()` incurs a
2705:    data transfer.

2707: .seealso: [](chapter_vectors), `Vec`, `VecHIPGetArray()`, `VecHIPGetArrayRead()`, `VecHIPGetArrayWrite()`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayRead()`
2708: @*/
2709: PETSC_EXTERN PetscErrorCode VecHIPRestoreArray(Vec v, PetscScalar **a)
2710: {
2712: #if defined(PETSC_HAVE_HIP)
2713:   v->offloadmask = PETSC_OFFLOAD_GPU;
2714: #endif

2716:   PetscObjectStateIncrease((PetscObject)v);
2717:   return 0;
2718: }

2720: /*@C
2721:    VecHIPGetArrayRead - Provides read access to the HIP buffer inside a vector.

2723:    Not Collective; No Fortran Support

2725:    Input Parameter:
2726: .  v - the vector

2728:    Output Parameter:
2729: .  a - the HIP pointer.

2731:    Level: intermediate

2733:    This function is analogous to `VecGetArrayRead()`:  The pointer
2734:    returned by this function points to a consistent view of the vector
2735:    data.  This may involve a copy operation of data from the host to the
2736:    device if the data on the device is out of date.  If the device
2737:    memory hasn't been allocated previously it will be allocated as part
2738:    of this function call.  `VecHIPGetArrayRead()` assumes that the
2739:    user will not modify the vector data.  This is analgogous to
2740:    intent(in) in Fortran.

2742:    The HIP device pointer has to be released by calling
2743:    `VecHIPRestoreArrayRead()`.  If the data on the host side was
2744:    previously up to date it will remain so, i.e. data on both the device
2745:    and the host is up to date.  Accessing data on the host side does not
2746:    incur a device to host data transfer.

2748: .seealso: [](chapter_vectors), `Vec`, `VecHIPRestoreArrayRead()`, `VecHIPGetArray()`, `VecHIPGetArrayWrite()`, `VecGetArray()`, `VecGetArrayRead()`
2749: @*/
2750: PETSC_EXTERN PetscErrorCode VecHIPGetArrayRead(Vec v, const PetscScalar **a)
2751: {
2753: #if defined(PETSC_HAVE_HIP)
2754:   *a = 0;
2755:   VecHIPCopyToGPU(v);
2756:   *a = ((Vec_HIP *)v->spptr)->GPUarray;
2757: #endif
2758:   return 0;
2759: }

2761: /*@C
2762:    VecHIPRestoreArrayRead - Restore a HIP device pointer previously acquired with `VecHIPGetArrayRead()`.

2764:    Not Collective; No Fortran Support

2766:    Input Parameters:
2767: +  v - the vector
2768: -  a - the HIP device pointer.  This pointer is invalid after
2769:        `VecHIPRestoreArrayRead()` returns.

2771:    Level: intermediate

2773:    Note:
2774:    If the data on the host side was previously up to date it will remain
2775:    so, i.e. data on both the device and the host is up to date.
2776:    Accessing data on the host side e.g. with `VecGetArray()` does not
2777:    incur a device to host data transfer.

2779: .seealso: [](chapter_vectors), `Vec`, `VecHIPGetArrayRead()`, `VecHIPGetArrayWrite()`, `VecHIPGetArray()`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayRead()`
2780: @*/
2781: PETSC_EXTERN PetscErrorCode VecHIPRestoreArrayRead(Vec v, const PetscScalar **a)
2782: {
2784:   *a = NULL;
2785:   return 0;
2786: }

2788: /*@C
2789:    VecHIPGetArrayWrite - Provides write access to the HIP buffer inside a vector.

2791:    Not Collective; No Fortran Support

2793:    Input Parameter:
2794: .  v - the vector

2796:    Output Parameter:
2797: .  a - the HIP pointer

2799:    Level: advanced

2801:    Notes:
2802:    The data pointed to by the device pointer is uninitialized.  The user
2803:    may not read from this data.  Furthermore, the entire array needs to
2804:    be filled by the user to obtain well-defined behaviour.  The device
2805:    memory will be allocated by this function if it hasn't been allocated
2806:    previously.  This is analogous to intent(out) in Fortran.

2808:    The device pointer needs to be released with
2809:    `VecHIPRestoreArrayWrite()`.  When the pointer is released the
2810:    host data of the vector is marked as out of data.  Subsequent access
2811:    of the host data with e.g. `VecGetArray()` incurs a device to host data
2812:    transfer.

2814: .seealso: [](chapter_vectors), `Vec`, `VecHIPRestoreArrayWrite()`, `VecHIPGetArray()`, `VecHIPGetArrayRead()`, `VecHIPGetArrayWrite()`, `VecGetArray()`, `VecGetArrayRead()`
2815: @*/
2816: PETSC_EXTERN PetscErrorCode VecHIPGetArrayWrite(Vec v, PetscScalar **a)
2817: {
2819: #if defined(PETSC_HAVE_HIP)
2820:   *a = 0;
2821:   VecHIPAllocateCheck(v);
2822:   *a = ((Vec_HIP *)v->spptr)->GPUarray;
2823: #endif
2824:   return 0;
2825: }

2827: /*@C
2828:    VecHIPRestoreArrayWrite - Restore a HIP device pointer previously acquired with `VecHIPGetArrayWrite()`.

2830:    Not Collective; No Fortran Support

2832:    Input Parameters:
2833: +  v - the vector
2834: -  a - the HIP device pointer.  This pointer is invalid after
2835:        `VecHIPRestoreArrayWrite()` returns.

2837:    Level: intermediate

2839:    Note:
2840:    Data on the host will be marked as out of date.  Subsequent access of
2841:    the data on the host side e.g. with `VecGetArray()` will incur a device
2842:    to host data transfer.

2844: .seealso: [](chapter_vectors), `Vec`, `VecHIPGetArrayWrite()`, `VecHIPGetArray()`, `VecHIPGetArrayRead()`, `VecHIPGetArrayWrite()`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayRead()`
2845: @*/
2846: PETSC_EXTERN PetscErrorCode VecHIPRestoreArrayWrite(Vec v, PetscScalar **a)
2847: {
2849: #if defined(PETSC_HAVE_HIP)
2850:   v->offloadmask = PETSC_OFFLOAD_GPU;
2851: #endif

2853:   PetscObjectStateIncrease((PetscObject)v);
2854:   return 0;
2855: }

2857: /*@C
2858:    VecHIPPlaceArray - Allows one to replace the GPU array in a vector with a
2859:    GPU array provided by the user. This is useful to avoid copying an
2860:    array into a vector.

2862:    Not Collective; No Fortran Support

2864:    Input Parameters:
2865: +  vec - the vector
2866: -  array - the GPU array

2868:    Level: developer

2870:    Notes:
2871:    You can return to the original GPU array with a call to `VecHIPResetArray()`
2872:    It is not possible to use `VecHIPPlaceArray()` and `VecPlaceArray()` at the
2873:    same time on the same vector.

2875:    `vec` does not take ownership of `array` in any way. The user must free `array` themselves
2876:    but be careful not to do so before the vector has either been destroyed, had its original
2877:    array restored with `VecHIPResetArray()` or permanently replaced with
2878:    `VecHIPReplaceArray()`.

2880: .seealso: [](chapter_vectors), `Vec`, `VecPlaceArray()`, `VecGetArray()`, `VecRestoreArray()`, `VecReplaceArray()`, `VecResetArray()`, `VecHIPResetArray()`, `VecHIPReplaceArray()`
2881: @*/
2882: PetscErrorCode VecHIPPlaceArray(Vec vin, const PetscScalar a[])
2883: {
2885: #if defined(PETSC_HAVE_HIP)
2886:   VecHIPCopyToGPU(vin);
2888:   ((Vec_Seq *)vin->data)->unplacedarray = (PetscScalar *)((Vec_HIP *)vin->spptr)->GPUarray; /* save previous GPU array so reset can bring it back */
2889:   ((Vec_HIP *)vin->spptr)->GPUarray     = (PetscScalar *)a;
2890:   vin->offloadmask                      = PETSC_OFFLOAD_GPU;
2891: #endif
2892:   PetscObjectStateIncrease((PetscObject)vin);
2893:   return 0;
2894: }

2896: /*@C
2897:    VecHIPReplaceArray - Allows one to replace the GPU array in a vector
2898:    with a GPU array provided by the user. This is useful to avoid copying
2899:    a GPU array into a vector.

2901:    Not Collective; No Fortran Support

2903:    Input Parameters:
2904: +  vec - the vector
2905: -  array - the GPU array

2907:    Level: developer

2909:    Notes:
2910:    This permanently replaces the GPU array and frees the memory associated
2911:    with the old GPU array.

2913:    The memory passed in CANNOT be freed by the user. It will be freed
2914:    when the vector is destroyed.

2916: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecPlaceArray()`, `VecResetArray()`, `VecHIPResetArray()`, `VecHIPPlaceArray()`, `VecReplaceArray()`
2917: @*/
2918: PetscErrorCode VecHIPReplaceArray(Vec vin, const PetscScalar a[])
2919: {
2921: #if defined(PETSC_HAVE_HIP)
2922:   hipFree(((Vec_HIP *)vin->spptr)->GPUarray);
2923:   ((Vec_HIP *)vin->spptr)->GPUarray = (PetscScalar *)a;
2924:   vin->offloadmask                  = PETSC_OFFLOAD_GPU;
2925: #endif
2926:   PetscObjectStateIncrease((PetscObject)vin);
2927:   return 0;
2928: }

2930: /*@C
2931:    VecHIPResetArray - Resets a vector to use its default memory. Call this
2932:    after the use of `VecHIPPlaceArray()`.

2934:    Not Collective; No Fortran Support

2936:    Input Parameters:
2937: .  vec - the vector

2939:    Level: developer

2941: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecReplaceArray()`, `VecPlaceArray()`, `VecResetArray()`, `VecHIPPlaceArray()`, `VecHIPReplaceArray()`
2942: @*/
2943: PetscErrorCode VecHIPResetArray(Vec vin)
2944: {
2946: #if defined(PETSC_HAVE_HIP)
2947:   VecHIPCopyToGPU(vin);
2948:   ((Vec_HIP *)vin->spptr)->GPUarray     = (PetscScalar *)((Vec_Seq *)vin->data)->unplacedarray;
2949:   ((Vec_Seq *)vin->data)->unplacedarray = 0;
2950:   vin->offloadmask                      = PETSC_OFFLOAD_GPU;
2951: #endif
2952:   PetscObjectStateIncrease((PetscObject)vin);
2953:   return 0;
2954: }

2956: /*MC
2957:     VecDuplicateVecsF90 - Creates several vectors of the same type as an existing vector
2958:     and makes them accessible via a Fortran pointer.

2960:     Synopsis:
2961:     VecDuplicateVecsF90(Vec x,PetscInt n,{Vec, pointer :: y(:)},integer ierr)

2963:     Collective

2965:     Input Parameters:
2966: +   x - a vector to mimic
2967: -   n - the number of vectors to obtain

2969:     Output Parameters:
2970: +   y - Fortran pointer to the array of vectors
2971: -   ierr - error code

2973:     Example of Usage:
2974: .vb
2975: #include <petsc/finclude/petscvec.h>
2976:     use petscvec

2978:     Vec x
2979:     Vec, pointer :: y(:)
2980:     ....
2981:     call VecDuplicateVecsF90(x,2,y,ierr)
2982:     call VecSet(y(2),alpha,ierr)
2983:     call VecSet(y(2),alpha,ierr)
2984:     ....
2985:     call VecDestroyVecsF90(2,y,ierr)
2986: .ve

2988:     Level: beginner

2990:     Note:
2991:     Use `VecDestroyVecsF90()` to free the space.

2993: .seealso: [](chapter_vectors), `Vec`, `VecDestroyVecsF90()`, `VecDuplicateVecs()`
2994: M*/

2996: /*MC
2997:     VecRestoreArrayF90 - Restores a vector to a usable state after a call to
2998:     `VecGetArrayF90()`.

3000:     Synopsis:
3001:     VecRestoreArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)

3003:     Logically Collective

3005:     Input Parameters:
3006: +   x - vector
3007: -   xx_v - the Fortran pointer to the array

3009:     Output Parameter:
3010: .   ierr - error code

3012:     Example of Usage:
3013: .vb
3014: #include <petsc/finclude/petscvec.h>
3015:     use petscvec

3017:     PetscScalar, pointer :: xx_v(:)
3018:     ....
3019:     call VecGetArrayF90(x,xx_v,ierr)
3020:     xx_v(3) = a
3021:     call VecRestoreArrayF90(x,xx_v,ierr)
3022: .ve

3024:     Level: beginner

3026: .seealso: [](chapter_vectors), `Vec`, `VecGetArrayF90()`, `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrayReadF90()`
3027: M*/

3029: /*MC
3030:     VecDestroyVecsF90 - Frees a block of vectors obtained with `VecDuplicateVecsF90()`.

3032:     Synopsis:
3033:     VecDestroyVecsF90(PetscInt n,{Vec, pointer :: x(:)},PetscErrorCode ierr)

3035:     Collective

3037:     Input Parameters:
3038: +   n - the number of vectors previously obtained
3039: -   x - pointer to array of vector pointers

3041:     Output Parameter:
3042: .   ierr - error code

3044:     Level: beginner

3046: .seealso: [](chapter_vectors), `Vec`, `VecDestroyVecs()`, `VecDuplicateVecsF90()`
3047: M*/

3049: /*MC
3050:     VecGetArrayF90 - Accesses a vector array from Fortran. For default PETSc
3051:     vectors, `VecGetArrayF90()` returns a pointer to the local data array. Otherwise,
3052:     this routine is implementation dependent. You MUST call `VecRestoreArrayF90()`
3053:     when you no longer need access to the array.

3055:     Synopsis:
3056:     VecGetArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)

3058:     Logically Collective

3060:     Input Parameter:
3061: .   x - vector

3063:     Output Parameters:
3064: +   xx_v - the Fortran pointer to the array
3065: -   ierr - error code

3067:     Example of Usage:
3068: .vb
3069: #include <petsc/finclude/petscvec.h>
3070:     use petscvec

3072:     PetscScalar, pointer :: xx_v(:)
3073:     ....
3074:     call VecGetArrayF90(x,xx_v,ierr)
3075:     xx_v(3) = a
3076:     call VecRestoreArrayF90(x,xx_v,ierr)
3077: .ve

3079:      Level: beginner

3081:     Note:
3082:     If you ONLY intend to read entries from the array and not change any entries you should use `VecGetArrayReadF90()`.

3084: .seealso: [](chapter_vectors), `Vec`, `VecRestoreArrayF90()`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayReadF90()`
3085: M*/

3087: /*MC
3088:     VecGetArrayReadF90 - Accesses a read only array from Fortran. For default PETSc
3089:     vectors, `VecGetArrayF90()` returns a pointer to the local data array. Otherwise,
3090:     this routine is implementation dependent. You MUST call `VecRestoreArrayReadF90()`
3091:     when you no longer need access to the array.

3093:     Synopsis:
3094:     VecGetArrayReadF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)

3096:     Logically Collective

3098:     Input Parameter:
3099: .   x - vector

3101:     Output Parameters:
3102: +   xx_v - the Fortran pointer to the array
3103: -   ierr - error code

3105:     Example of Usage:
3106: .vb
3107: #include <petsc/finclude/petscvec.h>
3108:     use petscvec

3110:     PetscScalar, pointer :: xx_v(:)
3111:     ....
3112:     call VecGetArrayReadF90(x,xx_v,ierr)
3113:     a = xx_v(3)
3114:     call VecRestoreArrayReadF90(x,xx_v,ierr)
3115: .ve

3117:     Level: beginner

3119:     Note:
3120:     If you intend to write entries into the array you must use `VecGetArrayF90()`.

3122: .seealso: [](chapter_vectors), `Vec`, `VecRestoreArrayReadF90()`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayRead()`, `VecRestoreArrayRead()`, `VecGetArrayF90()`
3123: M*/

3125: /*MC
3126:     VecRestoreArrayReadF90 - Restores a readonly vector to a usable state after a call to
3127:     `VecGetArrayReadF90()`.

3129:     Synopsis:
3130:     VecRestoreArrayReadF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)

3132:     Logically Collective

3134:     Input Parameters:
3135: +   x - vector
3136: -   xx_v - the Fortran pointer to the array

3138:     Output Parameter:
3139: .   ierr - error code

3141:     Example of Usage:
3142: .vb
3143: #include <petsc/finclude/petscvec.h>
3144:     use petscvec

3146:     PetscScalar, pointer :: xx_v(:)
3147:     ....
3148:     call VecGetArrayReadF90(x,xx_v,ierr)
3149:     a = xx_v(3)
3150:     call VecRestoreArrayReadF90(x,xx_v,ierr)
3151: .ve

3153:     Level: beginner

3155: .seealso: [](chapter_vectors), `Vec`, `VecGetArrayReadF90()`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrayRead()`, `VecRestoreArrayRead()`, `VecRestoreArrayF90()`
3156: M*/

3158: /*@C
3159:    VecGetArray2d - Returns a pointer to a 2d contiguous array that contains this
3160:    processor's portion of the vector data.  You MUST call `VecRestoreArray2d()`
3161:    when you no longer need access to the array.

3163:    Logically Collective

3165:    Input Parameters:
3166: +  x - the vector
3167: .  m - first dimension of two dimensional array
3168: .  n - second dimension of two dimensional array
3169: .  mstart - first index you will use in first coordinate direction (often 0)
3170: -  nstart - first index in the second coordinate direction (often 0)

3172:    Output Parameter:
3173: .  a - location to put pointer to the array

3175:    Level: developer

3177:   Notes:
3178:    For a vector obtained from `DMCreateLocalVector()` mstart and nstart are likely
3179:    obtained from the corner indices obtained from `DMDAGetGhostCorners()` while for
3180:    `DMCreateGlobalVector()` they are the corner indices from `DMDAGetCorners()`. In both cases
3181:    the arguments from `DMDAGet[Ghost]Corners()` are reversed in the call to `VecGetArray2d()`.

3183:    For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.

3185: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
3186:           `VecRestoreArray2d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
3187:           `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3188: @*/
3189: PetscErrorCode VecGetArray2d(Vec x, PetscInt m, PetscInt n, PetscInt mstart, PetscInt nstart, PetscScalar **a[])
3190: {
3191:   PetscInt     i, N;
3192:   PetscScalar *aa;

3197:   VecGetLocalSize(x, &N);
3199:   VecGetArray(x, &aa);

3201:   PetscMalloc1(m, a);
3202:   for (i = 0; i < m; i++) (*a)[i] = aa + i * n - nstart;
3203:   *a -= mstart;
3204:   return 0;
3205: }

3207: /*@C
3208:    VecGetArray2dWrite - Returns a pointer to a 2d contiguous array that will contain this
3209:    processor's portion of the vector data.  You MUST call `VecRestoreArray2dWrite()`
3210:    when you no longer need access to the array.

3212:    Logically Collective

3214:    Input Parameters:
3215: +  x - the vector
3216: .  m - first dimension of two dimensional array
3217: .  n - second dimension of two dimensional array
3218: .  mstart - first index you will use in first coordinate direction (often 0)
3219: -  nstart - first index in the second coordinate direction (often 0)

3221:    Output Parameter:
3222: .  a - location to put pointer to the array

3224:    Level: developer

3226:   Notes:
3227:    For a vector obtained from `DMCreateLocalVector()` mstart and nstart are likely
3228:    obtained from the corner indices obtained from `DMDAGetGhostCorners()` while for
3229:    `DMCreateGlobalVector()` they are the corner indices from `DMDAGetCorners()`. In both cases
3230:    the arguments from `DMDAGet[Ghost]Corners()` are reversed in the call to `VecGetArray2d()`.

3232:    For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.

3234: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
3235:           `VecRestoreArray2d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
3236:           `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3237: @*/
3238: PetscErrorCode VecGetArray2dWrite(Vec x, PetscInt m, PetscInt n, PetscInt mstart, PetscInt nstart, PetscScalar **a[])
3239: {
3240:   PetscInt     i, N;
3241:   PetscScalar *aa;

3246:   VecGetLocalSize(x, &N);
3248:   VecGetArrayWrite(x, &aa);

3250:   PetscMalloc1(m, a);
3251:   for (i = 0; i < m; i++) (*a)[i] = aa + i * n - nstart;
3252:   *a -= mstart;
3253:   return 0;
3254: }

3256: /*@C
3257:    VecRestoreArray2d - Restores a vector after `VecGetArray2d()` has been called.

3259:    Logically Collective

3261:    Input Parameters:
3262: +  x - the vector
3263: .  m - first dimension of two dimensional array
3264: .  n - second dimension of the two dimensional array
3265: .  mstart - first index you will use in first coordinate direction (often 0)
3266: .  nstart - first index in the second coordinate direction (often 0)
3267: -  a - location of pointer to array obtained from `VecGetArray2d()`

3269:    Level: developer

3271:    Notes:
3272:    For regular PETSc vectors this routine does not involve any copies. For
3273:    any special vectors that do not store local vector data in a contiguous
3274:    array, this routine will copy the data back into the underlying
3275:    vector data structure from the array obtained with `VecGetArray()`.

3277:    This routine actually zeros out the a pointer.

3279: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
3280:           `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
3281:           `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3282: @*/
3283: PetscErrorCode VecRestoreArray2d(Vec x, PetscInt m, PetscInt n, PetscInt mstart, PetscInt nstart, PetscScalar **a[])
3284: {
3285:   void *dummy;

3290:   dummy = (void *)(*a + mstart);
3291:   PetscFree(dummy);
3292:   VecRestoreArray(x, NULL);
3293:   return 0;
3294: }

3296: /*@C
3297:    VecRestoreArray2dWrite - Restores a vector after VecGetArray2dWrite`()` has been called.

3299:    Logically Collective

3301:    Input Parameters:
3302: +  x - the vector
3303: .  m - first dimension of two dimensional array
3304: .  n - second dimension of the two dimensional array
3305: .  mstart - first index you will use in first coordinate direction (often 0)
3306: .  nstart - first index in the second coordinate direction (often 0)
3307: -  a - location of pointer to array obtained from `VecGetArray2d()`

3309:    Level: developer

3311:    Notes:
3312:    For regular PETSc vectors this routine does not involve any copies. For
3313:    any special vectors that do not store local vector data in a contiguous
3314:    array, this routine will copy the data back into the underlying
3315:    vector data structure from the array obtained with `VecGetArray()`.

3317:    This routine actually zeros out the a pointer.

3319: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
3320:           `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
3321:           `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3322: @*/
3323: PetscErrorCode VecRestoreArray2dWrite(Vec x, PetscInt m, PetscInt n, PetscInt mstart, PetscInt nstart, PetscScalar **a[])
3324: {
3325:   void *dummy;

3330:   dummy = (void *)(*a + mstart);
3331:   PetscFree(dummy);
3332:   VecRestoreArrayWrite(x, NULL);
3333:   return 0;
3334: }

3336: /*@C
3337:    VecGetArray1d - Returns a pointer to a 1d contiguous array that contains this
3338:    processor's portion of the vector data.  You MUST call `VecRestoreArray1d()`
3339:    when you no longer need access to the array.

3341:    Logically Collective

3343:    Input Parameters:
3344: +  x - the vector
3345: .  m - first dimension of two dimensional array
3346: -  mstart - first index you will use in first coordinate direction (often 0)

3348:    Output Parameter:
3349: .  a - location to put pointer to the array

3351:    Level: developer

3353:   Notes:
3354:    For a vector obtained from `DMCreateLocalVector()` mstart are likely
3355:    obtained from the corner indices obtained from `DMDAGetGhostCorners()` while for
3356:    `DMCreateGlobalVector()` they are the corner indices from `DMDAGetCorners()`.

3358:    For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.

3360: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
3361:           `VecRestoreArray2d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
3362:           `VecGetArray2d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3363: @*/
3364: PetscErrorCode VecGetArray1d(Vec x, PetscInt m, PetscInt mstart, PetscScalar *a[])
3365: {
3366:   PetscInt N;

3371:   VecGetLocalSize(x, &N);
3373:   VecGetArray(x, a);
3374:   *a -= mstart;
3375:   return 0;
3376: }

3378: /*@C
3379:    VecGetArray1dWrite - Returns a pointer to a 1d contiguous array that will contain this
3380:    processor's portion of the vector data.  You MUST call `VecRestoreArray1dWrite()`
3381:    when you no longer need access to the array.

3383:    Logically Collective

3385:    Input Parameters:
3386: +  x - the vector
3387: .  m - first dimension of two dimensional array
3388: -  mstart - first index you will use in first coordinate direction (often 0)

3390:    Output Parameter:
3391: .  a - location to put pointer to the array

3393:    Level: developer

3395:   Notes:
3396:    For a vector obtained from `DMCreateLocalVector()` mstart are likely
3397:    obtained from the corner indices obtained from `DMDAGetGhostCorners()` while for
3398:    `DMCreateGlobalVector()` they are the corner indices from `DMDAGetCorners()`.

3400:    For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.

3402: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
3403:           `VecRestoreArray2d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
3404:           `VecGetArray2d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3405: @*/
3406: PetscErrorCode VecGetArray1dWrite(Vec x, PetscInt m, PetscInt mstart, PetscScalar *a[])
3407: {
3408:   PetscInt N;

3413:   VecGetLocalSize(x, &N);
3415:   VecGetArrayWrite(x, a);
3416:   *a -= mstart;
3417:   return 0;
3418: }

3420: /*@C
3421:    VecRestoreArray1d - Restores a vector after `VecGetArray1d()` has been called.

3423:    Logically Collective

3425:    Input Parameters:
3426: +  x - the vector
3427: .  m - first dimension of two dimensional array
3428: .  mstart - first index you will use in first coordinate direction (often 0)
3429: -  a - location of pointer to array obtained from `VecGetArray1d()`

3431:    Level: developer

3433:    Notes:
3434:    For regular PETSc vectors this routine does not involve any copies. For
3435:    any special vectors that do not store local vector data in a contiguous
3436:    array, this routine will copy the data back into the underlying
3437:    vector data structure from the array obtained with `VecGetArray1d()`.

3439:    This routine actually zeros out the a pointer.

3441: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
3442:           `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
3443:           `VecGetArray1d()`, `VecRestoreArray2d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3444: @*/
3445: PetscErrorCode VecRestoreArray1d(Vec x, PetscInt m, PetscInt mstart, PetscScalar *a[])
3446: {
3449:   VecRestoreArray(x, NULL);
3450:   return 0;
3451: }

3453: /*@C
3454:    VecRestoreArray1dWrite - Restores a vector after `VecGetArray1dWrite()` has been called.

3456:    Logically Collective

3458:    Input Parameters:
3459: +  x - the vector
3460: .  m - first dimension of two dimensional array
3461: .  mstart - first index you will use in first coordinate direction (often 0)
3462: -  a - location of pointer to array obtained from `VecGetArray1d()`

3464:    Level: developer

3466:    Notes:
3467:    For regular PETSc vectors this routine does not involve any copies. For
3468:    any special vectors that do not store local vector data in a contiguous
3469:    array, this routine will copy the data back into the underlying
3470:    vector data structure from the array obtained with `VecGetArray1d()`.

3472:    This routine actually zeros out the a pointer.

3474: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
3475:           `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
3476:           `VecGetArray1d()`, `VecRestoreArray2d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3477: @*/
3478: PetscErrorCode VecRestoreArray1dWrite(Vec x, PetscInt m, PetscInt mstart, PetscScalar *a[])
3479: {
3482:   VecRestoreArrayWrite(x, NULL);
3483:   return 0;
3484: }

3486: /*@C
3487:    VecGetArray3d - Returns a pointer to a 3d contiguous array that contains this
3488:    processor's portion of the vector data.  You MUST call `VecRestoreArray3d()`
3489:    when you no longer need access to the array.

3491:    Logically Collective

3493:    Input Parameters:
3494: +  x - the vector
3495: .  m - first dimension of three dimensional array
3496: .  n - second dimension of three dimensional array
3497: .  p - third dimension of three dimensional array
3498: .  mstart - first index you will use in first coordinate direction (often 0)
3499: .  nstart - first index in the second coordinate direction (often 0)
3500: -  pstart - first index in the third coordinate direction (often 0)

3502:    Output Parameter:
3503: .  a - location to put pointer to the array

3505:    Level: developer

3507:   Notes:
3508:    For a vector obtained from `DMCreateLocalVector()` mstart, nstart, and pstart are likely
3509:    obtained from the corner indices obtained from `DMDAGetGhostCorners()` while for
3510:    `DMCreateGlobalVector()` they are the corner indices from `DMDAGetCorners()`. In both cases
3511:    the arguments from `DMDAGet[Ghost]Corners()` are reversed in the call to `VecGetArray3d()`.

3513:    For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.

3515: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
3516:           `VecRestoreArray2d()`, `DMDAVecGetarray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
3517:           `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3518: @*/
3519: PetscErrorCode VecGetArray3d(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscScalar ***a[])
3520: {
3521:   PetscInt     i, N, j;
3522:   PetscScalar *aa, **b;

3527:   VecGetLocalSize(x, &N);
3529:   VecGetArray(x, &aa);

3531:   PetscMalloc(m * sizeof(PetscScalar **) + m * n * sizeof(PetscScalar *), a);
3532:   b = (PetscScalar **)((*a) + m);
3533:   for (i = 0; i < m; i++) (*a)[i] = b + i * n - nstart;
3534:   for (i = 0; i < m; i++)
3535:     for (j = 0; j < n; j++) b[i * n + j] = aa + i * n * p + j * p - pstart;
3536:   *a -= mstart;
3537:   return 0;
3538: }

3540: /*@C
3541:    VecGetArray3dWrite - Returns a pointer to a 3d contiguous array that will contain this
3542:    processor's portion of the vector data.  You MUST call `VecRestoreArray3dWrite()`
3543:    when you no longer need access to the array.

3545:    Logically Collective

3547:    Input Parameters:
3548: +  x - the vector
3549: .  m - first dimension of three dimensional array
3550: .  n - second dimension of three dimensional array
3551: .  p - third dimension of three dimensional array
3552: .  mstart - first index you will use in first coordinate direction (often 0)
3553: .  nstart - first index in the second coordinate direction (often 0)
3554: -  pstart - first index in the third coordinate direction (often 0)

3556:    Output Parameter:
3557: .  a - location to put pointer to the array

3559:    Level: developer

3561:   Notes:
3562:    For a vector obtained from `DMCreateLocalVector()` mstart, nstart, and pstart are likely
3563:    obtained from the corner indices obtained from `DMDAGetGhostCorners()` while for
3564:    `DMCreateGlobalVector()` they are the corner indices from `DMDAGetCorners()`. In both cases
3565:    the arguments from `DMDAGet[Ghost]Corners()` are reversed in the call to `VecGetArray3d()`.

3567:    For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.

3569: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
3570:           `VecRestoreArray2d()`, `DMDAVecGetarray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
3571:           `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3572: @*/
3573: PetscErrorCode VecGetArray3dWrite(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscScalar ***a[])
3574: {
3575:   PetscInt     i, N, j;
3576:   PetscScalar *aa, **b;

3581:   VecGetLocalSize(x, &N);
3583:   VecGetArrayWrite(x, &aa);

3585:   PetscMalloc(m * sizeof(PetscScalar **) + m * n * sizeof(PetscScalar *), a);
3586:   b = (PetscScalar **)((*a) + m);
3587:   for (i = 0; i < m; i++) (*a)[i] = b + i * n - nstart;
3588:   for (i = 0; i < m; i++)
3589:     for (j = 0; j < n; j++) b[i * n + j] = aa + i * n * p + j * p - pstart;

3591:   *a -= mstart;
3592:   return 0;
3593: }

3595: /*@C
3596:    VecRestoreArray3d - Restores a vector after `VecGetArray3d()` has been called.

3598:    Logically Collective

3600:    Input Parameters:
3601: +  x - the vector
3602: .  m - first dimension of three dimensional array
3603: .  n - second dimension of the three dimensional array
3604: .  p - third dimension of the three dimensional array
3605: .  mstart - first index you will use in first coordinate direction (often 0)
3606: .  nstart - first index in the second coordinate direction (often 0)
3607: .  pstart - first index in the third coordinate direction (often 0)
3608: -  a - location of pointer to array obtained from VecGetArray3d()

3610:    Level: developer

3612:    Notes:
3613:    For regular PETSc vectors this routine does not involve any copies. For
3614:    any special vectors that do not store local vector data in a contiguous
3615:    array, this routine will copy the data back into the underlying
3616:    vector data structure from the array obtained with `VecGetArray()`.

3618:    This routine actually zeros out the a pointer.

3620: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
3621:           `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
3622:           `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`, `VecGet`
3623: @*/
3624: PetscErrorCode VecRestoreArray3d(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscScalar ***a[])
3625: {
3626:   void *dummy;

3631:   dummy = (void *)(*a + mstart);
3632:   PetscFree(dummy);
3633:   VecRestoreArray(x, NULL);
3634:   return 0;
3635: }

3637: /*@C
3638:    VecRestoreArray3dWrite - Restores a vector after `VecGetArray3dWrite()` has been called.

3640:    Logically Collective

3642:    Input Parameters:
3643: +  x - the vector
3644: .  m - first dimension of three dimensional array
3645: .  n - second dimension of the three dimensional array
3646: .  p - third dimension of the three dimensional array
3647: .  mstart - first index you will use in first coordinate direction (often 0)
3648: .  nstart - first index in the second coordinate direction (often 0)
3649: .  pstart - first index in the third coordinate direction (often 0)
3650: -  a - location of pointer to array obtained from VecGetArray3d()

3652:    Level: developer

3654:    Notes:
3655:    For regular PETSc vectors this routine does not involve any copies. For
3656:    any special vectors that do not store local vector data in a contiguous
3657:    array, this routine will copy the data back into the underlying
3658:    vector data structure from the array obtained with `VecGetArray()`.

3660:    This routine actually zeros out the a pointer.

3662: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
3663:           `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
3664:           `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`, `VecGet`
3665: @*/
3666: PetscErrorCode VecRestoreArray3dWrite(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscScalar ***a[])
3667: {
3668:   void *dummy;

3673:   dummy = (void *)(*a + mstart);
3674:   PetscFree(dummy);
3675:   VecRestoreArrayWrite(x, NULL);
3676:   return 0;
3677: }

3679: /*@C
3680:    VecGetArray4d - Returns a pointer to a 4d contiguous array that contains this
3681:    processor's portion of the vector data.  You MUST call `VecRestoreArray4d()`
3682:    when you no longer need access to the array.

3684:    Logically Collective

3686:    Input Parameters:
3687: +  x - the vector
3688: .  m - first dimension of four dimensional array
3689: .  n - second dimension of four dimensional array
3690: .  p - third dimension of four dimensional array
3691: .  q - fourth dimension of four dimensional array
3692: .  mstart - first index you will use in first coordinate direction (often 0)
3693: .  nstart - first index in the second coordinate direction (often 0)
3694: .  pstart - first index in the third coordinate direction (often 0)
3695: -  qstart - first index in the fourth coordinate direction (often 0)

3697:    Output Parameter:
3698: .  a - location to put pointer to the array

3700:    Level: beginner

3702:   Notes:
3703:    For a vector obtained from `DMCreateLocalVector()` mstart, nstart, and pstart are likely
3704:    obtained from the corner indices obtained from `DMDAGetGhostCorners()` while for
3705:    `DMCreateGlobalVector()` they are the corner indices from `DMDAGetCorners()`. In both cases
3706:    the arguments from `DMDAGet[Ghost]Corners()` are reversed in the call to `VecGetArray3d()`.

3708:    For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.

3710: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
3711:           `VecRestoreArray2d()`, `DMDAVecGetarray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
3712:           `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3713: @*/
3714: PetscErrorCode VecGetArray4d(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt q, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscInt qstart, PetscScalar ****a[])
3715: {
3716:   PetscInt     i, N, j, k;
3717:   PetscScalar *aa, ***b, **c;

3722:   VecGetLocalSize(x, &N);
3724:   VecGetArray(x, &aa);

3726:   PetscMalloc(m * sizeof(PetscScalar ***) + m * n * sizeof(PetscScalar **) + m * n * p * sizeof(PetscScalar *), a);
3727:   b = (PetscScalar ***)((*a) + m);
3728:   c = (PetscScalar **)(b + m * n);
3729:   for (i = 0; i < m; i++) (*a)[i] = b + i * n - nstart;
3730:   for (i = 0; i < m; i++)
3731:     for (j = 0; j < n; j++) b[i * n + j] = c + i * n * p + j * p - pstart;
3732:   for (i = 0; i < m; i++)
3733:     for (j = 0; j < n; j++)
3734:       for (k = 0; k < p; k++) c[i * n * p + j * p + k] = aa + i * n * p * q + j * p * q + k * q - qstart;
3735:   *a -= mstart;
3736:   return 0;
3737: }

3739: /*@C
3740:    VecGetArray4dWrite - Returns a pointer to a 4d contiguous array that will contain this
3741:    processor's portion of the vector data.  You MUST call `VecRestoreArray4dWrite()`
3742:    when you no longer need access to the array.

3744:    Logically Collective

3746:    Input Parameters:
3747: +  x - the vector
3748: .  m - first dimension of four dimensional array
3749: .  n - second dimension of four dimensional array
3750: .  p - third dimension of four dimensional array
3751: .  q - fourth dimension of four dimensional array
3752: .  mstart - first index you will use in first coordinate direction (often 0)
3753: .  nstart - first index in the second coordinate direction (often 0)
3754: .  pstart - first index in the third coordinate direction (often 0)
3755: -  qstart - first index in the fourth coordinate direction (often 0)

3757:    Output Parameter:
3758: .  a - location to put pointer to the array

3760:    Level: beginner

3762:   Notes:
3763:    For a vector obtained from `DMCreateLocalVector()` mstart, nstart, and pstart are likely
3764:    obtained from the corner indices obtained from `DMDAGetGhostCorners()` while for
3765:    `DMCreateGlobalVector()` they are the corner indices from `DMDAGetCorners()`. In both cases
3766:    the arguments from `DMDAGet[Ghost]Corners()` are reversed in the call to `VecGetArray3d()`.

3768:    For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.

3770: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
3771:           `VecRestoreArray2d()`, `DMDAVecGetarray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
3772:           `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3773: @*/
3774: PetscErrorCode VecGetArray4dWrite(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt q, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscInt qstart, PetscScalar ****a[])
3775: {
3776:   PetscInt     i, N, j, k;
3777:   PetscScalar *aa, ***b, **c;

3782:   VecGetLocalSize(x, &N);
3784:   VecGetArrayWrite(x, &aa);

3786:   PetscMalloc(m * sizeof(PetscScalar ***) + m * n * sizeof(PetscScalar **) + m * n * p * sizeof(PetscScalar *), a);
3787:   b = (PetscScalar ***)((*a) + m);
3788:   c = (PetscScalar **)(b + m * n);
3789:   for (i = 0; i < m; i++) (*a)[i] = b + i * n - nstart;
3790:   for (i = 0; i < m; i++)
3791:     for (j = 0; j < n; j++) b[i * n + j] = c + i * n * p + j * p - pstart;
3792:   for (i = 0; i < m; i++)
3793:     for (j = 0; j < n; j++)
3794:       for (k = 0; k < p; k++) c[i * n * p + j * p + k] = aa + i * n * p * q + j * p * q + k * q - qstart;
3795:   *a -= mstart;
3796:   return 0;
3797: }

3799: /*@C
3800:    VecRestoreArray4d - Restores a vector after `VecGetArray4d()` has been called.

3802:    Logically Collective

3804:    Input Parameters:
3805: +  x - the vector
3806: .  m - first dimension of four dimensional array
3807: .  n - second dimension of the four dimensional array
3808: .  p - third dimension of the four dimensional array
3809: .  q - fourth dimension of the four dimensional array
3810: .  mstart - first index you will use in first coordinate direction (often 0)
3811: .  nstart - first index in the second coordinate direction (often 0)
3812: .  pstart - first index in the third coordinate direction (often 0)
3813: .  qstart - first index in the fourth coordinate direction (often 0)
3814: -  a - location of pointer to array obtained from VecGetArray4d()

3816:    Level: beginner

3818:    Notes:
3819:    For regular PETSc vectors this routine does not involve any copies. For
3820:    any special vectors that do not store local vector data in a contiguous
3821:    array, this routine will copy the data back into the underlying
3822:    vector data structure from the array obtained with `VecGetArray()`.

3824:    This routine actually zeros out the a pointer.

3826: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
3827:           `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
3828:           `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`, `VecGet`
3829: @*/
3830: PetscErrorCode VecRestoreArray4d(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt q, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscInt qstart, PetscScalar ****a[])
3831: {
3832:   void *dummy;

3837:   dummy = (void *)(*a + mstart);
3838:   PetscFree(dummy);
3839:   VecRestoreArray(x, NULL);
3840:   return 0;
3841: }

3843: /*@C
3844:    VecRestoreArray4dWrite - Restores a vector after `VecGetArray4dWrite()` has been called.

3846:    Logically Collective

3848:    Input Parameters:
3849: +  x - the vector
3850: .  m - first dimension of four dimensional array
3851: .  n - second dimension of the four dimensional array
3852: .  p - third dimension of the four dimensional array
3853: .  q - fourth dimension of the four dimensional array
3854: .  mstart - first index you will use in first coordinate direction (often 0)
3855: .  nstart - first index in the second coordinate direction (often 0)
3856: .  pstart - first index in the third coordinate direction (often 0)
3857: .  qstart - first index in the fourth coordinate direction (often 0)
3858: -  a - location of pointer to array obtained from `VecGetArray4d()`

3860:    Level: beginner

3862:    Notes:
3863:    For regular PETSc vectors this routine does not involve any copies. For
3864:    any special vectors that do not store local vector data in a contiguous
3865:    array, this routine will copy the data back into the underlying
3866:    vector data structure from the array obtained with `VecGetArray()`.

3868:    This routine actually zeros out the a pointer.

3870: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
3871:           `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
3872:           `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`, `VecGet`
3873: @*/
3874: PetscErrorCode VecRestoreArray4dWrite(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt q, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscInt qstart, PetscScalar ****a[])
3875: {
3876:   void *dummy;

3881:   dummy = (void *)(*a + mstart);
3882:   PetscFree(dummy);
3883:   VecRestoreArrayWrite(x, NULL);
3884:   return 0;
3885: }

3887: /*@C
3888:    VecGetArray2dRead - Returns a pointer to a 2d contiguous array that contains this
3889:    processor's portion of the vector data.  You MUST call `VecRestoreArray2dRead()`
3890:    when you no longer need access to the array.

3892:    Logically Collective

3894:    Input Parameters:
3895: +  x - the vector
3896: .  m - first dimension of two dimensional array
3897: .  n - second dimension of two dimensional array
3898: .  mstart - first index you will use in first coordinate direction (often 0)
3899: -  nstart - first index in the second coordinate direction (often 0)

3901:    Output Parameter:
3902: .  a - location to put pointer to the array

3904:    Level: developer

3906:   Notes:
3907:    For a vector obtained from `DMCreateLocalVector()` mstart and nstart are likely
3908:    obtained from the corner indices obtained from `DMDAGetGhostCorners()` while for
3909:    `DMCreateGlobalVector()` they are the corner indices from `DMDAGetCorners()`. In both cases
3910:    the arguments from `DMDAGet[Ghost]Corners()` are reversed in the call to `VecGetArray2d()`.

3912:    For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.

3914: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
3915:           `VecRestoreArray2d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
3916:           `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3917: @*/
3918: PetscErrorCode VecGetArray2dRead(Vec x, PetscInt m, PetscInt n, PetscInt mstart, PetscInt nstart, PetscScalar **a[])
3919: {
3920:   PetscInt           i, N;
3921:   const PetscScalar *aa;

3926:   VecGetLocalSize(x, &N);
3928:   VecGetArrayRead(x, &aa);

3930:   PetscMalloc1(m, a);
3931:   for (i = 0; i < m; i++) (*a)[i] = (PetscScalar *)aa + i * n - nstart;
3932:   *a -= mstart;
3933:   return 0;
3934: }

3936: /*@C
3937:    VecRestoreArray2dRead - Restores a vector after `VecGetArray2dRead()` has been called.

3939:    Logically Collective

3941:    Input Parameters:
3942: +  x - the vector
3943: .  m - first dimension of two dimensional array
3944: .  n - second dimension of the two dimensional array
3945: .  mstart - first index you will use in first coordinate direction (often 0)
3946: .  nstart - first index in the second coordinate direction (often 0)
3947: -  a - location of pointer to array obtained from VecGetArray2d()

3949:    Level: developer

3951:    Notes:
3952:    For regular PETSc vectors this routine does not involve any copies. For
3953:    any special vectors that do not store local vector data in a contiguous
3954:    array, this routine will copy the data back into the underlying
3955:    vector data structure from the array obtained with `VecGetArray()`.

3957:    This routine actually zeros out the a pointer.

3959: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
3960:           `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
3961:           `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
3962: @*/
3963: PetscErrorCode VecRestoreArray2dRead(Vec x, PetscInt m, PetscInt n, PetscInt mstart, PetscInt nstart, PetscScalar **a[])
3964: {
3965:   void *dummy;

3970:   dummy = (void *)(*a + mstart);
3971:   PetscFree(dummy);
3972:   VecRestoreArrayRead(x, NULL);
3973:   return 0;
3974: }

3976: /*@C
3977:    VecGetArray1dRead - Returns a pointer to a 1d contiguous array that contains this
3978:    processor's portion of the vector data.  You MUST call `VecRestoreArray1dRead()`
3979:    when you no longer need access to the array.

3981:    Logically Collective

3983:    Input Parameters:
3984: +  x - the vector
3985: .  m - first dimension of two dimensional array
3986: -  mstart - first index you will use in first coordinate direction (often 0)

3988:    Output Parameter:
3989: .  a - location to put pointer to the array

3991:    Level: developer

3993:   Notes:
3994:    For a vector obtained from `DMCreateLocalVector()` mstart are likely
3995:    obtained from the corner indices obtained from `DMDAGetGhostCorners()` while for
3996:    `DMCreateGlobalVector()` they are the corner indices from `DMDAGetCorners()`.

3998:    For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.

4000: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
4001:           `VecRestoreArray2d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
4002:           `VecGetArray2d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
4003: @*/
4004: PetscErrorCode VecGetArray1dRead(Vec x, PetscInt m, PetscInt mstart, PetscScalar *a[])
4005: {
4006:   PetscInt N;

4011:   VecGetLocalSize(x, &N);
4013:   VecGetArrayRead(x, (const PetscScalar **)a);
4014:   *a -= mstart;
4015:   return 0;
4016: }

4018: /*@C
4019:    VecRestoreArray1dRead - Restores a vector after `VecGetArray1dRead()` has been called.

4021:    Logically Collective

4023:    Input Parameters:
4024: +  x - the vector
4025: .  m - first dimension of two dimensional array
4026: .  mstart - first index you will use in first coordinate direction (often 0)
4027: -  a - location of pointer to array obtained from `VecGetArray1dRead()`

4029:    Level: developer

4031:    Notes:
4032:    For regular PETSc vectors this routine does not involve any copies. For
4033:    any special vectors that do not store local vector data in a contiguous
4034:    array, this routine will copy the data back into the underlying
4035:    vector data structure from the array obtained with `VecGetArray1dRead()`.

4037:    This routine actually zeros out the a pointer.

4039: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
4040:           `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
4041:           `VecGetArray1d()`, `VecRestoreArray2d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
4042: @*/
4043: PetscErrorCode VecRestoreArray1dRead(Vec x, PetscInt m, PetscInt mstart, PetscScalar *a[])
4044: {
4047:   VecRestoreArrayRead(x, NULL);
4048:   return 0;
4049: }

4051: /*@C
4052:    VecGetArray3dRead - Returns a pointer to a 3d contiguous array that contains this
4053:    processor's portion of the vector data.  You MUST call `VecRestoreArray3dRead()`
4054:    when you no longer need access to the array.

4056:    Logically Collective

4058:    Input Parameters:
4059: +  x - the vector
4060: .  m - first dimension of three dimensional array
4061: .  n - second dimension of three dimensional array
4062: .  p - third dimension of three dimensional array
4063: .  mstart - first index you will use in first coordinate direction (often 0)
4064: .  nstart - first index in the second coordinate direction (often 0)
4065: -  pstart - first index in the third coordinate direction (often 0)

4067:    Output Parameter:
4068: .  a - location to put pointer to the array

4070:    Level: developer

4072:   Notes:
4073:    For a vector obtained from `DMCreateLocalVector()` mstart, nstart, and pstart are likely
4074:    obtained from the corner indices obtained from `DMDAGetGhostCorners()` while for
4075:    `DMCreateGlobalVector()` they are the corner indices from `DMDAGetCorners()`. In both cases
4076:    the arguments from `DMDAGet[Ghost]Corners()` are reversed in the call to `VecGetArray3dRead()`.

4078:    For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.

4080: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
4081:           `VecRestoreArray2d()`, `DMDAVecGetarray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
4082:           `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
4083: @*/
4084: PetscErrorCode VecGetArray3dRead(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscScalar ***a[])
4085: {
4086:   PetscInt           i, N, j;
4087:   const PetscScalar *aa;
4088:   PetscScalar      **b;

4093:   VecGetLocalSize(x, &N);
4095:   VecGetArrayRead(x, &aa);

4097:   PetscMalloc(m * sizeof(PetscScalar **) + m * n * sizeof(PetscScalar *), a);
4098:   b = (PetscScalar **)((*a) + m);
4099:   for (i = 0; i < m; i++) (*a)[i] = b + i * n - nstart;
4100:   for (i = 0; i < m; i++)
4101:     for (j = 0; j < n; j++) b[i * n + j] = (PetscScalar *)aa + i * n * p + j * p - pstart;
4102:   *a -= mstart;
4103:   return 0;
4104: }

4106: /*@C
4107:    VecRestoreArray3dRead - Restores a vector after `VecGetArray3dRead()` has been called.

4109:    Logically Collective

4111:    Input Parameters:
4112: +  x - the vector
4113: .  m - first dimension of three dimensional array
4114: .  n - second dimension of the three dimensional array
4115: .  p - third dimension of the three dimensional array
4116: .  mstart - first index you will use in first coordinate direction (often 0)
4117: .  nstart - first index in the second coordinate direction (often 0)
4118: .  pstart - first index in the third coordinate direction (often 0)
4119: -  a - location of pointer to array obtained from `VecGetArray3dRead()`

4121:    Level: developer

4123:    Notes:
4124:    For regular PETSc vectors this routine does not involve any copies. For
4125:    any special vectors that do not store local vector data in a contiguous
4126:    array, this routine will copy the data back into the underlying
4127:    vector data structure from the array obtained with `VecGetArray()`.

4129:    This routine actually zeros out the a pointer.

4131: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
4132:           `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
4133:           `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`, `VecGet`
4134: @*/
4135: PetscErrorCode VecRestoreArray3dRead(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscScalar ***a[])
4136: {
4137:   void *dummy;

4142:   dummy = (void *)(*a + mstart);
4143:   PetscFree(dummy);
4144:   VecRestoreArrayRead(x, NULL);
4145:   return 0;
4146: }

4148: /*@C
4149:    VecGetArray4dRead - Returns a pointer to a 4d contiguous array that contains this
4150:    processor's portion of the vector data.  You MUST call `VecRestoreArray4dRead()`
4151:    when you no longer need access to the array.

4153:    Logically Collective

4155:    Input Parameters:
4156: +  x - the vector
4157: .  m - first dimension of four dimensional array
4158: .  n - second dimension of four dimensional array
4159: .  p - third dimension of four dimensional array
4160: .  q - fourth dimension of four dimensional array
4161: .  mstart - first index you will use in first coordinate direction (often 0)
4162: .  nstart - first index in the second coordinate direction (often 0)
4163: .  pstart - first index in the third coordinate direction (often 0)
4164: -  qstart - first index in the fourth coordinate direction (often 0)

4166:    Output Parameter:
4167: .  a - location to put pointer to the array

4169:    Level: beginner

4171:   Notes:
4172:    For a vector obtained from `DMCreateLocalVector()` mstart, nstart, and pstart are likely
4173:    obtained from the corner indices obtained from `DMDAGetGhostCorners()` while for
4174:    `DMCreateGlobalVector()` they are the corner indices from `DMDAGetCorners()`. In both cases
4175:    the arguments from `DMDAGet[Ghost]Corners()` are reversed in the call to `VecGetArray3d()`.

4177:    For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.

4179: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecGetArrays()`, `VecGetArrayF90()`, `VecPlaceArray()`,
4180:           `VecRestoreArray2d()`, `DMDAVecGetarray()`, `DMDAVecRestoreArray()`, `VecGetArray3d()`, `VecRestoreArray3d()`,
4181:           `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`
4182: @*/
4183: PetscErrorCode VecGetArray4dRead(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt q, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscInt qstart, PetscScalar ****a[])
4184: {
4185:   PetscInt           i, N, j, k;
4186:   const PetscScalar *aa;
4187:   PetscScalar     ***b, **c;

4192:   VecGetLocalSize(x, &N);
4194:   VecGetArrayRead(x, &aa);

4196:   PetscMalloc(m * sizeof(PetscScalar ***) + m * n * sizeof(PetscScalar **) + m * n * p * sizeof(PetscScalar *), a);
4197:   b = (PetscScalar ***)((*a) + m);
4198:   c = (PetscScalar **)(b + m * n);
4199:   for (i = 0; i < m; i++) (*a)[i] = b + i * n - nstart;
4200:   for (i = 0; i < m; i++)
4201:     for (j = 0; j < n; j++) b[i * n + j] = c + i * n * p + j * p - pstart;
4202:   for (i = 0; i < m; i++)
4203:     for (j = 0; j < n; j++)
4204:       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;
4205:   *a -= mstart;
4206:   return 0;
4207: }

4209: /*@C
4210:    VecRestoreArray4dRead - Restores a vector after `VecGetArray4d()` has been called.

4212:    Logically Collective

4214:    Input Parameters:
4215: +  x - the vector
4216: .  m - first dimension of four dimensional array
4217: .  n - second dimension of the four dimensional array
4218: .  p - third dimension of the four dimensional array
4219: .  q - fourth dimension of the four dimensional array
4220: .  mstart - first index you will use in first coordinate direction (often 0)
4221: .  nstart - first index in the second coordinate direction (often 0)
4222: .  pstart - first index in the third coordinate direction (often 0)
4223: .  qstart - first index in the fourth coordinate direction (often 0)
4224: -  a - location of pointer to array obtained from `VecGetArray4dRead()`

4226:    Level: beginner

4228:    Notes:
4229:    For regular PETSc vectors this routine does not involve any copies. For
4230:    any special vectors that do not store local vector data in a contiguous
4231:    array, this routine will copy the data back into the underlying
4232:    vector data structure from the array obtained with `VecGetArray()`.

4234:    This routine actually zeros out the a pointer.

4236: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecRestoreArrays()`, `VecRestoreArrayF90()`, `VecPlaceArray()`,
4237:           `VecGetArray2d()`, `VecGetArray3d()`, `VecRestoreArray3d()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`
4238:           `VecGetArray1d()`, `VecRestoreArray1d()`, `VecGetArray4d()`, `VecRestoreArray4d()`, `VecGet`
4239: @*/
4240: PetscErrorCode VecRestoreArray4dRead(Vec x, PetscInt m, PetscInt n, PetscInt p, PetscInt q, PetscInt mstart, PetscInt nstart, PetscInt pstart, PetscInt qstart, PetscScalar ****a[])
4241: {
4242:   void *dummy;

4247:   dummy = (void *)(*a + mstart);
4248:   PetscFree(dummy);
4249:   VecRestoreArrayRead(x, NULL);
4250:   return 0;
4251: }

4253: #if defined(PETSC_USE_DEBUG)

4255: /*@
4256:    VecLockGet  - Gets the current lock status of a vector

4258:    Logically Collective

4260:    Input Parameter:
4261: .  x - the vector

4263:    Output Parameter:
4264: .  state - greater than zero indicates the vector is locked for read; less then zero indicates the vector is
4265:            locked for write; equal to zero means the vector is unlocked, that is, it is free to read or write.

4267:    Level: beginner

4269: .seealso: [](chapter_vectors), `Vec`, `VecRestoreArray()`, `VecGetArrayRead()`, `VecLockReadPush()`, `VecLockReadPop()`
4270: @*/
4271: PetscErrorCode VecLockGet(Vec x, PetscInt *state)
4272: {
4274:   *state = x->lock;
4275:   return 0;
4276: }

4278: PetscErrorCode VecLockGetLocation(Vec x, const char *file[], const char *func[], int *line)
4279: {
4284:   {
4285:     const int index = x->lockstack.currentsize - 1;

4288:     *file = x->lockstack.file[index];
4289:     *func = x->lockstack.function[index];
4290:     *line = x->lockstack.line[index];
4291:   }
4292:   return 0;
4293: }

4295: /*@
4296:    VecLockReadPush  - Pushes a read-only lock on a vector to prevent it from writing

4298:    Logically Collective

4300:    Input Parameter:
4301: .  x - the vector

4303:    Level: beginner

4305:    Notes:
4306:     If this is set then calls to `VecGetArray()` or `VecSetValues()` or any other routines that change the vectors values will fail.

4308:     The call can be nested, i.e., called multiple times on the same vector, but each `VecLockReadPush()` has to have one matching
4309:     `VecLockReadPop()`, which removes the latest read-only lock.

4311: .seealso: [](chapter_vectors), `Vec`, `VecRestoreArray()`, `VecGetArrayRead()`, `VecLockReadPop()`, `VecLockGet()`
4312: @*/
4313: PetscErrorCode VecLockReadPush(Vec x)
4314: {
4315:   const char *file, *func;
4316:   int         index, line;

4320:   if ((index = petscstack.currentsize - 2) == -1) {
4321:     // vec was locked "outside" of petsc, either in user-land or main. the error message will
4322:     // now show this function as the culprit, but it will include the stacktrace
4323:     file = "unknown user-file";
4324:     func = "unknown_user_function";
4325:     line = 0;
4326:   } else {
4328:     file = petscstack.file[index];
4329:     func = petscstack.function[index];
4330:     line = petscstack.line[index];
4331:   }
4332:   PetscStackPush_Private(x->lockstack, file, func, line, petscstack.petscroutine[index], PETSC_FALSE);
4333:   return 0;
4334: }

4336: /*@
4337:    VecLockReadPop  - Pops a read-only lock from a vector

4339:    Logically Collective

4341:    Input Parameter:
4342: .  x - the vector

4344:    Level: beginner

4346: .seealso: [](chapter_vectors), `Vec`, `VecRestoreArray()`, `VecGetArrayRead()`, `VecLockReadPush()`, `VecLockGet()`
4347: @*/
4348: PetscErrorCode VecLockReadPop(Vec x)
4349: {
4352:   {
4353:     const char *previous = x->lockstack.function[x->lockstack.currentsize - 1];

4355:     PetscStackPop_Private(x->lockstack, previous);
4356:   }
4357:   return 0;
4358: }

4360: /*@C
4361:    VecLockWriteSet  - Lock or unlock a vector for exclusive read/write access

4363:    Logically Collective

4365:    Input Parameters:
4366: +  x   - the vector
4367: -  flg - PETSC_TRUE to lock the vector for exclusive read/write access; PETSC_FALSE to unlock it.

4369:    Level: beginner

4371:    Notes:
4372:     The function is useful in split-phase computations, which usually have a begin phase and an end phase.
4373:     One can call `VecLockWriteSet`(x,`PETSC_TRUE`) in the begin phase to lock a vector for exclusive
4374:     access, and call `VecLockWriteSet`(x,`PETSC_FALSE`) in the end phase to unlock the vector from exclusive
4375:     access. In this way, one is ensured no other operations can access the vector in between. The code may like

4377: .vb
4378:        VecGetArray(x,&xdata); // begin phase
4379:        VecLockWriteSet(v,PETSC_TRUE);

4381:        Other operations, which can not access x anymore (they can access xdata, of course)

4383:        VecRestoreArray(x,&vdata); // end phase
4384:        VecLockWriteSet(v,PETSC_FALSE);
4385: .ve

4387:     The call can not be nested on the same vector, in other words, one can not call `VecLockWriteSet`(x,`PETSC_TRUE`)
4388:     again before calling `VecLockWriteSet`(v,`PETSC_FALSE`).

4390: .seealso: [](chapter_vectors), `Vec`, `VecRestoreArray()`, `VecGetArrayRead()`, `VecLockReadPush()`, `VecLockReadPop()`, `VecLockGet()`
4391: @*/
4392: PetscErrorCode VecLockWriteSet(Vec x, PetscBool flg)
4393: {
4395:   if (flg) {
4398:     x->lock = -1;
4399:   } else {
4401:     x->lock = 0;
4402:   }
4403:   return 0;
4404: }

4406: /*@
4407:    VecLockPush  - Pushes a read-only lock on a vector to prevent it from writing

4409:    Level: deprecated

4411: .seealso: [](chapter_vectors), `Vec`, `VecLockReadPush()`
4412: @*/
4413: PetscErrorCode VecLockPush(Vec x)
4414: {
4415:   VecLockReadPush(x);
4416:   return 0;
4417: }

4419: /*@
4420:    VecLockPop  - Pops a read-only lock from a vector

4422:    Level: deprecated

4424: .seealso: [](chapter_vectors), `Vec`, `VecLockReadPop()`
4425: @*/
4426: PetscErrorCode VecLockPop(Vec x)
4427: {
4428:   VecLockReadPop(x);
4429:   return 0;
4430: }

4432: #endif