Actual source code: vector.c
1: /*
2: Provides the interface functions for vector operations that do NOT have PetscScalar/PetscReal in the signature
3: These are the vector functions the user calls.
4: */
5: #include <petsc/private/vecimpl.h>
6: #include <petsc/private/deviceimpl.h>
8: /* Logging support */
9: PetscClassId VEC_CLASSID;
10: PetscLogEvent VEC_View, VEC_Max, VEC_Min, VEC_Dot, VEC_MDot, VEC_TDot;
11: PetscLogEvent VEC_Norm, VEC_Normalize, VEC_Scale, VEC_Copy, VEC_Set, VEC_AXPY, VEC_AYPX, VEC_WAXPY;
12: PetscLogEvent VEC_MTDot, VEC_MAXPY, VEC_Swap, VEC_AssemblyBegin, VEC_ScatterBegin, VEC_ScatterEnd;
13: PetscLogEvent VEC_AssemblyEnd, VEC_PointwiseMult, VEC_SetValues, VEC_Load, VEC_SetPreallocateCOO, VEC_SetValuesCOO;
14: PetscLogEvent VEC_SetRandom, VEC_ReduceArithmetic, VEC_ReduceCommunication, VEC_ReduceBegin, VEC_ReduceEnd, VEC_Ops;
15: PetscLogEvent VEC_DotNorm2, VEC_AXPBYPCZ;
16: PetscLogEvent VEC_ViennaCLCopyFromGPU, VEC_ViennaCLCopyToGPU;
17: PetscLogEvent VEC_CUDACopyFromGPU, VEC_CUDACopyToGPU;
18: PetscLogEvent VEC_HIPCopyFromGPU, VEC_HIPCopyToGPU;
20: /*@
21: VecStashGetInfo - Gets how many values are currently in the vector stash, i.e. need
22: to be communicated to other processors during the VecAssemblyBegin/End() process
24: Not collective
26: Input Parameter:
27: . vec - the vector
29: Output Parameters:
30: + nstash - the size of the stash
31: . reallocs - the number of additional mallocs incurred.
32: . bnstash - the size of the block stash
33: - breallocs - the number of additional mallocs incurred.in the block stash
35: Level: advanced
37: .seealso: [](chapter_vectors), `Vec`, `VecAssemblyBegin()`, `VecAssemblyEnd()`, `Vec`, `VecStashSetInitialSize()`, `VecStashView()`
38: @*/
39: PetscErrorCode VecStashGetInfo(Vec vec, PetscInt *nstash, PetscInt *reallocs, PetscInt *bnstash, PetscInt *breallocs)
40: {
41: PetscFunctionBegin;
42: PetscCall(VecStashGetInfo_Private(&vec->stash, nstash, reallocs));
43: PetscCall(VecStashGetInfo_Private(&vec->bstash, bnstash, breallocs));
44: PetscFunctionReturn(PETSC_SUCCESS);
45: }
47: /*@
48: VecSetLocalToGlobalMapping - Sets a local numbering to global numbering used
49: by the routine `VecSetValuesLocal()` to allow users to insert vector entries
50: using a local (per-processor) numbering.
52: Logically Collective
54: Input Parameters:
55: + x - vector
56: - mapping - mapping created with ISLocalToGlobalMappingCreate() or ISLocalToGlobalMappingCreateIS()
58: Level: intermediate
60: Note:
61: All vectors obtained with `VecDuplicate()` from this vector inherit the same mapping.
63: seealso: [](chapter_vectors), `Vec`, `VecAssemblyBegin()`, `VecAssemblyEnd()`, `VecSetValues()`, `VecSetValuesLocal()`,
64: `VecSetLocalToGlobalMapping()`, `VecSetValuesBlockedLocal()`
65: @*/
66: PetscErrorCode VecSetLocalToGlobalMapping(Vec x, ISLocalToGlobalMapping mapping)
67: {
68: PetscFunctionBegin;
71: if (x->ops->setlocaltoglobalmapping) PetscUseTypeMethod(x, setlocaltoglobalmapping, mapping);
72: else PetscCall(PetscLayoutSetISLocalToGlobalMapping(x->map, mapping));
73: PetscFunctionReturn(PETSC_SUCCESS);
74: }
76: /*@
77: VecGetLocalToGlobalMapping - Gets the local-to-global numbering set by `VecSetLocalToGlobalMapping()`
79: Not Collective
81: Input Parameter:
82: . X - the vector
84: Output Parameter:
85: . mapping - the mapping
87: Level: advanced
89: .seealso: [](chapter_vectors), `Vec`, `VecSetValuesLocal()`
90: @*/
91: PetscErrorCode VecGetLocalToGlobalMapping(Vec X, ISLocalToGlobalMapping *mapping)
92: {
93: PetscFunctionBegin;
97: *mapping = X->map->mapping;
98: PetscFunctionReturn(PETSC_SUCCESS);
99: }
101: /*@
102: VecAssemblyBegin - Begins assembling the vector. This routine should
103: be called after completing all calls to `VecSetValues()`.
105: Collective
107: Input Parameter:
108: . vec - the vector
110: Level: beginner
112: .seealso: [](chapter_vectors), `Vec`, `VecAssemblyEnd()`, `VecSetValues()`
113: @*/
114: PetscErrorCode VecAssemblyBegin(Vec vec)
115: {
116: PetscFunctionBegin;
119: PetscCall(VecStashViewFromOptions(vec, NULL, "-vec_view_stash"));
120: PetscCall(PetscLogEventBegin(VEC_AssemblyBegin, vec, 0, 0, 0));
121: PetscTryTypeMethod(vec, assemblybegin);
122: PetscCall(PetscLogEventEnd(VEC_AssemblyBegin, vec, 0, 0, 0));
123: PetscCall(PetscObjectStateIncrease((PetscObject)vec));
124: PetscFunctionReturn(PETSC_SUCCESS);
125: }
127: /*@
128: VecAssemblyEnd - Completes assembling the vector. This routine should
129: be called after `VecAssemblyBegin()`.
131: Collective
133: Input Parameter:
134: . vec - the vector
136: Options Database Keys:
137: + -vec_view - Prints vector in `PETSC_VIEWER_DEFAULT` format
138: . -vec_view ::ascii_matlab - Prints vector in `PETSC_VIEWER_ASCII_MATLAB` format to stdout
139: . -vec_view matlab:filename - Prints vector in MATLAB .mat file to filename (requires PETSc configured with --with-matlab)
140: . -vec_view draw - Activates vector viewing using drawing tools
141: . -display <name> - Sets display name (default is host)
142: . -draw_pause <sec> - Sets number of seconds to pause after display
143: - -vec_view socket - Activates vector viewing using a socket
145: Level: beginner
147: .seealso: [](chapter_vectors), `Vec`, `VecAssemblyBegin()`, `VecSetValues()`
148: @*/
149: PetscErrorCode VecAssemblyEnd(Vec vec)
150: {
151: PetscFunctionBegin;
153: PetscCall(PetscLogEventBegin(VEC_AssemblyEnd, vec, 0, 0, 0));
155: PetscTryTypeMethod(vec, assemblyend);
156: PetscCall(PetscLogEventEnd(VEC_AssemblyEnd, vec, 0, 0, 0));
157: PetscCall(VecViewFromOptions(vec, NULL, "-vec_view"));
158: PetscFunctionReturn(PETSC_SUCCESS);
159: }
161: /*@
162: VecSetPreallocationCOO - set preallocation for a vector using a coordinate format of the entries with global indices
164: Collective
166: Input Parameters:
167: + x - vector being preallocated
168: . ncoo - number of entries
169: - coo_i - entry indices
171: Level: beginner
173: Notes:
174: Entries can be repeated, see `VecSetValuesCOO()`. Negative indices are not allowed unless vector option `VEC_IGNORE_NEGATIVE_INDICES` is set,
175: in which case they, along with the corresponding entries in `VecSetValuesCOO()`, are ignored. If vector option `VEC_NO_OFF_PROC_ENTRIES` is set,
176: remote entries are ignored, otherwise, they will be properly added or inserted to the vector.
178: The array coo_i[] may be freed immediately after calling this function.
180: .seealso: [](chapter_vectors), `Vec`, VecSetValuesCOO(), VecSetPreallocationCOOLocal()
181: @*/
182: PetscErrorCode VecSetPreallocationCOO(Vec x, PetscCount ncoo, const PetscInt coo_i[])
183: {
184: PetscFunctionBegin;
188: PetscCall(PetscLogEventBegin(VEC_SetPreallocateCOO, x, 0, 0, 0));
189: PetscCall(PetscLayoutSetUp(x->map));
190: if (x->ops->setpreallocationcoo) {
191: PetscUseTypeMethod(x, setpreallocationcoo, ncoo, coo_i);
192: } else {
193: IS is_coo_i;
194: /* The default implementation only supports ncoo within limit of PetscInt */
195: PetscCheck(ncoo <= PETSC_MAX_INT, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "ncoo %" PetscCount_FMT " overflowed PetscInt; configure --with-64-bit-indices or request support", ncoo);
196: PetscCall(ISCreateGeneral(PETSC_COMM_SELF, ncoo, coo_i, PETSC_COPY_VALUES, &is_coo_i));
197: PetscCall(PetscObjectCompose((PetscObject)x, "__PETSc_coo_i", (PetscObject)is_coo_i));
198: PetscCall(ISDestroy(&is_coo_i));
199: }
200: PetscCall(PetscLogEventEnd(VEC_SetPreallocateCOO, x, 0, 0, 0));
201: PetscFunctionReturn(PETSC_SUCCESS);
202: }
204: /*@
205: VecSetPreallocationCOOLocal - set preallocation for vectors using a coordinate format of the entries with local indices
207: Collective
209: Input Parameters:
210: + x - vector being preallocated
211: . ncoo - number of entries
212: - coo_i - row indices (local numbering; may be modified)
214: Level: beginner
216: Notes:
217: The local indices are translated using the local to global mapping, thus `VecSetLocalToGlobalMapping()` must have been
218: called prior to this function.
220: The indices coo_i may be modified within this function. They might be translated to corresponding global
221: indices, but the caller should not rely on them having any specific value after this function returns. The arrays
222: can be freed or reused immediately after this function returns.
224: Entries can be repeated. Negative indices and remote indices might be allowed. see `VecSetPreallocationCOO()`.
226: .seealso: [](chapter_vectors), `Vec`, VecSetPreallocationCOO(), VecSetValuesCOO()
227: @*/
228: PetscErrorCode VecSetPreallocationCOOLocal(Vec x, PetscCount ncoo, PetscInt coo_i[])
229: {
230: ISLocalToGlobalMapping ltog;
232: PetscFunctionBegin;
236: PetscCheck(ncoo <= PETSC_MAX_INT, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "ncoo %" PetscCount_FMT " overflowed PetscInt; configure --with-64-bit-indices or request support", ncoo);
237: PetscCall(PetscLayoutSetUp(x->map));
238: PetscCall(VecGetLocalToGlobalMapping(x, <og));
239: if (ltog) PetscCall(ISLocalToGlobalMappingApply(ltog, ncoo, coo_i, coo_i));
240: PetscCall(VecSetPreallocationCOO(x, ncoo, coo_i));
241: PetscFunctionReturn(PETSC_SUCCESS);
242: }
244: /*@
245: VecSetValuesCOO - set values at once in a vector preallocated using `VecSetPreallocationCOO()`
247: Collective
249: Input Parameters:
250: + x - vector being set
251: . coo_v - the value array
252: - imode - the insert mode
254: Level: beginner
256: Note:
257: The values must follow the order of the indices prescribed with `VecSetPreallocationCOO()` or `VecSetPreallocationCOOLocal()`.
258: When repeated entries are specified in the COO indices the coo_v values are first properly summed, regardless of the value of imode.
259: The imode flag indicates if coo_v must be added to the current values of the vector (`ADD_VALUES`) or overwritten (`INSERT_VALUES`).
260: `VecAssemblyBegin()` and `VecAssemblyEnd()` do not need to be called after this routine. It automatically handles the assembly process.
262: .seealso: [](chapter_vectors), `Vec`, VecSetPreallocationCOO(), VecSetPreallocationCOOLocal(), VecSetValues()
263: @*/
264: PetscErrorCode VecSetValuesCOO(Vec x, const PetscScalar coo_v[], InsertMode imode)
265: {
266: PetscFunctionBegin;
270: PetscCall(PetscLogEventBegin(VEC_SetValuesCOO, x, 0, 0, 0));
271: if (x->ops->setvaluescoo) {
272: PetscUseTypeMethod(x, setvaluescoo, coo_v, imode);
273: PetscCall(PetscObjectStateIncrease((PetscObject)x));
274: } else {
275: IS is_coo_i;
276: const PetscInt *coo_i;
277: PetscInt ncoo;
278: PetscMemType mtype;
280: PetscCall(PetscGetMemType(coo_v, &mtype));
281: PetscCheck(mtype == PETSC_MEMTYPE_HOST, PetscObjectComm((PetscObject)x), PETSC_ERR_ARG_WRONG, "The basic VecSetValuesCOO() only supports v[] on host");
282: PetscCall(PetscObjectQuery((PetscObject)x, "__PETSc_coo_i", (PetscObject *)&is_coo_i));
283: PetscCheck(is_coo_i, PetscObjectComm((PetscObject)x), PETSC_ERR_COR, "Missing coo_i IS");
284: PetscCall(ISGetLocalSize(is_coo_i, &ncoo));
285: PetscCall(ISGetIndices(is_coo_i, &coo_i));
286: if (imode != ADD_VALUES) PetscCall(VecZeroEntries(x));
287: PetscCall(VecSetValues(x, ncoo, coo_i, coo_v, ADD_VALUES));
288: PetscCall(ISRestoreIndices(is_coo_i, &coo_i));
289: PetscCall(VecAssemblyBegin(x));
290: PetscCall(VecAssemblyEnd(x));
291: }
292: PetscCall(PetscLogEventEnd(VEC_SetValuesCOO, x, 0, 0, 0));
293: PetscFunctionReturn(PETSC_SUCCESS);
294: }
296: static PetscErrorCode VecPointwiseApply_Private(Vec w, Vec x, Vec y, PetscLogEvent event, PetscErrorCode (*const pointwise_op)(Vec, Vec, Vec))
297: {
298: PetscFunctionBegin;
305: PetscCheckSameTypeAndComm(x, 2, y, 3);
306: PetscCheckSameTypeAndComm(y, 3, w, 1);
307: VecCheckSameSize(w, 1, x, 2);
308: VecCheckSameSize(w, 1, y, 3);
309: VecCheckAssembled(x);
310: VecCheckAssembled(y);
311: PetscCall(VecSetErrorIfLocked(w, 1));
314: if (event) PetscCall(PetscLogEventBegin(event, x, y, w, 0));
315: PetscCall((*pointwise_op)(w, x, y));
316: if (event) PetscCall(PetscLogEventEnd(event, x, y, w, 0));
317: PetscCall(PetscObjectStateIncrease((PetscObject)w));
318: PetscFunctionReturn(PETSC_SUCCESS);
319: }
321: /*@
322: VecPointwiseMax - Computes the component-wise maximum `w[i] = max(x[i], y[i])`.
324: Logically Collective
326: Input Parameters:
327: + x - the first input vector
328: - y - the second input vector
330: Output Parameter:
331: . w - the result
333: Level: advanced
335: Notes:
336: Any subset of the `x`, `y`, and `w` may be the same vector.
338: For complex numbers compares only the real part
340: .seealso: [](chapter_vectors), `Vec`, `VecPointwiseDivide()`, `VecPointwiseMult()`, `VecPointwiseMin()`, `VecPointwiseMaxAbs()`, `VecMaxPointwiseDivide()`
341: @*/
342: PetscErrorCode VecPointwiseMax(Vec w, Vec x, Vec y)
343: {
344: PetscFunctionBegin;
346: // REVIEW ME: no log event?
347: PetscCall(VecPointwiseApply_Private(w, x, y, 0, w->ops->pointwisemax));
348: PetscFunctionReturn(PETSC_SUCCESS);
349: }
351: /*@
352: VecPointwiseMin - Computes the component-wise minimum `w[i] = min(x[i], y[i])`.
354: Logically Collective
356: Input Parameters:
357: + x - the first input vector
358: - y - the second input vector
360: Output Parameter:
361: . w - the result
363: Level: advanced
365: Notes:
366: Any subset of the `x`, `y`, and `w` may be the same vector.
368: For complex numbers compares only the real part
370: .seealso: [](chapter_vectors), `Vec`, `VecPointwiseDivide()`, `VecPointwiseMult()`, `VecPointwiseMin()`, `VecPointwiseMaxAbs()`, `VecMaxPointwiseDivide()`
371: @*/
372: PetscErrorCode VecPointwiseMin(Vec w, Vec x, Vec y)
373: {
374: PetscFunctionBegin;
376: VecCheckAssembled(x);
377: // REVIEW ME: no log event?
378: PetscCall(VecPointwiseApply_Private(w, x, y, 0, w->ops->pointwisemin));
379: PetscFunctionReturn(PETSC_SUCCESS);
380: }
382: /*@
383: VecPointwiseMaxAbs - Computes the component-wise maximum of the absolute values `w[i] = max(abs(x[i]), abs(y[i]))`.
385: Logically Collective
387: Input Parameters:
388: + x - the first input vector
389: - y - the second input vector
391: Output Parameter:
392: . w - the result
394: Level: advanced
396: Notes:
397: Any subset of the `x`, `y`, and `w` may be the same vector.
399: .seealso: [](chapter_vectors), `Vec`, `VecPointwiseDivide()`, `VecPointwiseMult()`, `VecPointwiseMin()`, `VecPointwiseMax()`, `VecMaxPointwiseDivide()`
400: @*/
401: PetscErrorCode VecPointwiseMaxAbs(Vec w, Vec x, Vec y)
402: {
403: PetscFunctionBegin;
405: // REVIEW ME: no log event?
406: PetscCall(VecPointwiseApply_Private(w, x, y, 0, w->ops->pointwisemaxabs));
407: PetscFunctionReturn(PETSC_SUCCESS);
408: }
410: /*@
411: VecPointwiseDivide - Computes the component-wise division `w[i] = x[i] / y[i]`.
413: Logically Collective
415: Input Parameters:
416: + x - the numerator vector
417: - y - the denominator vector
419: Output Parameter:
420: . w - the result
422: Level: advanced
424: Note:
425: Any subset of the `x`, `y`, and `w` may be the same vector.
427: .seealso: [](chapter_vectors), `Vec`, `VecPointwiseMult()`, `VecPointwiseMax()`, `VecPointwiseMin()`, `VecPointwiseMaxAbs()`, `VecMaxPointwiseDivide()`
428: @*/
429: PetscErrorCode VecPointwiseDivide(Vec w, Vec x, Vec y)
430: {
431: PetscFunctionBegin;
433: // REVIEW ME: no log event?
434: PetscCall(VecPointwiseApply_Private(w, x, y, 0, w->ops->pointwisedivide));
435: PetscFunctionReturn(PETSC_SUCCESS);
436: }
438: /*@
439: VecPointwiseMult - Computes the component-wise multiplication `w[i] = x[i] * y[i]`.
441: Logically Collective
443: Input Parameters:
444: . x, y - the vectors
446: Output Parameter:
447: . w - the result
449: Level: advanced
451: Note:
452: Any subset of the `x`, `y`, and `w` may be the same vector.
454: .seealso: [](chapter_vectors), `Vec`, `VecPointwiseDivide()`, `VecPointwiseMax()`, `VecPointwiseMin()`, `VecPointwiseMaxAbs()`, `VecMaxPointwiseDivide()`
455: @*/
456: PetscErrorCode VecPointwiseMult(Vec w, Vec x, Vec y)
457: {
458: PetscFunctionBegin;
460: PetscCall(VecPointwiseApply_Private(w, x, y, VEC_PointwiseMult, w->ops->pointwisemult));
461: PetscFunctionReturn(PETSC_SUCCESS);
462: }
464: /*@
465: VecDuplicate - Creates a new vector of the same type as an existing vector.
467: Collective
469: Input Parameters:
470: . v - a vector to mimic
472: Output Parameter:
473: . newv - location to put new vector
475: Level: beginner
477: Notes:
478: VecDuplicate() DOES NOT COPY the vector entries, but rather allocates storage
479: for the new vector. Use `VecCopy()` to copy a vector.
481: Use `VecDestroy()` to free the space. Use `VecDuplicateVecs()` to get several
482: vectors.
484: .seealso: [](chapter_vectors), `Vec`, `VecDestroy()`, `VecDuplicateVecs()`, `VecCreate()`, `VecCopy()`
485: @*/
486: PetscErrorCode VecDuplicate(Vec v, Vec *newv)
487: {
488: PetscFunctionBegin;
492: PetscUseTypeMethod(v, duplicate, newv);
493: #if PetscDefined(HAVE_DEVICE)
494: if (v->boundtocpu && v->bindingpropagates) {
495: PetscCall(VecSetBindingPropagates(*newv, PETSC_TRUE));
496: PetscCall(VecBindToCPU(*newv, PETSC_TRUE));
497: }
498: #endif
499: PetscCall(PetscObjectStateIncrease((PetscObject)(*newv)));
500: PetscFunctionReturn(PETSC_SUCCESS);
501: }
503: /*@C
504: VecDestroy - Destroys a vector.
506: Collective
508: Input Parameters:
509: . v - the vector
511: Level: beginner
513: .seealso: [](chapter_vectors), `Vec`, `VecDuplicate()`, `VecDestroyVecs()`
514: @*/
515: PetscErrorCode VecDestroy(Vec *v)
516: {
517: PetscFunctionBegin;
519: if (!*v) PetscFunctionReturn(PETSC_SUCCESS);
521: if (--((PetscObject)(*v))->refct > 0) {
522: *v = NULL;
523: PetscFunctionReturn(PETSC_SUCCESS);
524: }
526: PetscCall(PetscObjectSAWsViewOff((PetscObject)*v));
527: /* destroy the internal part */
528: PetscTryTypeMethod(*v, destroy);
529: PetscCall(PetscFree((*v)->defaultrandtype));
530: /* destroy the external/common part */
531: PetscCall(PetscLayoutDestroy(&(*v)->map));
532: PetscCall(PetscHeaderDestroy(v));
533: PetscFunctionReturn(PETSC_SUCCESS);
534: }
536: /*@C
537: VecDuplicateVecs - Creates several vectors of the same type as an existing vector.
539: Collective
541: Input Parameters:
542: + m - the number of vectors to obtain
543: - v - a vector to mimic
545: Output Parameter:
546: . V - location to put pointer to array of vectors
548: Level: intermediate
550: Note:
551: Use `VecDestroyVecs()` to free the space. Use `VecDuplicate()` to form a single
552: vector.
554: Fortran Note:
555: The Fortran interface is slightly different from that given below, it
556: requires one to pass in V a `Vec` array of size at least m.
557: See the [](chapter_fortran) for details.
559: .seealso: [](chapter_vectors), `Vec`, [](chapter_fortran), `VecDestroyVecs()`, `VecDuplicate()`, `VecCreate()`, `VecDuplicateVecsF90()`
560: @*/
561: PetscErrorCode VecDuplicateVecs(Vec v, PetscInt m, Vec *V[])
562: {
563: PetscFunctionBegin;
567: PetscUseTypeMethod(v, duplicatevecs, m, V);
568: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA) || defined(PETSC_HAVE_HIP)
569: if (v->boundtocpu && v->bindingpropagates) {
570: PetscInt i;
572: for (i = 0; i < m; i++) {
573: /* Since ops->duplicatevecs might itself propagate the value of boundtocpu,
574: * avoid unnecessary overhead by only calling VecBindToCPU() if the vector isn't already bound. */
575: if (!(*V)[i]->boundtocpu) {
576: PetscCall(VecSetBindingPropagates((*V)[i], PETSC_TRUE));
577: PetscCall(VecBindToCPU((*V)[i], PETSC_TRUE));
578: }
579: }
580: }
581: #endif
582: PetscFunctionReturn(PETSC_SUCCESS);
583: }
585: /*@C
586: VecDestroyVecs - Frees a block of vectors obtained with `VecDuplicateVecs()`.
588: Collective
590: Input Parameters:
591: + vv - pointer to pointer to array of vector pointers, if `NULL` no vectors are destroyed
592: - m - the number of vectors previously obtained, if zero no vectors are destroyed
594: Level: intermediate
596: Fortran Note:
597: The Fortran interface is slightly different from that given below.
598: See the [](chapter_fortran) for details.
600: .seealso: [](chapter_vectors), `Vec`, [](chapter_fortran), `VecDuplicateVecs()`, `VecDestroyVecsf90()`
601: @*/
602: PetscErrorCode VecDestroyVecs(PetscInt m, Vec *vv[])
603: {
604: PetscFunctionBegin;
606: PetscCheck(m >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Trying to destroy negative number of vectors %" PetscInt_FMT, m);
607: if (!m || !*vv) {
608: *vv = NULL;
609: PetscFunctionReturn(PETSC_SUCCESS);
610: }
613: PetscCall((*(**vv)->ops->destroyvecs)(m, *vv));
614: *vv = NULL;
615: PetscFunctionReturn(PETSC_SUCCESS);
616: }
618: /*@C
619: VecViewFromOptions - View a vector based on values in the options database
621: Collective
623: Input Parameters:
624: + A - the vector
625: . obj - Optional object that provides the options prefix for this viewing
626: - name - command line option
628: Level: intermediate
630: .seealso: [](chapter_vectors), `Vec`, `VecView`, `PetscObjectViewFromOptions()`, `VecCreate()`
631: @*/
632: PetscErrorCode VecViewFromOptions(Vec A, PetscObject obj, const char name[])
633: {
634: PetscFunctionBegin;
636: PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name));
637: PetscFunctionReturn(PETSC_SUCCESS);
638: }
640: /*@C
641: VecView - Views a vector object.
643: Collective
645: Input Parameters:
646: + vec - the vector
647: - viewer - an optional visualization context
649: Level: beginner
651: Notes:
652: The available visualization contexts include
653: + `PETSC_VIEWER_STDOUT_SELF` - for sequential vectors
654: . `PETSC_VIEWER_STDOUT_WORLD` - for parallel vectors created on `PETSC_COMM_WORLD`
655: - `PETSC_VIEWER_STDOUT`_(comm) - for parallel vectors created on MPI communicator comm
657: You can change the format the vector is printed using the
658: option `PetscViewerPushFormat()`.
660: The user can open alternative viewers with
661: + `PetscViewerASCIIOpen()` - Outputs vector to a specified file
662: . `PetscViewerBinaryOpen()` - Outputs vector in binary to a
663: specified file; corresponding input uses `VecLoad()`
664: . `PetscViewerDrawOpen()` - Outputs vector to an X window display
665: . `PetscViewerSocketOpen()` - Outputs vector to Socket viewer
666: - `PetscViewerHDF5Open()` - Outputs vector to HDF5 file viewer
668: The user can call `PetscViewerPushFormat()` to specify the output
669: format of ASCII printed objects (when using `PETSC_VIEWER_STDOUT_SELF`,
670: `PETSC_VIEWER_STDOUT_WORLD` and `PetscViewerASCIIOpen()`). Available formats include
671: + `PETSC_VIEWER_DEFAULT` - default, prints vector contents
672: . `PETSC_VIEWER_ASCII_MATLAB` - prints vector contents in MATLAB format
673: . `PETSC_VIEWER_ASCII_INDEX` - prints vector contents, including indices of vector elements
674: - `PETSC_VIEWER_ASCII_COMMON` - prints vector contents, using a
675: format common among all vector types
677: You can pass any number of vector objects, or other PETSc objects to the same viewer.
679: In the debugger you can do call `VecView`(v,0) to display the vector. (The same holds for any PETSc object viewer).
681: Notes for binary viewer:
682: If you pass multiple vectors to a binary viewer you can read them back in in the same order
683: with `VecLoad()`.
685: If the blocksize of the vector is greater than one then you must provide a unique prefix to
686: the vector with `PetscObjectSetOptionsPrefix`((`PetscObject`)vec,"uniqueprefix"); BEFORE calling `VecView()` on the
687: vector to be stored and then set that same unique prefix on the vector that you pass to `VecLoad()`. The blocksize
688: information is stored in an ASCII file with the same name as the binary file plus a ".info" appended to the
689: filename. If you copy the binary file, make sure you copy the associated .info file with it.
691: See the manual page for `VecLoad()` on the exact format the binary viewer stores
692: the values in the file.
694: Notes for HDF5 Viewer:
695: The name of the `Vec` (given with `PetscObjectSetName()` is the name that is used
696: for the object in the HDF5 file. If you wish to store the same Vec into multiple
697: datasets in the same file (typically with different values), you must change its
698: name each time before calling the `VecView()`. To load the same vector,
699: the name of the Vec object passed to `VecLoad()` must be the same.
701: If the block size of the vector is greater than 1 then it is used as the first dimension in the HDF5 array.
702: If the function `PetscViewerHDF5SetBaseDimension2()`is called then even if the block size is one it will
703: be used as the first dimension in the HDF5 array (that is the HDF5 array will always be two dimensional)
704: See also `PetscViewerHDF5SetTimestep()` which adds an additional complication to reading and writing `Vec`
705: with the HDF5 viewer.
707: .seealso: [](chapter_vectors), `Vec`, `PetscViewerASCIIOpen()`, `PetscViewerDrawOpen()`, `PetscDrawLGCreate()`,
708: `PetscViewerSocketOpen()`, `PetscViewerBinaryOpen()`, `VecLoad()`, `PetscViewerCreate()`,
709: `PetscRealView()`, `PetscScalarView()`, `PetscIntView()`, `PetscViewerHDF5SetTimestep()`
710: @*/
711: PetscErrorCode VecView(Vec vec, PetscViewer viewer)
712: {
713: PetscBool iascii;
714: PetscViewerFormat format;
715: PetscMPIInt size;
717: PetscFunctionBegin;
720: VecCheckAssembled(vec);
721: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)vec), &viewer));
723: PetscCall(PetscViewerGetFormat(viewer, &format));
724: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)vec), &size));
725: if (size == 1 && format == PETSC_VIEWER_LOAD_BALANCE) PetscFunctionReturn(PETSC_SUCCESS);
727: PetscCheck(!vec->stash.n && !vec->bstash.n, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Must call VecAssemblyBegin/End() before viewing this vector");
729: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
730: if (iascii) {
731: PetscInt rows, bs;
733: PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)vec, viewer));
734: if (format == PETSC_VIEWER_ASCII_INFO || format == PETSC_VIEWER_ASCII_INFO_DETAIL) {
735: PetscCall(PetscViewerASCIIPushTab(viewer));
736: PetscCall(VecGetSize(vec, &rows));
737: PetscCall(VecGetBlockSize(vec, &bs));
738: if (bs != 1) {
739: PetscCall(PetscViewerASCIIPrintf(viewer, "length=%" PetscInt_FMT ", bs=%" PetscInt_FMT "\n", rows, bs));
740: } else {
741: PetscCall(PetscViewerASCIIPrintf(viewer, "length=%" PetscInt_FMT "\n", rows));
742: }
743: PetscCall(PetscViewerASCIIPopTab(viewer));
744: }
745: }
746: PetscCall(VecLockReadPush(vec));
747: PetscCall(PetscLogEventBegin(VEC_View, vec, viewer, 0, 0));
748: if ((format == PETSC_VIEWER_NATIVE || format == PETSC_VIEWER_LOAD_BALANCE) && vec->ops->viewnative) {
749: PetscUseTypeMethod(vec, viewnative, viewer);
750: } else {
751: PetscUseTypeMethod(vec, view, viewer);
752: }
753: PetscCall(VecLockReadPop(vec));
754: PetscCall(PetscLogEventEnd(VEC_View, vec, viewer, 0, 0));
755: PetscFunctionReturn(PETSC_SUCCESS);
756: }
758: #if defined(PETSC_USE_DEBUG)
759: #include <../src/sys/totalview/tv_data_display.h>
760: PETSC_UNUSED static int TV_display_type(const struct _p_Vec *v)
761: {
762: const PetscScalar *values;
763: char type[32];
765: TV_add_row("Local rows", "int", &v->map->n);
766: TV_add_row("Global rows", "int", &v->map->N);
767: TV_add_row("Typename", TV_ascii_string_type, ((PetscObject)v)->type_name);
768: PetscCall(VecGetArrayRead((Vec)v, &values));
769: PetscCall(PetscSNPrintf(type, 32, "double[%" PetscInt_FMT "]", v->map->n));
770: TV_add_row("values", type, values);
771: PetscCall(VecRestoreArrayRead((Vec)v, &values));
772: return TV_format_OK;
773: }
774: #endif
776: /*@C
777: VecViewNative - Views a vector object with the original type specific viewer
779: Collective
781: Input Parameters:
782: + vec - the vector
783: - viewer - an optional visualization context
785: Level: developer
787: .seealso: [](chapter_vectors), `Vec`, `PetscViewerASCIIOpen()`, `PetscViewerDrawOpen()`, `PetscDrawLGCreate()`, `VecView()`
788: `PetscViewerSocketOpen()`, `PetscViewerBinaryOpen()`, `VecLoad()`, `PetscViewerCreate()`,
789: `PetscRealView()`, `PetscScalarView()`, `PetscIntView()`, `PetscViewerHDF5SetTimestep()`
790: @*/
791: PetscErrorCode VecViewNative(Vec vec, PetscViewer viewer)
792: {
793: PetscFunctionBegin;
796: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)vec), &viewer));
798: PetscUseTypeMethod(vec, viewnative, viewer);
799: PetscFunctionReturn(PETSC_SUCCESS);
800: }
802: /*@
803: VecGetSize - Returns the global number of elements of the vector.
805: Not Collective
807: Input Parameter:
808: . x - the vector
810: Output Parameters:
811: . size - the global length of the vector
813: Level: beginner
815: .seealso: [](chapter_vectors), `Vec`, `VecGetLocalSize()`
816: @*/
817: PetscErrorCode VecGetSize(Vec x, PetscInt *size)
818: {
819: PetscFunctionBegin;
823: PetscUseTypeMethod(x, getsize, size);
824: PetscFunctionReturn(PETSC_SUCCESS);
825: }
827: /*@
828: VecGetLocalSize - Returns the number of elements of the vector stored
829: in local memory.
831: Not Collective
833: Input Parameter:
834: . x - the vector
836: Output Parameter:
837: . size - the length of the local piece of the vector
839: Level: beginner
841: .seealso: [](chapter_vectors), `Vec`, `VecGetSize()`
842: @*/
843: PetscErrorCode VecGetLocalSize(Vec x, PetscInt *size)
844: {
845: PetscFunctionBegin;
849: PetscUseTypeMethod(x, getlocalsize, size);
850: PetscFunctionReturn(PETSC_SUCCESS);
851: }
853: /*@C
854: VecGetOwnershipRange - Returns the range of indices owned by
855: this processor, assuming that the vectors are laid out with the
856: first n1 elements on the first processor, next n2 elements on the
857: second, etc. For certain parallel layouts this range may not be
858: well defined.
860: Not Collective
862: Input Parameter:
863: . x - the vector
865: Output Parameters:
866: + low - the first local element, pass in `NULL` if not interested
867: - high - one more than the last local element, pass in `NULL` if not interested
869: Level: beginner
871: Note:
872: The high argument is one more than the last element stored locally.
874: Fortran Note:
875: `PETSC_NULL_INTEGER` should be used instead of NULL
877: .seealso: [](chapter_vectors), `Vec`, `MatGetOwnershipRange()`, `MatGetOwnershipRanges()`, `VecGetOwnershipRanges()`
878: @*/
879: PetscErrorCode VecGetOwnershipRange(Vec x, PetscInt *low, PetscInt *high)
880: {
881: PetscFunctionBegin;
886: if (low) *low = x->map->rstart;
887: if (high) *high = x->map->rend;
888: PetscFunctionReturn(PETSC_SUCCESS);
889: }
891: /*@C
892: VecGetOwnershipRanges - Returns the range of indices owned by EACH processor,
893: assuming that the vectors are laid out with the
894: first n1 elements on the first processor, next n2 elements on the
895: second, etc. For certain parallel layouts this range may not be
896: well defined.
898: Not Collective
900: Input Parameter:
901: . x - the vector
903: Output Parameters:
904: . range - array of length size+1 with the start and end+1 for each process
906: Level: beginner
908: Notes:
909: The high argument is one more than the last element stored locally.
911: If the ranges are used after all vectors that share the ranges has been destroyed then the program will crash accessing ranges[].
913: Fortran Note:
914: You must PASS in an array of length size+1
916: .seealso: [](chapter_vectors), `Vec`, `MatGetOwnershipRange()`, `MatGetOwnershipRanges()`, `VecGetOwnershipRange()`
917: @*/
918: PetscErrorCode VecGetOwnershipRanges(Vec x, const PetscInt *ranges[])
919: {
920: PetscFunctionBegin;
923: PetscCall(PetscLayoutGetRanges(x->map, ranges));
924: PetscFunctionReturn(PETSC_SUCCESS);
925: }
927: /*@
928: VecSetOption - Sets an option for controlling a vector's behavior.
930: Collective
932: Input Parameters:
933: + x - the vector
934: . op - the option
935: - flag - turn the option on or off
937: Supported Options:
938: + `VEC_IGNORE_OFF_PROC_ENTRIES`, which causes `VecSetValues()` to ignore
939: entries destined to be stored on a separate processor. This can be used
940: to eliminate the global reduction in the `VecAssemblyBegin()` if you know
941: that you have only used `VecSetValues()` to set local elements
942: . `VEC_IGNORE_NEGATIVE_INDICES`, which means you can pass negative indices
943: in ix in calls to `VecSetValues()` or `VecGetValues()`. These rows are simply
944: ignored.
945: - `VEC_SUBSET_OFF_PROC_ENTRIES`, which causes `VecAssemblyBegin()` to assume that the off-process
946: entries will always be a subset (possibly equal) of the off-process entries set on the
947: first assembly which had a true `VEC_SUBSET_OFF_PROC_ENTRIES` and the vector has not
948: changed this flag afterwards. If this assembly is not such first assembly, then this
949: assembly can reuse the communication pattern setup in that first assembly, thus avoiding
950: a global reduction. Subsequent assemblies setting off-process values should use the same
951: InsertMode as the first assembly.
953: Level: intermediate
955: Developer Note:
956: The `InsertMode` restriction could be removed by packing the stash messages out of place.
958: .seealso: [](chapter_vectors), `Vec`, `VecSetValues()`
959: @*/
960: PetscErrorCode VecSetOption(Vec x, VecOption op, PetscBool flag)
961: {
962: PetscFunctionBegin;
965: PetscTryTypeMethod(x, setoption, op, flag);
966: PetscFunctionReturn(PETSC_SUCCESS);
967: }
969: /* Default routines for obtaining and releasing; */
970: /* may be used by any implementation */
971: PetscErrorCode VecDuplicateVecs_Default(Vec w, PetscInt m, Vec *V[])
972: {
973: PetscFunctionBegin;
976: PetscCheck(m > 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "m must be > 0: m = %" PetscInt_FMT, m);
977: PetscCall(PetscMalloc1(m, V));
978: for (PetscInt i = 0; i < m; i++) PetscCall(VecDuplicate(w, *V + i));
979: PetscFunctionReturn(PETSC_SUCCESS);
980: }
982: PetscErrorCode VecDestroyVecs_Default(PetscInt m, Vec v[])
983: {
984: PetscInt i;
986: PetscFunctionBegin;
988: for (i = 0; i < m; i++) PetscCall(VecDestroy(&v[i]));
989: PetscCall(PetscFree(v));
990: PetscFunctionReturn(PETSC_SUCCESS);
991: }
993: /*@
994: VecResetArray - Resets a vector to use its default memory. Call this
995: after the use of `VecPlaceArray()`.
997: Not Collective
999: Input Parameters:
1000: . vec - the vector
1002: Level: developer
1004: .seealso: [](chapter_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecReplaceArray()`, `VecPlaceArray()`
1005: @*/
1006: PetscErrorCode VecResetArray(Vec vec)
1007: {
1008: PetscFunctionBegin;
1011: PetscUseTypeMethod(vec, resetarray);
1012: PetscCall(PetscObjectStateIncrease((PetscObject)vec));
1013: PetscFunctionReturn(PETSC_SUCCESS);
1014: }
1016: /*@C
1017: VecLoad - Loads a vector that has been stored in binary or HDF5 format
1018: with `VecView()`.
1020: Collective
1022: Input Parameters:
1023: + vec - the newly loaded vector, this needs to have been created with `VecCreate()` or
1024: some related function before a call to VecLoad().
1025: - viewer - binary file viewer, obtained from PetscViewerBinaryOpen() or
1026: HDF5 file viewer, obtained from PetscViewerHDF5Open()
1028: Level: intermediate
1030: Notes:
1031: Defaults to the standard Seq or MPI Vec, if you want some other type of `Vec` call `VecSetFromOptions()`
1032: before calling this.
1034: The input file must contain the full global vector, as
1035: written by the routine `VecView()`.
1037: If the type or size of vec is not set before a call to `VecLoad()`, PETSc
1038: sets the type and the local and global sizes. If type and/or
1039: sizes are already set, then the same are used.
1041: If using the binary viewer and the blocksize of the vector is greater than one then you must provide a unique prefix to
1042: the vector with `PetscObjectSetOptionsPrefix`((`PetscObject`)vec,"uniqueprefix"); BEFORE calling `VecView()` on the
1043: vector to be stored and then set that same unique prefix on the vector that you pass to VecLoad(). The blocksize
1044: information is stored in an ASCII file with the same name as the binary file plus a ".info" appended to the
1045: filename. If you copy the binary file, make sure you copy the associated .info file with it.
1047: If using HDF5, you must assign the Vec the same name as was used in the Vec
1048: that was stored in the file using `PetscObjectSetName(). Otherwise you will
1049: get the error message: "Cannot H5DOpen2() with `Vec` name NAMEOFOBJECT".
1051: If the HDF5 file contains a two dimensional array the first dimension is treated as the block size
1052: in loading the vector. Hence, for example, using MATLAB notation h5create('vector.dat','/Test_Vec',[27 1]);
1053: will load a vector of size 27 and block size 27 thus resulting in all 27 entries being on the first process of
1054: vectors communicator and the rest of the processes having zero entries
1056: Notes for advanced users when using the binary viewer:
1057: Most users should not need to know the details of the binary storage
1058: format, since `VecLoad()` and `VecView()` completely hide these details.
1059: But for anyone who's interested, the standard binary vector storage
1060: format is
1061: .vb
1062: PetscInt VEC_FILE_CLASSID
1063: PetscInt number of rows
1064: PetscScalar *values of all entries
1065: .ve
1067: In addition, PETSc automatically uses byte swapping to work on all machines; the files
1068: are written ALWAYS using big-endian ordering. On small-endian machines the numbers
1069: are converted to the small-endian format when they are read in from the file.
1070: See PetscBinaryRead() and PetscBinaryWrite() to see how this may be done.
1072: .seealso: [](chapter_vectors), `Vec`, `PetscViewerBinaryOpen()`, `VecView()`, `MatLoad()`, `VecLoad()`
1073: @*/
1074: PetscErrorCode VecLoad(Vec vec, PetscViewer viewer)
1075: {
1076: PetscBool isbinary, ishdf5, isadios, isexodusii;
1077: PetscViewerFormat format;
1079: PetscFunctionBegin;
1082: PetscCheckSameComm(vec, 1, viewer, 2);
1083: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary));
1084: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERHDF5, &ishdf5));
1085: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERADIOS, &isadios));
1086: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWEREXODUSII, &isexodusii));
1087: PetscCheck(isbinary || ishdf5 || isadios || isexodusii, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid viewer; open viewer with PetscViewerBinaryOpen()");
1089: PetscCall(VecSetErrorIfLocked(vec, 1));
1090: if (!((PetscObject)vec)->type_name && !vec->ops->create) PetscCall(VecSetType(vec, VECSTANDARD));
1091: PetscCall(PetscLogEventBegin(VEC_Load, viewer, 0, 0, 0));
1092: PetscCall(PetscViewerGetFormat(viewer, &format));
1093: if (format == PETSC_VIEWER_NATIVE && vec->ops->loadnative) {
1094: PetscUseTypeMethod(vec, loadnative, viewer);
1095: } else {
1096: PetscUseTypeMethod(vec, load, viewer);
1097: }
1098: PetscCall(PetscLogEventEnd(VEC_Load, viewer, 0, 0, 0));
1099: PetscFunctionReturn(PETSC_SUCCESS);
1100: }
1102: /*@
1103: VecReciprocal - Replaces each component of a vector by its reciprocal.
1105: Logically Collective
1107: Input Parameter:
1108: . vec - the vector
1110: Output Parameter:
1111: . vec - the vector reciprocal
1113: Level: intermediate
1115: .seealso: [](chapter_vectors), `Vec`, `VecLog()`, `VecExp()`, `VecSqrtAbs()`
1116: @*/
1117: PetscErrorCode VecReciprocal(Vec vec)
1118: {
1119: PetscFunctionBegin;
1122: VecCheckAssembled(vec);
1123: PetscCall(VecSetErrorIfLocked(vec, 1));
1124: PetscUseTypeMethod(vec, reciprocal);
1125: PetscCall(PetscObjectStateIncrease((PetscObject)vec));
1126: PetscFunctionReturn(PETSC_SUCCESS);
1127: }
1129: /*@C
1130: VecSetOperation - Allows the user to override a particular vector operation.
1132: Logically Collective; No Fortran Support
1134: Input Parameters:
1135: + vec - The vector to modify
1136: . op - The name of the operation
1137: - f - The function that provides the operation.
1139: Notes:
1140: `f` may be `NULL` to remove the operation from `vec`. Depending on the operation this may be
1141: allowed, however some always expect a valid function. In these cases an error will be raised
1142: when calling the interface routine in question.
1144: See `VecOperation` for an up-to-date list of override-able operations. The operations listed
1145: there have the form `VECOP_<OPERATION>`, where `<OPERATION>` is the suffix (in all capital
1146: letters) of the public interface routine (e.g., `VecView()` -> `VECOP_VIEW`).
1148: Overriding a particular `Vec`'s operation has no affect on any other `Vec`s past, present,
1149: or future. The user should also note that overriding a method is "destructive"; the previous
1150: method is not retained in any way.
1152: Level: advanced
1154: Example Usage:
1155: .vb
1156: // some new VecView() implementation, must have the same signature as the function it seeks
1157: // to replace
1158: PetscErrorCode UserVecView(Vec x, PetscViewer viewer)
1159: {
1160: PetscFunctionBeginUser;
1161: // ...
1162: PetscFunctionReturn(PETSC_SUCCESS);
1163: }
1165: // Create a VECMPI which has a pre-defined VecView() implementation
1166: VecCreateMPI(comm, n, N, &x);
1167: // Calls the VECMPI implementation for VecView()
1168: VecView(x, viewer);
1170: VecSetOperation(x, VECOP_VIEW, (void (*)(void))UserVecView);
1171: // Now calls UserVecView()
1172: VecView(x, viewer);
1173: .ve
1175: .seealso: [](chapter_vectors), `Vec`, `VecCreate()`, `MatShellSetOperation()`
1176: @*/
1177: PetscErrorCode VecSetOperation(Vec vec, VecOperation op, void (*f)(void))
1178: {
1179: PetscFunctionBegin;
1181: if (op == VECOP_VIEW && !vec->ops->viewnative) {
1182: vec->ops->viewnative = vec->ops->view;
1183: } else if (op == VECOP_LOAD && !vec->ops->loadnative) {
1184: vec->ops->loadnative = vec->ops->load;
1185: }
1186: ((void (**)(void))vec->ops)[(int)op] = f;
1187: PetscFunctionReturn(PETSC_SUCCESS);
1188: }
1190: /*@
1191: VecStashSetInitialSize - sets the sizes of the vec-stash, that is
1192: used during the assembly process to store values that belong to
1193: other processors.
1195: Not Collective, different processes can have different size stashes
1197: Input Parameters:
1198: + vec - the vector
1199: . size - the initial size of the stash.
1200: - bsize - the initial size of the block-stash(if used).
1202: Options Database Keys:
1203: + -vecstash_initial_size <size> or <size0,size1,...sizep-1>
1204: - -vecstash_block_initial_size <bsize> or <bsize0,bsize1,...bsizep-1>
1206: Level: intermediate
1208: Notes:
1209: The block-stash is used for values set with `VecSetValuesBlocked()` while
1210: the stash is used for values set with `VecSetValues()`
1212: Run with the option -info and look for output of the form
1213: VecAssemblyBegin_MPIXXX:Stash has MM entries, uses nn mallocs.
1214: to determine the appropriate value, MM, to use for size and
1215: VecAssemblyBegin_MPIXXX:Block-Stash has BMM entries, uses nn mallocs.
1216: to determine the value, BMM to use for bsize
1218: .seealso: [](chapter_vectors), `Vec`, `VecSetBlockSize()`, `VecSetValues()`, `VecSetValuesBlocked()`, `VecStashView()`
1219: @*/
1220: PetscErrorCode VecStashSetInitialSize(Vec vec, PetscInt size, PetscInt bsize)
1221: {
1222: PetscFunctionBegin;
1224: PetscCall(VecStashSetInitialSize_Private(&vec->stash, size));
1225: PetscCall(VecStashSetInitialSize_Private(&vec->bstash, bsize));
1226: PetscFunctionReturn(PETSC_SUCCESS);
1227: }
1229: /*@
1230: VecConjugate - Conjugates a vector.
1232: Logically Collective
1234: Input Parameters:
1235: . x - the vector
1237: Level: intermediate
1239: .seealso: [](chapter_vectors), `Vec`, `VecSet()`
1240: @*/
1241: PetscErrorCode VecConjugate(Vec x)
1242: {
1243: PetscFunctionBegin;
1246: VecCheckAssembled(x);
1247: PetscCall(VecSetErrorIfLocked(x, 1));
1248: if (PetscDefined(USE_COMPLEX)) {
1249: PetscUseTypeMethod(x, conjugate);
1250: /* we need to copy norms here */
1251: PetscCall(PetscObjectStateIncrease((PetscObject)x));
1252: }
1253: PetscFunctionReturn(PETSC_SUCCESS);
1254: }
1256: /*@
1257: VecSetRandom - Sets all components of a vector to random numbers.
1259: Logically Collective
1261: Input Parameters:
1262: + x - the vector
1263: - rctx - the random number context, formed by PetscRandomCreate(), or NULL and
1264: it will create one internally.
1266: Output Parameter:
1267: . x - the vector
1269: Example of Usage:
1270: .vb
1271: PetscRandomCreate(PETSC_COMM_WORLD,&rctx);
1272: VecSetRandom(x,rctx);
1273: PetscRandomDestroy(&rctx);
1274: .ve
1276: Level: intermediate
1278: .seealso: [](chapter_vectors), `Vec`, `VecSet()`, `VecSetValues()`, `PetscRandomCreate()`, `PetscRandomDestroy()`
1279: @*/
1280: PetscErrorCode VecSetRandom(Vec x, PetscRandom rctx)
1281: {
1282: PetscRandom randObj = NULL;
1284: PetscFunctionBegin;
1288: VecCheckAssembled(x);
1289: PetscCall(VecSetErrorIfLocked(x, 1));
1291: if (!rctx) {
1292: PetscCall(PetscRandomCreate(PetscObjectComm((PetscObject)x), &randObj));
1293: PetscCall(PetscRandomSetType(randObj, x->defaultrandtype));
1294: PetscCall(PetscRandomSetFromOptions(randObj));
1295: rctx = randObj;
1296: }
1298: PetscCall(PetscLogEventBegin(VEC_SetRandom, x, rctx, 0, 0));
1299: PetscUseTypeMethod(x, setrandom, rctx);
1300: PetscCall(PetscLogEventEnd(VEC_SetRandom, x, rctx, 0, 0));
1302: PetscCall(PetscRandomDestroy(&randObj));
1303: PetscCall(PetscObjectStateIncrease((PetscObject)x));
1304: PetscFunctionReturn(PETSC_SUCCESS);
1305: }
1307: /*@
1308: VecZeroEntries - puts a `0.0` in each element of a vector
1310: Logically Collective
1312: Input Parameter:
1313: . vec - The vector
1315: Level: beginner
1317: .seealso: [](chapter_vectors), `Vec`, `VecCreate()`, `VecSetOptionsPrefix()`, `VecSet()`, `VecSetValues()`
1318: @*/
1319: PetscErrorCode VecZeroEntries(Vec vec)
1320: {
1321: PetscFunctionBegin;
1322: PetscCall(VecSet(vec, 0));
1323: PetscFunctionReturn(PETSC_SUCCESS);
1324: }
1326: /*
1327: VecSetTypeFromOptions_Private - Sets the type of vector from user options. Defaults to a PETSc sequential vector on one
1328: processor and a PETSc MPI vector on more than one processor.
1330: Collective
1332: Input Parameter:
1333: . vec - The vector
1335: Level: intermediate
1337: .seealso: [](chapter_vectors), `Vec`, `VecSetFromOptions()`, `VecSetType()`
1338: */
1339: static PetscErrorCode VecSetTypeFromOptions_Private(Vec vec, PetscOptionItems *PetscOptionsObject)
1340: {
1341: PetscBool opt;
1342: VecType defaultType;
1343: char typeName[256];
1344: PetscMPIInt size;
1346: PetscFunctionBegin;
1347: if (((PetscObject)vec)->type_name) defaultType = ((PetscObject)vec)->type_name;
1348: else {
1349: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)vec), &size));
1350: if (size > 1) defaultType = VECMPI;
1351: else defaultType = VECSEQ;
1352: }
1354: PetscCall(VecRegisterAll());
1355: PetscCall(PetscOptionsFList("-vec_type", "Vector type", "VecSetType", VecList, defaultType, typeName, 256, &opt));
1356: if (opt) {
1357: PetscCall(VecSetType(vec, typeName));
1358: } else {
1359: PetscCall(VecSetType(vec, defaultType));
1360: }
1361: PetscFunctionReturn(PETSC_SUCCESS);
1362: }
1364: /*@
1365: VecSetFromOptions - Configures the vector from the options database.
1367: Collective
1369: Input Parameter:
1370: . vec - The vector
1372: Level: beginner
1374: Notes:
1375: To see all options, run your program with the -help option, or consult the users manual.
1377: Must be called after `VecCreate()` but before the vector is used.
1379: .seealso: [](chapter_vectors), `Vec`, `VecCreate()`, `VecSetOptionsPrefix()`
1380: @*/
1381: PetscErrorCode VecSetFromOptions(Vec vec)
1382: {
1383: PetscBool flg;
1384: PetscInt bind_below = 0;
1386: PetscFunctionBegin;
1389: PetscObjectOptionsBegin((PetscObject)vec);
1390: /* Handle vector type options */
1391: PetscCall(VecSetTypeFromOptions_Private(vec, PetscOptionsObject));
1393: /* Handle specific vector options */
1394: PetscTryTypeMethod(vec, setfromoptions, PetscOptionsObject);
1396: /* Bind to CPU if below a user-specified size threshold.
1397: * This perhaps belongs in the options for the GPU Vec types, but VecBindToCPU() does nothing when called on non-GPU types,
1398: * and putting it here makes is more maintainable than duplicating this for all. */
1399: PetscCall(PetscOptionsInt("-vec_bind_below", "Set the size threshold (in local entries) below which the Vec is bound to the CPU", "VecBindToCPU", bind_below, &bind_below, &flg));
1400: if (flg && vec->map->n < bind_below) PetscCall(VecBindToCPU(vec, PETSC_TRUE));
1402: /* process any options handlers added with PetscObjectAddOptionsHandler() */
1403: PetscCall(PetscObjectProcessOptionsHandlers((PetscObject)vec, PetscOptionsObject));
1404: PetscOptionsEnd();
1405: PetscFunctionReturn(PETSC_SUCCESS);
1406: }
1408: /*@
1409: VecSetSizes - Sets the local and global sizes, and checks to determine compatibility
1411: Collective
1413: Input Parameters:
1414: + v - the vector
1415: . n - the local size (or `PETSC_DECIDE` to have it set)
1416: - N - the global size (or `PETSC_DETERMINE` to have it set)
1418: Level: intermediate
1420: Notes:
1421: N cannot be `PETSC_DETERMINE` if n is `PETSC_DECIDE`
1423: If one processor calls this with N of `PETSC_DETERMINE` then all processors must, otherwise the program will hang.
1425: .seealso: [](chapter_vectors), `Vec`, `VecGetSize()`, `PetscSplitOwnership()`
1426: @*/
1427: PetscErrorCode VecSetSizes(Vec v, PetscInt n, PetscInt N)
1428: {
1429: PetscFunctionBegin;
1431: if (N >= 0) {
1433: PetscCheck(n <= N, PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Local size %" PetscInt_FMT " cannot be larger than global size %" PetscInt_FMT, n, N);
1434: }
1435: PetscCheck(!(v->map->n >= 0 || v->map->N >= 0) || !(v->map->n != n || v->map->N != N), PETSC_COMM_SELF, PETSC_ERR_SUP, "Cannot change/reset vector sizes to %" PetscInt_FMT " local %" PetscInt_FMT " global after previously setting them to %" PetscInt_FMT " local %" PetscInt_FMT " global", n, N,
1436: v->map->n, v->map->N);
1437: v->map->n = n;
1438: v->map->N = N;
1439: PetscTryTypeMethod(v, create);
1440: v->ops->create = NULL;
1441: PetscFunctionReturn(PETSC_SUCCESS);
1442: }
1444: /*@
1445: VecSetBlockSize - Sets the block size for future calls to `VecSetValuesBlocked()`
1446: and `VecSetValuesBlockedLocal()`.
1448: Logically Collective
1450: Input Parameters:
1451: + v - the vector
1452: - bs - the blocksize
1454: Level: advanced
1456: Note:
1457: All vectors obtained by `VecDuplicate()` inherit the same blocksize.
1459: .seealso: [](chapter_vectors), `Vec`, `VecSetValuesBlocked()`, `VecSetLocalToGlobalMapping()`, `VecGetBlockSize()`
1460: @*/
1461: PetscErrorCode VecSetBlockSize(Vec v, PetscInt bs)
1462: {
1463: PetscFunctionBegin;
1466: PetscCall(PetscLayoutSetBlockSize(v->map, bs));
1467: v->bstash.bs = bs; /* use the same blocksize for the vec's block-stash */
1468: PetscFunctionReturn(PETSC_SUCCESS);
1469: }
1471: /*@
1472: VecGetBlockSize - Gets the blocksize for the vector, i.e. what is used for `VecSetValuesBlocked()`
1473: and `VecSetValuesBlockedLocal()`.
1475: Not Collective
1477: Input Parameter:
1478: . v - the vector
1480: Output Parameter:
1481: . bs - the blocksize
1483: Level: advanced
1485: Note:
1486: All vectors obtained by `VecDuplicate()` inherit the same blocksize.
1488: .seealso: [](chapter_vectors), `Vec`, `VecSetValuesBlocked()`, `VecSetLocalToGlobalMapping()`, `VecSetBlockSize()`
1489: @*/
1490: PetscErrorCode VecGetBlockSize(Vec v, PetscInt *bs)
1491: {
1492: PetscFunctionBegin;
1495: PetscCall(PetscLayoutGetBlockSize(v->map, bs));
1496: PetscFunctionReturn(PETSC_SUCCESS);
1497: }
1499: /*@C
1500: VecSetOptionsPrefix - Sets the prefix used for searching for all
1501: `Vec` options in the database.
1503: Logically Collective
1505: Input Parameters:
1506: + v - the `Vec` context
1507: - prefix - the prefix to prepend to all option names
1509: Level: advanced
1511: Note:
1512: A hyphen (-) must NOT be given at the beginning of the prefix name.
1513: The first character of all runtime options is AUTOMATICALLY the hyphen.
1515: .seealso: [](chapter_vectors), `Vec`, `VecSetFromOptions()`
1516: @*/
1517: PetscErrorCode VecSetOptionsPrefix(Vec v, const char prefix[])
1518: {
1519: PetscFunctionBegin;
1521: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)v, prefix));
1522: PetscFunctionReturn(PETSC_SUCCESS);
1523: }
1525: /*@C
1526: VecAppendOptionsPrefix - Appends to the prefix used for searching for all
1527: `Vec` options in the database.
1529: Logically Collective
1531: Input Parameters:
1532: + v - the `Vec` context
1533: - prefix - the prefix to prepend to all option names
1535: Level: advanced
1537: Note:
1538: A hyphen (-) must NOT be given at the beginning of the prefix name.
1539: The first character of all runtime options is AUTOMATICALLY the hyphen.
1541: .seealso: [](chapter_vectors), `Vec`, `VecGetOptionsPrefix()`
1542: @*/
1543: PetscErrorCode VecAppendOptionsPrefix(Vec v, const char prefix[])
1544: {
1545: PetscFunctionBegin;
1547: PetscCall(PetscObjectAppendOptionsPrefix((PetscObject)v, prefix));
1548: PetscFunctionReturn(PETSC_SUCCESS);
1549: }
1551: /*@C
1552: VecGetOptionsPrefix - Sets the prefix used for searching for all
1553: Vec options in the database.
1555: Not Collective
1557: Input Parameter:
1558: . v - the `Vec` context
1560: Output Parameter:
1561: . prefix - pointer to the prefix string used
1563: Level: advanced
1565: Fortran Note:
1566: The user must pass in a string `prefix` of
1567: sufficient length to hold the prefix.
1569: .seealso: [](chapter_vectors), `Vec`, `VecAppendOptionsPrefix()`
1570: @*/
1571: PetscErrorCode VecGetOptionsPrefix(Vec v, const char *prefix[])
1572: {
1573: PetscFunctionBegin;
1575: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)v, prefix));
1576: PetscFunctionReturn(PETSC_SUCCESS);
1577: }
1579: /*@
1580: VecSetUp - Sets up the internal vector data structures for the later use.
1582: Collective
1584: Input Parameters:
1585: . v - the `Vec` context
1587: Level: advanced
1589: Notes:
1590: For basic use of the `Vec` classes the user need not explicitly call
1591: `VecSetUp()`, since these actions will happen automatically.
1593: .seealso: [](chapter_vectors), `Vec`, `VecCreate()`, `VecDestroy()`
1594: @*/
1595: PetscErrorCode VecSetUp(Vec v)
1596: {
1597: PetscMPIInt size;
1599: PetscFunctionBegin;
1601: PetscCheck(v->map->n >= 0 || v->map->N >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Sizes not set");
1602: if (!((PetscObject)v)->type_name) {
1603: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)v), &size));
1604: if (size == 1) {
1605: PetscCall(VecSetType(v, VECSEQ));
1606: } else {
1607: PetscCall(VecSetType(v, VECMPI));
1608: }
1609: }
1610: PetscFunctionReturn(PETSC_SUCCESS);
1611: }
1613: /*
1614: These currently expose the PetscScalar/PetscReal in updating the
1615: cached norm. If we push those down into the implementation these
1616: will become independent of PetscScalar/PetscReal
1617: */
1619: /*@
1620: VecCopy - Copies a vector `y = x`
1622: Logically Collective
1624: Input Parameter:
1625: . x - the vector
1627: Output Parameter:
1628: . y - the copy
1630: Level: beginner
1632: Note:
1633: For default parallel PETSc vectors, both `x` and `y` must be distributed in
1634: the same manner; local copies are done.
1636: Developer Note:
1637: `PetscCheckSameTypeAndComm`(x,1,y,2) is not used on these vectors because we allow one
1638: of the vectors to be sequential and one to be parallel so long as both have the same
1639: local sizes. This is used in some internal functions in PETSc.
1641: .seealso: [](chapter_vectors), `Vec`, `VecDuplicate()`
1642: @*/
1643: PetscErrorCode VecCopy(Vec x, Vec y)
1644: {
1645: PetscBool flgs[4];
1646: PetscReal norms[4] = {0.0, 0.0, 0.0, 0.0};
1648: PetscFunctionBegin;
1653: if (x == y) PetscFunctionReturn(PETSC_SUCCESS);
1654: VecCheckSameLocalSize(x, 1, y, 2);
1655: VecCheckAssembled(x);
1656: PetscCall(VecSetErrorIfLocked(y, 2));
1658: #if !defined(PETSC_USE_MIXED_PRECISION)
1659: for (PetscInt i = 0; i < 4; i++) PetscCall(PetscObjectComposedDataGetReal((PetscObject)x, NormIds[i], norms[i], flgs[i]));
1660: #endif
1662: PetscCall(PetscLogEventBegin(VEC_Copy, x, y, 0, 0));
1663: #if defined(PETSC_USE_MIXED_PRECISION)
1664: extern PetscErrorCode VecGetArray(Vec, double **);
1665: extern PetscErrorCode VecRestoreArray(Vec, double **);
1666: extern PetscErrorCode VecGetArray(Vec, float **);
1667: extern PetscErrorCode VecRestoreArray(Vec, float **);
1668: extern PetscErrorCode VecGetArrayRead(Vec, const double **);
1669: extern PetscErrorCode VecRestoreArrayRead(Vec, const double **);
1670: extern PetscErrorCode VecGetArrayRead(Vec, const float **);
1671: extern PetscErrorCode VecRestoreArrayRead(Vec, const float **);
1672: if ((((PetscObject)x)->precision == PETSC_PRECISION_SINGLE) && (((PetscObject)y)->precision == PETSC_PRECISION_DOUBLE)) {
1673: PetscInt i, n;
1674: const float *xx;
1675: double *yy;
1676: PetscCall(VecGetArrayRead(x, &xx));
1677: PetscCall(VecGetArray(y, &yy));
1678: PetscCall(VecGetLocalSize(x, &n));
1679: for (i = 0; i < n; i++) yy[i] = xx[i];
1680: PetscCall(VecRestoreArrayRead(x, &xx));
1681: PetscCall(VecRestoreArray(y, &yy));
1682: } else if ((((PetscObject)x)->precision == PETSC_PRECISION_DOUBLE) && (((PetscObject)y)->precision == PETSC_PRECISION_SINGLE)) {
1683: PetscInt i, n;
1684: float *yy;
1685: const double *xx;
1686: PetscCall(VecGetArrayRead(x, &xx));
1687: PetscCall(VecGetArray(y, &yy));
1688: PetscCall(VecGetLocalSize(x, &n));
1689: for (i = 0; i < n; i++) yy[i] = (float)xx[i];
1690: PetscCall(VecRestoreArrayRead(x, &xx));
1691: PetscCall(VecRestoreArray(y, &yy));
1692: } else PetscUseTypeMethod(x, copy, y);
1693: #else
1694: PetscUseTypeMethod(x, copy, y);
1695: #endif
1697: PetscCall(PetscObjectStateIncrease((PetscObject)y));
1698: #if !defined(PETSC_USE_MIXED_PRECISION)
1699: for (PetscInt i = 0; i < 4; i++) {
1700: if (flgs[i]) PetscCall(PetscObjectComposedDataSetReal((PetscObject)y, NormIds[i], norms[i]));
1701: }
1702: #endif
1704: PetscCall(PetscLogEventEnd(VEC_Copy, x, y, 0, 0));
1705: PetscFunctionReturn(PETSC_SUCCESS);
1706: }
1708: /*@
1709: VecSwap - Swaps the vectors `x` and `y`.
1711: Logically Collective
1713: Input Parameters:
1714: . x, y - the vectors
1716: Level: advanced
1718: .seealso: [](chapter_vectors), `Vec`, `VecSet()`
1719: @*/
1720: PetscErrorCode VecSwap(Vec x, Vec y)
1721: {
1722: PetscReal normxs[4], normys[4];
1723: PetscBool flgxs[4], flgys[4];
1725: PetscFunctionBegin;
1730: PetscCheckSameTypeAndComm(x, 1, y, 2);
1731: VecCheckSameSize(x, 1, y, 2);
1732: VecCheckAssembled(x);
1733: VecCheckAssembled(y);
1734: PetscCall(VecSetErrorIfLocked(x, 1));
1735: PetscCall(VecSetErrorIfLocked(y, 2));
1737: for (PetscInt i = 0; i < 4; i++) {
1738: PetscCall(PetscObjectComposedDataGetReal((PetscObject)x, NormIds[i], normxs[i], flgxs[i]));
1739: PetscCall(PetscObjectComposedDataGetReal((PetscObject)y, NormIds[i], normys[i], flgys[i]));
1740: }
1742: PetscCall(PetscLogEventBegin(VEC_Swap, x, y, 0, 0));
1743: PetscUseTypeMethod(x, swap, y);
1744: PetscCall(PetscLogEventEnd(VEC_Swap, x, y, 0, 0));
1746: PetscCall(PetscObjectStateIncrease((PetscObject)x));
1747: PetscCall(PetscObjectStateIncrease((PetscObject)y));
1748: for (PetscInt i = 0; i < 4; i++) {
1749: if (flgxs[i]) PetscCall(PetscObjectComposedDataSetReal((PetscObject)y, NormIds[i], normxs[i]));
1750: if (flgys[i]) PetscCall(PetscObjectComposedDataSetReal((PetscObject)x, NormIds[i], normys[i]));
1751: }
1752: PetscFunctionReturn(PETSC_SUCCESS);
1753: }
1755: /*@C
1756: VecStashViewFromOptions - Processes command line options to determine if/how a `VecStash` object is to be viewed.
1758: Collective
1760: Input Parameters:
1761: + obj - the `VecStash` object
1762: . bobj - optional other object that provides the prefix
1763: - optionname - option to activate viewing
1765: Level: intermediate
1767: Developer Note:
1768: This cannot use `PetscObjectViewFromOptions()` because it takes a `Vec` as an argument but does not use VecView
1770: .seealso: [](chapter_vectors), `Vec`, `VecStashSetInitialSize()`
1771: @*/
1772: PetscErrorCode VecStashViewFromOptions(Vec obj, PetscObject bobj, const char optionname[])
1773: {
1774: PetscViewer viewer;
1775: PetscBool flg;
1776: PetscViewerFormat format;
1777: char *prefix;
1779: PetscFunctionBegin;
1780: prefix = bobj ? bobj->prefix : ((PetscObject)obj)->prefix;
1781: PetscCall(PetscOptionsGetViewer(PetscObjectComm((PetscObject)obj), ((PetscObject)obj)->options, prefix, optionname, &viewer, &format, &flg));
1782: if (flg) {
1783: PetscCall(PetscViewerPushFormat(viewer, format));
1784: PetscCall(VecStashView(obj, viewer));
1785: PetscCall(PetscViewerPopFormat(viewer));
1786: PetscCall(PetscViewerDestroy(&viewer));
1787: }
1788: PetscFunctionReturn(PETSC_SUCCESS);
1789: }
1791: /*@
1792: VecStashView - Prints the entries in the vector stash and block stash.
1794: Collective
1796: Input Parameters:
1797: + v - the vector
1798: - viewer - the viewer
1800: Level: advanced
1802: .seealso: [](chapter_vectors), `Vec`, `VecSetBlockSize()`, `VecSetValues()`, `VecSetValuesBlocked()`
1803: @*/
1804: PetscErrorCode VecStashView(Vec v, PetscViewer viewer)
1805: {
1806: PetscMPIInt rank;
1807: PetscInt i, j;
1808: PetscBool match;
1809: VecStash *s;
1810: PetscScalar val;
1812: PetscFunctionBegin;
1815: PetscCheckSameComm(v, 1, viewer, 2);
1817: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &match));
1818: PetscCheck(match, PETSC_COMM_SELF, PETSC_ERR_SUP, "Stash viewer only works with ASCII viewer not %s", ((PetscObject)v)->type_name);
1819: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE));
1820: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)v), &rank));
1821: s = &v->bstash;
1823: /* print block stash */
1824: PetscCall(PetscViewerASCIIPushSynchronized(viewer));
1825: PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "[%d]Vector Block stash size %" PetscInt_FMT " block size %" PetscInt_FMT "\n", rank, s->n, s->bs));
1826: for (i = 0; i < s->n; i++) {
1827: PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "[%d] Element %" PetscInt_FMT " ", rank, s->idx[i]));
1828: for (j = 0; j < s->bs; j++) {
1829: val = s->array[i * s->bs + j];
1830: #if defined(PETSC_USE_COMPLEX)
1831: PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "(%18.16e %18.16e) ", (double)PetscRealPart(val), (double)PetscImaginaryPart(val)));
1832: #else
1833: PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "%18.16e ", (double)val));
1834: #endif
1835: }
1836: PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "\n"));
1837: }
1838: PetscCall(PetscViewerFlush(viewer));
1840: s = &v->stash;
1842: /* print basic stash */
1843: PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "[%d]Vector stash size %" PetscInt_FMT "\n", rank, s->n));
1844: for (i = 0; i < s->n; i++) {
1845: val = s->array[i];
1846: #if defined(PETSC_USE_COMPLEX)
1847: PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "[%d] Element %" PetscInt_FMT " (%18.16e %18.16e) ", rank, s->idx[i], (double)PetscRealPart(val), (double)PetscImaginaryPart(val)));
1848: #else
1849: PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "[%d] Element %" PetscInt_FMT " %18.16e\n", rank, s->idx[i], (double)val));
1850: #endif
1851: }
1852: PetscCall(PetscViewerFlush(viewer));
1853: PetscCall(PetscViewerASCIIPopSynchronized(viewer));
1854: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE));
1855: PetscFunctionReturn(PETSC_SUCCESS);
1856: }
1858: PetscErrorCode PetscOptionsGetVec(PetscOptions options, const char prefix[], const char key[], Vec v, PetscBool *set)
1859: {
1860: PetscInt i, N, rstart, rend;
1861: PetscScalar *xx;
1862: PetscReal *xreal;
1863: PetscBool iset;
1865: PetscFunctionBegin;
1866: PetscCall(VecGetOwnershipRange(v, &rstart, &rend));
1867: PetscCall(VecGetSize(v, &N));
1868: PetscCall(PetscCalloc1(N, &xreal));
1869: PetscCall(PetscOptionsGetRealArray(options, prefix, key, xreal, &N, &iset));
1870: if (iset) {
1871: PetscCall(VecGetArray(v, &xx));
1872: for (i = rstart; i < rend; i++) xx[i - rstart] = xreal[i];
1873: PetscCall(VecRestoreArray(v, &xx));
1874: }
1875: PetscCall(PetscFree(xreal));
1876: if (set) *set = iset;
1877: PetscFunctionReturn(PETSC_SUCCESS);
1878: }
1880: /*@
1881: VecGetLayout - get `PetscLayout` describing vector layout
1883: Not Collective
1885: Input Parameter:
1886: . x - the vector
1888: Output Parameter:
1889: . map - the layout
1891: Level: developer
1893: .seealso: [](chapter_vectors), `Vec`, `VecGetSizes()`, `VecGetOwnershipRange()`, `VecGetOwnershipRanges()`
1894: @*/
1895: PetscErrorCode VecGetLayout(Vec x, PetscLayout *map)
1896: {
1897: PetscFunctionBegin;
1900: *map = x->map;
1901: PetscFunctionReturn(PETSC_SUCCESS);
1902: }
1904: /*@
1905: VecSetLayout - set `PetscLayout` describing vector layout
1907: Not Collective
1909: Input Parameters:
1910: + x - the vector
1911: - map - the layout
1913: Level: developer
1915: Note:
1916: It is normally only valid to replace the layout with a layout known to be equivalent.
1918: .seealso: [](chapter_vectors), `Vec`, `PetscLayout`, `VecGetLayout()`, `VecGetSizes()`, `VecGetOwnershipRange()`, `VecGetOwnershipRanges()`
1919: @*/
1920: PetscErrorCode VecSetLayout(Vec x, PetscLayout map)
1921: {
1922: PetscFunctionBegin;
1924: PetscCall(PetscLayoutReference(map, &x->map));
1925: PetscFunctionReturn(PETSC_SUCCESS);
1926: }
1928: PetscErrorCode VecSetInf(Vec xin)
1929: {
1930: // use of variables one and zero over just doing 1.0/0.0 is deliberate. MSVC complains that
1931: // we are dividing by zero in the latter case (ostensibly because dividing by 0 is UB, but
1932: // only for *integers* not floats).
1933: const PetscScalar one = 1.0, zero = 0.0, inf = one / zero;
1935: PetscFunctionBegin;
1936: if (xin->ops->set) {
1937: PetscUseTypeMethod(xin, set, inf);
1938: } else {
1939: PetscInt n;
1940: PetscScalar *xx;
1942: PetscCall(VecGetLocalSize(xin, &n));
1943: PetscCall(VecGetArrayWrite(xin, &xx));
1944: for (PetscInt i = 0; i < n; ++i) xx[i] = inf;
1945: PetscCall(VecRestoreArrayWrite(xin, &xx));
1946: }
1947: PetscFunctionReturn(PETSC_SUCCESS);
1948: }
1950: /*@
1951: VecBindToCPU - marks a vector to temporarily stay on the CPU and perform computations on the CPU
1953: Logically collective
1955: Input Parameters:
1956: + v - the vector
1957: - flg - bind to the CPU if value of `PETSC_TRUE`
1959: Level: intermediate
1961: .seelaso: `VecBoundToCPU()`
1962: @*/
1963: PetscErrorCode VecBindToCPU(Vec v, PetscBool flg)
1964: {
1965: PetscFunctionBegin;
1968: #if defined(PETSC_HAVE_DEVICE)
1969: if (v->boundtocpu == flg) PetscFunctionReturn(PETSC_SUCCESS);
1970: v->boundtocpu = flg;
1971: PetscTryTypeMethod(v, bindtocpu, flg);
1972: #endif
1973: PetscFunctionReturn(PETSC_SUCCESS);
1974: }
1976: /*@
1977: VecBoundToCPU - query if a vector is bound to the CPU
1979: Not collective
1981: Input Parameter:
1982: . v - the vector
1984: Output Parameter:
1985: . flg - the logical flag
1987: Level: intermediate
1989: .seealso: [](chapter_vectors), `Vec`, `VecBindToCPU()`
1990: @*/
1991: PetscErrorCode VecBoundToCPU(Vec v, PetscBool *flg)
1992: {
1993: PetscFunctionBegin;
1996: #if defined(PETSC_HAVE_DEVICE)
1997: *flg = v->boundtocpu;
1998: #else
1999: *flg = PETSC_TRUE;
2000: #endif
2001: PetscFunctionReturn(PETSC_SUCCESS);
2002: }
2004: /*@
2005: VecSetBindingPropagates - Sets whether the state of being bound to the CPU for a GPU vector type propagates to child and some other associated objects
2007: Input Parameters:
2008: + v - the vector
2009: - flg - flag indicating whether the boundtocpu flag should be propagated
2011: Level: developer
2013: Notes:
2014: If the value of flg is set to true, then `VecDuplicate()` and `VecDuplicateVecs()` will bind created vectors to GPU if the input vector is bound to the CPU.
2015: The created vectors will also have their bindingpropagates flag set to true.
2017: Developer Note:
2018: If a `DMDA` has the `-dm_bind_below option` set to true, then vectors created by `DMCreateGlobalVector()` will have `VecSetBindingPropagates()` called on them to
2019: set their bindingpropagates flag to true.
2021: .seealso: [](chapter_vectors), `Vec`, `MatSetBindingPropagates()`, `VecGetBindingPropagates()`
2022: @*/
2023: PetscErrorCode VecSetBindingPropagates(Vec v, PetscBool flg)
2024: {
2025: PetscFunctionBegin;
2027: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA) || defined(PETSC_HAVE_HIP)
2028: v->bindingpropagates = flg;
2029: #endif
2030: PetscFunctionReturn(PETSC_SUCCESS);
2031: }
2033: /*@
2034: VecGetBindingPropagates - Gets whether the state of being bound to the CPU for a GPU vector type propagates to child and some other associated objects
2036: Input Parameter:
2037: . v - the vector
2039: Output Parameter:
2040: . flg - flag indicating whether the boundtocpu flag will be propagated
2042: Level: developer
2044: .seealso: [](chapter_vectors), `Vec`, `VecSetBindingPropagates()`
2045: @*/
2046: PetscErrorCode VecGetBindingPropagates(Vec v, PetscBool *flg)
2047: {
2048: PetscFunctionBegin;
2051: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA) || defined(PETSC_HAVE_HIP)
2052: *flg = v->bindingpropagates;
2053: #else
2054: *flg = PETSC_FALSE;
2055: #endif
2056: PetscFunctionReturn(PETSC_SUCCESS);
2057: }
2059: /*@C
2060: VecSetPinnedMemoryMin - Set the minimum data size for which pinned memory will be used for host (CPU) allocations.
2062: Logically Collective
2064: Input Parameters:
2065: + v - the vector
2066: - mbytes - minimum data size in bytes
2068: Options Database Key:
2069: . -vec_pinned_memory_min <size> - minimum size (in bytes) for an allocation to use pinned memory on host.
2070: Note that this takes a PetscScalar, to accommodate large values;
2071: specifying -1 ensures that pinned memory will never be used.
2073: Level: developer
2075: .seealso: [](chapter_vectors), `Vec`, `VecGetPinnedMemoryMin()`
2076: @*/
2077: PetscErrorCode VecSetPinnedMemoryMin(Vec v, size_t mbytes)
2078: {
2079: PetscFunctionBegin;
2081: #if PetscDefined(HAVE_DEVICE)
2082: v->minimum_bytes_pinned_memory = mbytes;
2083: #endif
2084: PetscFunctionReturn(PETSC_SUCCESS);
2085: }
2087: /*@C
2088: VecGetPinnedMemoryMin - Get the minimum data size for which pinned memory will be used for host (CPU) allocations.
2090: Logically Collective
2092: Input Parameters:
2093: . v - the vector
2095: Output Parameters:
2096: . mbytes - minimum data size in bytes
2098: Level: developer
2100: .seealso: [](chapter_vectors), `Vec`, `VecSetPinnedMemoryMin()`
2101: @*/
2102: PetscErrorCode VecGetPinnedMemoryMin(Vec v, size_t *mbytes)
2103: {
2104: PetscFunctionBegin;
2107: #if PetscDefined(HAVE_DEVICE)
2108: *mbytes = v->minimum_bytes_pinned_memory;
2109: #endif
2110: PetscFunctionReturn(PETSC_SUCCESS);
2111: }
2113: /*@
2114: VecGetOffloadMask - Get the offload mask of a `Vec`
2116: Not Collective
2118: Input Parameters:
2119: . v - the vector
2121: Output Parameters:
2122: . mask - corresponding `PetscOffloadMask` enum value.
2124: Level: intermediate
2126: .seealso: [](chapter_vectors), `Vec`, `VecCreateSeqCUDA()`, `VecCreateSeqViennaCL()`, `VecGetArray()`, `VecGetType()`
2127: @*/
2128: PetscErrorCode VecGetOffloadMask(Vec v, PetscOffloadMask *mask)
2129: {
2130: PetscFunctionBegin;
2133: *mask = v->offloadmask;
2134: PetscFunctionReturn(PETSC_SUCCESS);
2135: }
2137: #if !defined(PETSC_HAVE_VIENNACL)
2138: PETSC_EXTERN PetscErrorCode VecViennaCLGetCLContext(Vec v, PETSC_UINTPTR_T *ctx)
2139: {
2140: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "PETSc must be configured with --with-opencl to get a Vec's cl_context");
2141: }
2143: PETSC_EXTERN PetscErrorCode VecViennaCLGetCLQueue(Vec v, PETSC_UINTPTR_T *queue)
2144: {
2145: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "PETSc must be configured with --with-opencl to get a Vec's cl_command_queue");
2146: }
2148: PETSC_EXTERN PetscErrorCode VecViennaCLGetCLMem(Vec v, PETSC_UINTPTR_T *queue)
2149: {
2150: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "PETSc must be configured with --with-opencl to get a Vec's cl_mem");
2151: }
2153: PETSC_EXTERN PetscErrorCode VecViennaCLGetCLMemRead(Vec v, PETSC_UINTPTR_T *queue)
2154: {
2155: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "PETSc must be configured with --with-opencl to get a Vec's cl_mem");
2156: }
2158: PETSC_EXTERN PetscErrorCode VecViennaCLGetCLMemWrite(Vec v, PETSC_UINTPTR_T *queue)
2159: {
2160: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "PETSc must be configured with --with-opencl to get a Vec's cl_mem");
2161: }
2163: PETSC_EXTERN PetscErrorCode VecViennaCLRestoreCLMemWrite(Vec v)
2164: {
2165: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "PETSc must be configured with --with-opencl to restore a Vec's cl_mem");
2166: }
2167: #endif