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()`/`VecAssemblyEnd()` 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 in building the stash
32: . bnstash - the size of the block stash
33: - breallocs - the number of additional mallocs incurred in building the block stash (from `VecSetValuesBlocked()`)
35: Level: advanced
37: .seealso: [](ch_vectors), `Vec`, `VecAssemblyBegin()`, `VecAssemblyEnd()`, `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: Notes:
61: All vectors obtained with `VecDuplicate()` from this vector inherit the same mapping.
63: Vectors obtained with `DMCreateGlobaVector()` will often have this attribute attached to the vector so this call is not needed
65: .seealso: [](ch_vectors), `Vec`, `VecAssemblyBegin()`, `VecAssemblyEnd()`, `VecSetValues()`, `VecSetValuesLocal()`,
66: `VecGetLocalToGlobalMapping()`, `VecSetValuesBlockedLocal()`
67: @*/
68: PetscErrorCode VecSetLocalToGlobalMapping(Vec x, ISLocalToGlobalMapping mapping)
69: {
70: PetscFunctionBegin;
73: if (x->ops->setlocaltoglobalmapping) PetscUseTypeMethod(x, setlocaltoglobalmapping, mapping);
74: else PetscCall(PetscLayoutSetISLocalToGlobalMapping(x->map, mapping));
75: PetscFunctionReturn(PETSC_SUCCESS);
76: }
78: /*@
79: VecGetLocalToGlobalMapping - Gets the local-to-global numbering set by `VecSetLocalToGlobalMapping()`
81: Not Collective
83: Input Parameter:
84: . X - the vector
86: Output Parameter:
87: . mapping - the mapping
89: Level: advanced
91: .seealso: [](ch_vectors), `Vec`, `VecSetValuesLocal()`, `VecSetLocalToGlobalMapping()`
92: @*/
93: PetscErrorCode VecGetLocalToGlobalMapping(Vec X, ISLocalToGlobalMapping *mapping)
94: {
95: PetscFunctionBegin;
98: PetscAssertPointer(mapping, 2);
99: *mapping = X->map->mapping;
100: PetscFunctionReturn(PETSC_SUCCESS);
101: }
103: /*@
104: VecAssemblyBegin - Begins assembling the vector; that is ensuring all the vector's entries are stored on the correct MPI process. This routine should
105: be called after completing all calls to `VecSetValues()`.
107: Collective
109: Input Parameter:
110: . vec - the vector
112: Level: beginner
114: .seealso: [](ch_vectors), `Vec`, `VecAssemblyEnd()`, `VecSetValues()`
115: @*/
116: PetscErrorCode VecAssemblyBegin(Vec vec)
117: {
118: PetscFunctionBegin;
121: PetscCall(VecStashViewFromOptions(vec, NULL, "-vec_view_stash"));
122: PetscCall(PetscLogEventBegin(VEC_AssemblyBegin, vec, 0, 0, 0));
123: PetscTryTypeMethod(vec, assemblybegin);
124: PetscCall(PetscLogEventEnd(VEC_AssemblyBegin, vec, 0, 0, 0));
125: PetscCall(PetscObjectStateIncrease((PetscObject)vec));
126: PetscFunctionReturn(PETSC_SUCCESS);
127: }
129: /*@
130: VecAssemblyEnd - Completes assembling the vector. This routine should be called after `VecAssemblyBegin()`.
132: Collective
134: Input Parameter:
135: . vec - the vector
137: Options Database Keys:
138: + -vec_view - Prints vector in `PETSC_VIEWER_DEFAULT` format
139: . -vec_view ::ascii_matlab - Prints vector in `PETSC_VIEWER_ASCII_MATLAB` format to stdout
140: . -vec_view matlab:filename - Prints vector in MATLAB .mat file to filename (requires PETSc configured with --with-matlab)
141: . -vec_view draw - Activates vector viewing using drawing tools
142: . -display <name> - Sets display name (default is host)
143: . -draw_pause <sec> - Sets number of seconds to pause after display
144: - -vec_view socket - Activates vector viewing using a socket
146: Level: beginner
148: .seealso: [](ch_vectors), `Vec`, `VecAssemblyBegin()`, `VecSetValues()`
149: @*/
150: PetscErrorCode VecAssemblyEnd(Vec vec)
151: {
152: PetscFunctionBegin;
154: PetscCall(PetscLogEventBegin(VEC_AssemblyEnd, vec, 0, 0, 0));
156: PetscTryTypeMethod(vec, assemblyend);
157: PetscCall(PetscLogEventEnd(VEC_AssemblyEnd, vec, 0, 0, 0));
158: PetscCall(VecViewFromOptions(vec, NULL, "-vec_view"));
159: PetscFunctionReturn(PETSC_SUCCESS);
160: }
162: /*@
163: VecSetPreallocationCOO - set preallocation for a vector using a coordinate format of the entries with global indices
165: Collective
167: Input Parameters:
168: + x - vector being preallocated
169: . ncoo - number of entries
170: - coo_i - entry indices
172: Level: beginner
174: Notes:
175: This and `VecSetValuesCOO()` provide an alternative API to using `VecSetValues()` to provide vector values.
177: This API is particularly efficient for use on GPUs.
179: Entries can be repeated, see `VecSetValuesCOO()`. Negative indices are not allowed unless vector option `VEC_IGNORE_NEGATIVE_INDICES` is set,
180: in which case they, along with the corresponding entries in `VecSetValuesCOO()`, are ignored. If vector option `VEC_NO_OFF_PROC_ENTRIES` is set,
181: remote entries are ignored, otherwise, they will be properly added or inserted to the vector.
183: The array coo_i[] may be freed immediately after calling this function.
185: .seealso: [](ch_vectors), `Vec`, `VecSetValuesCOO()`, `VecSetPreallocationCOOLocal()`
186: @*/
187: PetscErrorCode VecSetPreallocationCOO(Vec x, PetscCount ncoo, const PetscInt coo_i[])
188: {
189: PetscFunctionBegin;
192: if (ncoo) PetscAssertPointer(coo_i, 3);
193: PetscCall(PetscLogEventBegin(VEC_SetPreallocateCOO, x, 0, 0, 0));
194: PetscCall(PetscLayoutSetUp(x->map));
195: if (x->ops->setpreallocationcoo) {
196: PetscUseTypeMethod(x, setpreallocationcoo, ncoo, coo_i);
197: } else {
198: IS is_coo_i;
199: /* The default implementation only supports ncoo within limit of PetscInt */
200: 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);
201: PetscCall(ISCreateGeneral(PETSC_COMM_SELF, ncoo, coo_i, PETSC_COPY_VALUES, &is_coo_i));
202: PetscCall(PetscObjectCompose((PetscObject)x, "__PETSc_coo_i", (PetscObject)is_coo_i));
203: PetscCall(ISDestroy(&is_coo_i));
204: }
205: PetscCall(PetscLogEventEnd(VEC_SetPreallocateCOO, x, 0, 0, 0));
206: PetscFunctionReturn(PETSC_SUCCESS);
207: }
209: /*@
210: VecSetPreallocationCOOLocal - set preallocation for vectors using a coordinate format of the entries with local indices
212: Collective
214: Input Parameters:
215: + x - vector being preallocated
216: . ncoo - number of entries
217: - coo_i - row indices (local numbering; may be modified)
219: Level: beginner
221: Notes:
222: This and `VecSetValuesCOO()` provide an alternative API to using `VecSetValuesLocal()` to provide vector values.
224: This API is particularly efficient for use on GPUs.
226: The local indices are translated using the local to global mapping, thus `VecSetLocalToGlobalMapping()` must have been
227: called prior to this function.
229: The indices coo_i may be modified within this function. They might be translated to corresponding global
230: indices, but the caller should not rely on them having any specific value after this function returns. The arrays
231: can be freed or reused immediately after this function returns.
233: Entries can be repeated. Negative indices and remote indices might be allowed. see `VecSetPreallocationCOO()`.
235: .seealso: [](ch_vectors), `Vec`, `VecSetPreallocationCOO()`, `VecSetValuesCOO()`
236: @*/
237: PetscErrorCode VecSetPreallocationCOOLocal(Vec x, PetscCount ncoo, PetscInt coo_i[])
238: {
239: ISLocalToGlobalMapping ltog;
241: PetscFunctionBegin;
244: if (ncoo) PetscAssertPointer(coo_i, 3);
245: 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);
246: PetscCall(PetscLayoutSetUp(x->map));
247: PetscCall(VecGetLocalToGlobalMapping(x, <og));
248: if (ltog) PetscCall(ISLocalToGlobalMappingApply(ltog, ncoo, coo_i, coo_i));
249: PetscCall(VecSetPreallocationCOO(x, ncoo, coo_i));
250: PetscFunctionReturn(PETSC_SUCCESS);
251: }
253: /*@
254: VecSetValuesCOO - set values at once in a vector preallocated using `VecSetPreallocationCOO()`
256: Collective
258: Input Parameters:
259: + x - vector being set
260: . coo_v - the value array
261: - imode - the insert mode
263: Level: beginner
265: Note:
266: This and `VecSetPreallocationCOO() or ``VecSetPreallocationCOOLocal()` provide an alternative API to using `VecSetValues()` to provide vector values.
268: This API is particularly efficient for use on GPUs.
270: The values must follow the order of the indices prescribed with `VecSetPreallocationCOO()` or `VecSetPreallocationCOOLocal()`.
271: When repeated entries are specified in the COO indices the `coo_v` values are first properly summed, regardless of the value of `imode`.
272: The imode flag indicates if `coo_v` must be added to the current values of the vector (`ADD_VALUES`) or overwritten (`INSERT_VALUES`).
273: `VecAssemblyBegin()` and `VecAssemblyEnd()` do not need to be called after this routine. It automatically handles the assembly process.
275: .seealso: [](ch_vectors), `Vec`, `VecSetPreallocationCOO()`, `VecSetPreallocationCOOLocal()`, `VecSetValues()`
276: @*/
277: PetscErrorCode VecSetValuesCOO(Vec x, const PetscScalar coo_v[], InsertMode imode)
278: {
279: PetscFunctionBegin;
283: PetscCall(PetscLogEventBegin(VEC_SetValuesCOO, x, 0, 0, 0));
284: if (x->ops->setvaluescoo) {
285: PetscUseTypeMethod(x, setvaluescoo, coo_v, imode);
286: PetscCall(PetscObjectStateIncrease((PetscObject)x));
287: } else {
288: IS is_coo_i;
289: const PetscInt *coo_i;
290: PetscInt ncoo;
291: PetscMemType mtype;
293: PetscCall(PetscGetMemType(coo_v, &mtype));
294: PetscCheck(mtype == PETSC_MEMTYPE_HOST, PetscObjectComm((PetscObject)x), PETSC_ERR_ARG_WRONG, "The basic VecSetValuesCOO() only supports v[] on host");
295: PetscCall(PetscObjectQuery((PetscObject)x, "__PETSc_coo_i", (PetscObject *)&is_coo_i));
296: PetscCheck(is_coo_i, PetscObjectComm((PetscObject)x), PETSC_ERR_COR, "Missing coo_i IS");
297: PetscCall(ISGetLocalSize(is_coo_i, &ncoo));
298: PetscCall(ISGetIndices(is_coo_i, &coo_i));
299: if (imode != ADD_VALUES) PetscCall(VecZeroEntries(x));
300: PetscCall(VecSetValues(x, ncoo, coo_i, coo_v, ADD_VALUES));
301: PetscCall(ISRestoreIndices(is_coo_i, &coo_i));
302: PetscCall(VecAssemblyBegin(x));
303: PetscCall(VecAssemblyEnd(x));
304: }
305: PetscCall(PetscLogEventEnd(VEC_SetValuesCOO, x, 0, 0, 0));
306: PetscFunctionReturn(PETSC_SUCCESS);
307: }
309: static PetscErrorCode VecPointwiseApply_Private(Vec w, Vec x, Vec y, PetscDeviceContext dctx, PetscLogEvent event, const char async_name[], PetscErrorCode (*const pointwise_op)(Vec, Vec, Vec))
310: {
311: PetscErrorCode (*async_fn)(Vec, Vec, Vec, PetscDeviceContext) = NULL;
313: PetscFunctionBegin;
320: PetscCheckSameTypeAndComm(x, 2, y, 3);
321: PetscCheckSameTypeAndComm(y, 3, w, 1);
322: VecCheckSameSize(w, 1, x, 2);
323: VecCheckSameSize(w, 1, y, 3);
324: VecCheckAssembled(x);
325: VecCheckAssembled(y);
326: PetscCall(VecSetErrorIfLocked(w, 1));
329: if (dctx) PetscCall(PetscObjectQueryFunction((PetscObject)w, async_name, &async_fn));
330: if (event) PetscCall(PetscLogEventBegin(event, x, y, w, 0));
331: if (async_fn) {
332: PetscCall((*async_fn)(w, x, y, dctx));
333: } else {
334: PetscCall((*pointwise_op)(w, x, y));
335: }
336: if (event) PetscCall(PetscLogEventEnd(event, x, y, w, 0));
337: PetscCall(PetscObjectStateIncrease((PetscObject)w));
338: PetscFunctionReturn(PETSC_SUCCESS);
339: }
341: PetscErrorCode VecPointwiseMaxAsync_Private(Vec w, Vec x, Vec y, PetscDeviceContext dctx)
342: {
343: PetscFunctionBegin;
344: // REVIEW ME: no log event?
345: PetscCall(VecPointwiseApply_Private(w, x, y, dctx, 0, VecAsyncFnName(PointwiseMax), w->ops->pointwisemax));
346: PetscFunctionReturn(PETSC_SUCCESS);
347: }
349: /*@
350: VecPointwiseMax - Computes the component-wise maximum `w[i] = max(x[i], y[i])`.
352: Logically Collective
354: Input Parameters:
355: + x - the first input vector
356: - y - the second input vector
358: Output Parameter:
359: . w - the result
361: Level: advanced
363: Notes:
364: Any subset of the `x`, `y`, and `w` may be the same vector.
366: For complex numbers compares only the real part
368: .seealso: [](ch_vectors), `Vec`, `VecPointwiseDivide()`, `VecPointwiseMult()`, `VecPointwiseMin()`, `VecPointwiseMaxAbs()`, `VecMaxPointwiseDivide()`
369: @*/
370: PetscErrorCode VecPointwiseMax(Vec w, Vec x, Vec y)
371: {
372: PetscFunctionBegin;
373: PetscCall(VecPointwiseMaxAsync_Private(w, x, y, NULL));
374: PetscFunctionReturn(PETSC_SUCCESS);
375: }
377: PetscErrorCode VecPointwiseMinAsync_Private(Vec w, Vec x, Vec y, PetscDeviceContext dctx)
378: {
379: PetscFunctionBegin;
380: // REVIEW ME: no log event?
381: PetscCall(VecPointwiseApply_Private(w, x, y, dctx, 0, VecAsyncFnName(PointwiseMin), w->ops->pointwisemin));
382: PetscFunctionReturn(PETSC_SUCCESS);
383: }
385: /*@
386: VecPointwiseMin - Computes the component-wise minimum `w[i] = min(x[i], y[i])`.
388: Logically Collective
390: Input Parameters:
391: + x - the first input vector
392: - y - the second input vector
394: Output Parameter:
395: . w - the result
397: Level: advanced
399: Notes:
400: Any subset of the `x`, `y`, and `w` may be the same vector.
402: For complex numbers compares only the real part
404: .seealso: [](ch_vectors), `Vec`, `VecPointwiseDivide()`, `VecPointwiseMult()`, `VecPointwiseMaxAbs()`, `VecMaxPointwiseDivide()`
405: @*/
406: PetscErrorCode VecPointwiseMin(Vec w, Vec x, Vec y)
407: {
408: PetscFunctionBegin;
409: PetscCall(VecPointwiseMinAsync_Private(w, x, y, NULL));
410: PetscFunctionReturn(PETSC_SUCCESS);
411: }
413: PetscErrorCode VecPointwiseMaxAbsAsync_Private(Vec w, Vec x, Vec y, PetscDeviceContext dctx)
414: {
415: PetscFunctionBegin;
416: // REVIEW ME: no log event?
417: PetscCall(VecPointwiseApply_Private(w, x, y, dctx, 0, VecAsyncFnName(PointwiseMaxAbs), w->ops->pointwisemaxabs));
418: PetscFunctionReturn(PETSC_SUCCESS);
419: }
421: /*@
422: VecPointwiseMaxAbs - Computes the component-wise maximum of the absolute values `w[i] = max(abs(x[i]), abs(y[i]))`.
424: Logically Collective
426: Input Parameters:
427: + x - the first input vector
428: - y - the second input vector
430: Output Parameter:
431: . w - the result
433: Level: advanced
435: Notes:
436: Any subset of the `x`, `y`, and `w` may be the same vector.
438: .seealso: [](ch_vectors), `Vec`, `VecPointwiseDivide()`, `VecPointwiseMult()`, `VecPointwiseMin()`, `VecPointwiseMax()`, `VecMaxPointwiseDivide()`
439: @*/
440: PetscErrorCode VecPointwiseMaxAbs(Vec w, Vec x, Vec y)
441: {
442: PetscFunctionBegin;
443: PetscCall(VecPointwiseMaxAbsAsync_Private(w, x, y, NULL));
444: PetscFunctionReturn(PETSC_SUCCESS);
445: }
447: PetscErrorCode VecPointwiseDivideAsync_Private(Vec w, Vec x, Vec y, PetscDeviceContext dctx)
448: {
449: PetscFunctionBegin;
450: // REVIEW ME: no log event?
451: PetscCall(VecPointwiseApply_Private(w, x, y, dctx, 0, VecAsyncFnName(PointwiseDivide), w->ops->pointwisedivide));
452: PetscFunctionReturn(PETSC_SUCCESS);
453: }
455: /*@
456: VecPointwiseDivide - Computes the component-wise division `w[i] = x[i] / y[i]`.
458: Logically Collective
460: Input Parameters:
461: + x - the numerator vector
462: - y - the denominator vector
464: Output Parameter:
465: . w - the result
467: Level: advanced
469: Note:
470: Any subset of the `x`, `y`, and `w` may be the same vector.
472: .seealso: [](ch_vectors), `Vec`, `VecPointwiseMult()`, `VecPointwiseMax()`, `VecPointwiseMin()`, `VecPointwiseMaxAbs()`, `VecMaxPointwiseDivide()`
473: @*/
474: PetscErrorCode VecPointwiseDivide(Vec w, Vec x, Vec y)
475: {
476: PetscFunctionBegin;
477: PetscCall(VecPointwiseDivideAsync_Private(w, x, y, NULL));
478: PetscFunctionReturn(PETSC_SUCCESS);
479: }
481: PetscErrorCode VecPointwiseMultAsync_Private(Vec w, Vec x, Vec y, PetscDeviceContext dctx)
482: {
483: PetscFunctionBegin;
485: PetscCall(VecPointwiseApply_Private(w, x, y, dctx, VEC_PointwiseMult, VecAsyncFnName(PointwiseMult), w->ops->pointwisemult));
486: PetscFunctionReturn(PETSC_SUCCESS);
487: }
489: /*@
490: VecPointwiseMult - Computes the component-wise multiplication `w[i] = x[i] * y[i]`.
492: Logically Collective
494: Input Parameters:
495: + x - the first vector
496: - y - the second vector
498: Output Parameter:
499: . w - the result
501: Level: advanced
503: Note:
504: Any subset of the `x`, `y`, and `w` may be the same vector.
506: .seealso: [](ch_vectors), `Vec`, `VecPointwiseDivide()`, `VecPointwiseMax()`, `VecPointwiseMin()`, `VecPointwiseMaxAbs()`, `VecMaxPointwiseDivide()`
507: @*/
508: PetscErrorCode VecPointwiseMult(Vec w, Vec x, Vec y)
509: {
510: PetscFunctionBegin;
511: PetscCall(VecPointwiseMultAsync_Private(w, x, y, NULL));
512: PetscFunctionReturn(PETSC_SUCCESS);
513: }
515: /*@
516: VecDuplicate - Creates a new vector of the same type as an existing vector.
518: Collective
520: Input Parameter:
521: . v - a vector to mimic
523: Output Parameter:
524: . newv - location to put new vector
526: Level: beginner
528: Notes:
529: `VecDuplicate()` DOES NOT COPY the vector entries, but rather allocates storage
530: for the new vector. Use `VecCopy()` to copy a vector.
532: Use `VecDestroy()` to free the space. Use `VecDuplicateVecs()` to get several
533: vectors.
535: .seealso: [](ch_vectors), `Vec`, `VecDestroy()`, `VecDuplicateVecs()`, `VecCreate()`, `VecCopy()`
536: @*/
537: PetscErrorCode VecDuplicate(Vec v, Vec *newv)
538: {
539: PetscFunctionBegin;
541: PetscAssertPointer(newv, 2);
543: PetscUseTypeMethod(v, duplicate, newv);
544: #if PetscDefined(HAVE_DEVICE)
545: if (v->boundtocpu && v->bindingpropagates) {
546: PetscCall(VecSetBindingPropagates(*newv, PETSC_TRUE));
547: PetscCall(VecBindToCPU(*newv, PETSC_TRUE));
548: }
549: #endif
550: PetscCall(PetscObjectStateIncrease((PetscObject)(*newv)));
551: PetscFunctionReturn(PETSC_SUCCESS);
552: }
554: /*@C
555: VecDestroy - Destroys a vector.
557: Collective
559: Input Parameter:
560: . v - the vector
562: Level: beginner
564: .seealso: [](ch_vectors), `Vec`, `VecCreate()`, `VecDuplicate()`, `VecDestroyVecs()`
565: @*/
566: PetscErrorCode VecDestroy(Vec *v)
567: {
568: PetscFunctionBegin;
569: PetscAssertPointer(v, 1);
570: if (!*v) PetscFunctionReturn(PETSC_SUCCESS);
572: if (--((PetscObject)(*v))->refct > 0) {
573: *v = NULL;
574: PetscFunctionReturn(PETSC_SUCCESS);
575: }
577: PetscCall(PetscObjectSAWsViewOff((PetscObject)*v));
578: /* destroy the internal part */
579: PetscTryTypeMethod(*v, destroy);
580: PetscCall(PetscFree((*v)->defaultrandtype));
581: /* destroy the external/common part */
582: PetscCall(PetscLayoutDestroy(&(*v)->map));
583: PetscCall(PetscHeaderDestroy(v));
584: PetscFunctionReturn(PETSC_SUCCESS);
585: }
587: /*@C
588: VecDuplicateVecs - Creates several vectors of the same type as an existing vector.
590: Collective
592: Input Parameters:
593: + m - the number of vectors to obtain
594: - v - a vector to mimic
596: Output Parameter:
597: . V - location to put pointer to array of vectors
599: Level: intermediate
601: Note:
602: Use `VecDestroyVecs()` to free the space. Use `VecDuplicate()` to form a single
603: vector.
605: Fortran Notes:
606: The Fortran interface is slightly different from that given below, it
607: requires one to pass in `V` a `Vec` array of size at least `m`.
608: See the [](ch_fortran) for details.
610: .seealso: [](ch_vectors), `Vec`, [](ch_fortran), `VecDestroyVecs()`, `VecDuplicate()`, `VecCreate()`, `VecDuplicateVecsF90()`
611: @*/
612: PetscErrorCode VecDuplicateVecs(Vec v, PetscInt m, Vec *V[])
613: {
614: PetscFunctionBegin;
616: PetscAssertPointer(V, 3);
618: PetscUseTypeMethod(v, duplicatevecs, m, V);
619: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA) || defined(PETSC_HAVE_HIP)
620: if (v->boundtocpu && v->bindingpropagates) {
621: PetscInt i;
623: for (i = 0; i < m; i++) {
624: /* Since ops->duplicatevecs might itself propagate the value of boundtocpu,
625: * avoid unnecessary overhead by only calling VecBindToCPU() if the vector isn't already bound. */
626: if (!(*V)[i]->boundtocpu) {
627: PetscCall(VecSetBindingPropagates((*V)[i], PETSC_TRUE));
628: PetscCall(VecBindToCPU((*V)[i], PETSC_TRUE));
629: }
630: }
631: }
632: #endif
633: PetscFunctionReturn(PETSC_SUCCESS);
634: }
636: /*@C
637: VecDestroyVecs - Frees a block of vectors obtained with `VecDuplicateVecs()`.
639: Collective
641: Input Parameters:
642: + m - the number of vectors previously obtained, if zero no vectors are destroyed
643: - vv - pointer to pointer to array of vector pointers, if `NULL` no vectors are destroyed
645: Level: intermediate
647: Fortran Notes:
648: The Fortran interface is slightly different from that given below.
649: See the [](ch_fortran) for details.
651: .seealso: [](ch_vectors), `Vec`, [](ch_fortran), `VecDuplicateVecs()`, `VecDestroyVecsf90()`
652: @*/
653: PetscErrorCode VecDestroyVecs(PetscInt m, Vec *vv[])
654: {
655: PetscFunctionBegin;
656: PetscAssertPointer(vv, 2);
657: PetscCheck(m >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Trying to destroy negative number of vectors %" PetscInt_FMT, m);
658: if (!m || !*vv) {
659: *vv = NULL;
660: PetscFunctionReturn(PETSC_SUCCESS);
661: }
664: PetscCall((*(**vv)->ops->destroyvecs)(m, *vv));
665: *vv = NULL;
666: PetscFunctionReturn(PETSC_SUCCESS);
667: }
669: /*@C
670: VecViewFromOptions - View a vector based on values in the options database
672: Collective
674: Input Parameters:
675: + A - the vector
676: . obj - Optional object that provides the options prefix for this viewing
677: - name - command line option
679: Level: intermediate
681: Note:
682: See `PetscObjectViewFromOptions()` to see the `PetscViewer` and PetscViewerFormat` available
684: .seealso: [](ch_vectors), `Vec`, `VecView`, `PetscObjectViewFromOptions()`, `VecCreate()`
685: @*/
686: PetscErrorCode VecViewFromOptions(Vec A, PetscObject obj, const char name[])
687: {
688: PetscFunctionBegin;
690: PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name));
691: PetscFunctionReturn(PETSC_SUCCESS);
692: }
694: /*@C
695: VecView - Views a vector object.
697: Collective
699: Input Parameters:
700: + vec - the vector
701: - viewer - an optional `PetscViewer` visualization context
703: Level: beginner
705: Notes:
706: The available visualization contexts include
707: + `PETSC_VIEWER_STDOUT_SELF` - for sequential vectors
708: . `PETSC_VIEWER_STDOUT_WORLD` - for parallel vectors created on `PETSC_COMM_WORLD`
709: - `PETSC_VIEWER_STDOUT`_(comm) - for parallel vectors created on MPI communicator comm
711: You can change the format the vector is printed using the
712: option `PetscViewerPushFormat()`.
714: The user can open alternative viewers with
715: + `PetscViewerASCIIOpen()` - Outputs vector to a specified file
716: . `PetscViewerBinaryOpen()` - Outputs vector in binary to a
717: specified file; corresponding input uses `VecLoad()`
718: . `PetscViewerDrawOpen()` - Outputs vector to an X window display
719: . `PetscViewerSocketOpen()` - Outputs vector to Socket viewer
720: - `PetscViewerHDF5Open()` - Outputs vector to HDF5 file viewer
722: The user can call `PetscViewerPushFormat()` to specify the output
723: format of ASCII printed objects (when using `PETSC_VIEWER_STDOUT_SELF`,
724: `PETSC_VIEWER_STDOUT_WORLD` and `PetscViewerASCIIOpen()`). Available formats include
725: + `PETSC_VIEWER_DEFAULT` - default, prints vector contents
726: . `PETSC_VIEWER_ASCII_MATLAB` - prints vector contents in MATLAB format
727: . `PETSC_VIEWER_ASCII_INDEX` - prints vector contents, including indices of vector elements
728: - `PETSC_VIEWER_ASCII_COMMON` - prints vector contents, using a
729: format common among all vector types
731: You can pass any number of vector objects, or other PETSc objects to the same viewer.
733: In the debugger you can do call `VecView`(v,0) to display the vector. (The same holds for any PETSc object viewer).
735: Notes for binary viewer:
736: If you pass multiple vectors to a binary viewer you can read them back in in the same order
737: with `VecLoad()`.
739: If the blocksize of the vector is greater than one then you must provide a unique prefix to
740: the vector with `PetscObjectSetOptionsPrefix`((`PetscObject`)vec,"uniqueprefix"); BEFORE calling `VecView()` on the
741: vector to be stored and then set that same unique prefix on the vector that you pass to `VecLoad()`. The blocksize
742: information is stored in an ASCII file with the same name as the binary file plus a ".info" appended to the
743: filename. If you copy the binary file, make sure you copy the associated .info file with it.
745: See the manual page for `VecLoad()` on the exact format the binary viewer stores
746: the values in the file.
748: Notes for HDF5 Viewer:
749: The name of the `Vec` (given with `PetscObjectSetName()` is the name that is used
750: for the object in the HDF5 file. If you wish to store the same Vec into multiple
751: datasets in the same file (typically with different values), you must change its
752: name each time before calling the `VecView()`. To load the same vector,
753: the name of the Vec object passed to `VecLoad()` must be the same.
755: If the block size of the vector is greater than 1 then it is used as the first dimension in the HDF5 array.
756: If the function `PetscViewerHDF5SetBaseDimension2()`is called then even if the block size is one it will
757: be used as the first dimension in the HDF5 array (that is the HDF5 array will always be two dimensional)
758: See also `PetscViewerHDF5SetTimestep()` which adds an additional complication to reading and writing `Vec`
759: with the HDF5 viewer.
761: .seealso: [](ch_vectors), `Vec`, `VecViewFromOptions()`, `PetscViewerASCIIOpen()`, `PetscViewerDrawOpen()`, `PetscDrawLGCreate()`,
762: `PetscViewerSocketOpen()`, `PetscViewerBinaryOpen()`, `VecLoad()`, `PetscViewerCreate()`,
763: `PetscRealView()`, `PetscScalarView()`, `PetscIntView()`, `PetscViewerHDF5SetTimestep()`
764: @*/
765: PetscErrorCode VecView(Vec vec, PetscViewer viewer)
766: {
767: PetscBool iascii;
768: PetscViewerFormat format;
769: PetscMPIInt size;
771: PetscFunctionBegin;
774: VecCheckAssembled(vec);
775: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)vec), &viewer));
777: PetscCall(PetscViewerGetFormat(viewer, &format));
778: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)vec), &size));
779: if (size == 1 && format == PETSC_VIEWER_LOAD_BALANCE) PetscFunctionReturn(PETSC_SUCCESS);
781: PetscCheck(!vec->stash.n && !vec->bstash.n, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Must call VecAssemblyBegin/End() before viewing this vector");
783: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
784: if (iascii) {
785: PetscInt rows, bs;
787: PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)vec, viewer));
788: if (format == PETSC_VIEWER_ASCII_INFO || format == PETSC_VIEWER_ASCII_INFO_DETAIL) {
789: PetscCall(PetscViewerASCIIPushTab(viewer));
790: PetscCall(VecGetSize(vec, &rows));
791: PetscCall(VecGetBlockSize(vec, &bs));
792: if (bs != 1) {
793: PetscCall(PetscViewerASCIIPrintf(viewer, "length=%" PetscInt_FMT ", bs=%" PetscInt_FMT "\n", rows, bs));
794: } else {
795: PetscCall(PetscViewerASCIIPrintf(viewer, "length=%" PetscInt_FMT "\n", rows));
796: }
797: PetscCall(PetscViewerASCIIPopTab(viewer));
798: }
799: }
800: PetscCall(VecLockReadPush(vec));
801: PetscCall(PetscLogEventBegin(VEC_View, vec, viewer, 0, 0));
802: if ((format == PETSC_VIEWER_NATIVE || format == PETSC_VIEWER_LOAD_BALANCE) && vec->ops->viewnative) {
803: PetscUseTypeMethod(vec, viewnative, viewer);
804: } else {
805: PetscUseTypeMethod(vec, view, viewer);
806: }
807: PetscCall(VecLockReadPop(vec));
808: PetscCall(PetscLogEventEnd(VEC_View, vec, viewer, 0, 0));
809: PetscFunctionReturn(PETSC_SUCCESS);
810: }
812: #if defined(PETSC_USE_DEBUG)
813: #include <../src/sys/totalview/tv_data_display.h>
814: PETSC_UNUSED static int TV_display_type(const struct _p_Vec *v)
815: {
816: const PetscScalar *values;
817: char type[32];
819: TV_add_row("Local rows", "int", &v->map->n);
820: TV_add_row("Global rows", "int", &v->map->N);
821: TV_add_row("Typename", TV_ascii_string_type, ((PetscObject)v)->type_name);
822: PetscCall(VecGetArrayRead((Vec)v, &values));
823: PetscCall(PetscSNPrintf(type, 32, "double[%" PetscInt_FMT "]", v->map->n));
824: TV_add_row("values", type, values);
825: PetscCall(VecRestoreArrayRead((Vec)v, &values));
826: return TV_format_OK;
827: }
828: #endif
830: /*@C
831: VecViewNative - Views a vector object with the original type specific viewer
833: Collective
835: Input Parameters:
836: + vec - the vector
837: - viewer - an optional `PetscViewer` visualization context
839: Level: developer
841: Note:
842: This can be used with, for example, vectors obtained with `DMCreateGlobalVector()` for a `DMDA` to display the vector
843: in the PETSc storage format (each MPI process values follow the previous MPI processes) instead of the "natural" grid
844: ordering.
846: .seealso: [](ch_vectors), `Vec`, `PetscViewerASCIIOpen()`, `PetscViewerDrawOpen()`, `PetscDrawLGCreate()`, `VecView()`
847: `PetscViewerSocketOpen()`, `PetscViewerBinaryOpen()`, `VecLoad()`, `PetscViewerCreate()`,
848: `PetscRealView()`, `PetscScalarView()`, `PetscIntView()`, `PetscViewerHDF5SetTimestep()`
849: @*/
850: PetscErrorCode VecViewNative(Vec vec, PetscViewer viewer)
851: {
852: PetscFunctionBegin;
855: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)vec), &viewer));
857: PetscUseTypeMethod(vec, viewnative, viewer);
858: PetscFunctionReturn(PETSC_SUCCESS);
859: }
861: /*@
862: VecGetSize - Returns the global number of elements of the vector.
864: Not Collective
866: Input Parameter:
867: . x - the vector
869: Output Parameter:
870: . size - the global length of the vector
872: Level: beginner
874: .seealso: [](ch_vectors), `Vec`, `VecGetLocalSize()`
875: @*/
876: PetscErrorCode VecGetSize(Vec x, PetscInt *size)
877: {
878: PetscFunctionBegin;
880: PetscAssertPointer(size, 2);
882: PetscUseTypeMethod(x, getsize, size);
883: PetscFunctionReturn(PETSC_SUCCESS);
884: }
886: /*@
887: VecGetLocalSize - Returns the number of elements of the vector stored
888: in local memory (that is on this MPI process)
890: Not Collective
892: Input Parameter:
893: . x - the vector
895: Output Parameter:
896: . size - the length of the local piece of the vector
898: Level: beginner
900: .seealso: [](ch_vectors), `Vec`, `VecGetSize()`
901: @*/
902: PetscErrorCode VecGetLocalSize(Vec x, PetscInt *size)
903: {
904: PetscFunctionBegin;
906: PetscAssertPointer(size, 2);
908: PetscUseTypeMethod(x, getlocalsize, size);
909: PetscFunctionReturn(PETSC_SUCCESS);
910: }
912: /*@C
913: VecGetOwnershipRange - Returns the range of indices owned by
914: this process. The vector is laid out with the
915: first n1 elements on the first processor, next n2 elements on the
916: second, etc. For certain parallel layouts this range may not be
917: well defined.
919: Not Collective
921: Input Parameter:
922: . x - the vector
924: Output Parameters:
925: + low - the first local element, pass in `NULL` if not interested
926: - high - one more than the last local element, pass in `NULL` if not interested
928: Level: beginner
930: Note:
931: The high argument is one more than the last element stored locally.
933: Fortran Notes:
934: `PETSC_NULL_INTEGER` should be used instead of NULL
936: .seealso: [](ch_vectors), `Vec`, `MatGetOwnershipRange()`, `MatGetOwnershipRanges()`, `VecGetOwnershipRanges()`
937: @*/
938: PetscErrorCode VecGetOwnershipRange(Vec x, PetscInt *low, PetscInt *high)
939: {
940: PetscFunctionBegin;
943: if (low) PetscAssertPointer(low, 2);
944: if (high) PetscAssertPointer(high, 3);
945: if (low) *low = x->map->rstart;
946: if (high) *high = x->map->rend;
947: PetscFunctionReturn(PETSC_SUCCESS);
948: }
950: /*@C
951: VecGetOwnershipRanges - Returns the range of indices owned by EACH processor,
952: The vector is laid out with the
953: first n1 elements on the first processor, next n2 elements on the
954: second, etc. For certain parallel layouts this range may not be
955: well defined.
957: Not Collective
959: Input Parameter:
960: . x - the vector
962: Output Parameter:
963: . ranges - array of length size+1 with the start and end+1 for each process
965: Level: beginner
967: Notes:
968: The high argument is one more than the last element stored locally.
970: If the ranges are used after all vectors that share the ranges has been destroyed then the program will crash accessing ranges[].
972: Fortran Notes:
973: You must PASS in an array of length size+1
975: .seealso: [](ch_vectors), `Vec`, `MatGetOwnershipRange()`, `MatGetOwnershipRanges()`, `VecGetOwnershipRange()`
976: @*/
977: PetscErrorCode VecGetOwnershipRanges(Vec x, const PetscInt *ranges[])
978: {
979: PetscFunctionBegin;
982: PetscCall(PetscLayoutGetRanges(x->map, ranges));
983: PetscFunctionReturn(PETSC_SUCCESS);
984: }
986: // PetscClangLinter pragma disable: -fdoc-section-header-unknown
987: /*@
988: VecSetOption - Sets an option for controlling a vector's behavior.
990: Collective
992: Input Parameters:
993: + x - the vector
994: . op - the option
995: - flag - turn the option on or off
997: Supported Options:
998: + `VEC_IGNORE_OFF_PROC_ENTRIES` - which causes `VecSetValues()` to ignore
999: entries destined to be stored on a separate processor. This can be used
1000: to eliminate the global reduction in the `VecAssemblyBegin()` if you know
1001: that you have only used `VecSetValues()` to set local elements
1002: . `VEC_IGNORE_NEGATIVE_INDICES` - which means you can pass negative indices
1003: in ix in calls to `VecSetValues()` or `VecGetValues()`. These rows are simply
1004: ignored.
1005: - `VEC_SUBSET_OFF_PROC_ENTRIES` - which causes `VecAssemblyBegin()` to assume that the off-process
1006: entries will always be a subset (possibly equal) of the off-process entries set on the
1007: first assembly which had a true `VEC_SUBSET_OFF_PROC_ENTRIES` and the vector has not
1008: changed this flag afterwards. If this assembly is not such first assembly, then this
1009: assembly can reuse the communication pattern setup in that first assembly, thus avoiding
1010: a global reduction. Subsequent assemblies setting off-process values should use the same
1011: InsertMode as the first assembly.
1013: Level: intermediate
1015: Developer Notes:
1016: The `InsertMode` restriction could be removed by packing the stash messages out of place.
1018: .seealso: [](ch_vectors), `Vec`, `VecSetValues()`
1019: @*/
1020: PetscErrorCode VecSetOption(Vec x, VecOption op, PetscBool flag)
1021: {
1022: PetscFunctionBegin;
1025: PetscTryTypeMethod(x, setoption, op, flag);
1026: PetscFunctionReturn(PETSC_SUCCESS);
1027: }
1029: /* Default routines for obtaining and releasing; */
1030: /* may be used by any implementation */
1031: PetscErrorCode VecDuplicateVecs_Default(Vec w, PetscInt m, Vec *V[])
1032: {
1033: PetscFunctionBegin;
1035: PetscAssertPointer(V, 3);
1036: PetscCheck(m > 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "m must be > 0: m = %" PetscInt_FMT, m);
1037: PetscCall(PetscMalloc1(m, V));
1038: for (PetscInt i = 0; i < m; i++) PetscCall(VecDuplicate(w, *V + i));
1039: PetscFunctionReturn(PETSC_SUCCESS);
1040: }
1042: PetscErrorCode VecDestroyVecs_Default(PetscInt m, Vec v[])
1043: {
1044: PetscInt i;
1046: PetscFunctionBegin;
1047: PetscAssertPointer(v, 2);
1048: for (i = 0; i < m; i++) PetscCall(VecDestroy(&v[i]));
1049: PetscCall(PetscFree(v));
1050: PetscFunctionReturn(PETSC_SUCCESS);
1051: }
1053: /*@
1054: VecResetArray - Resets a vector to use its default memory. Call this
1055: after the use of `VecPlaceArray()`.
1057: Not Collective
1059: Input Parameter:
1060: . vec - the vector
1062: Level: developer
1064: .seealso: [](ch_vectors), `Vec`, `VecGetArray()`, `VecRestoreArray()`, `VecReplaceArray()`, `VecPlaceArray()`
1065: @*/
1066: PetscErrorCode VecResetArray(Vec vec)
1067: {
1068: PetscFunctionBegin;
1071: PetscUseTypeMethod(vec, resetarray);
1072: PetscCall(PetscObjectStateIncrease((PetscObject)vec));
1073: PetscFunctionReturn(PETSC_SUCCESS);
1074: }
1076: /*@C
1077: VecLoad - Loads a vector that has been stored in binary or HDF5 format
1078: with `VecView()`.
1080: Collective
1082: Input Parameters:
1083: + vec - the newly loaded vector, this needs to have been created with `VecCreate()` or
1084: some related function before the call to `VecLoad()`.
1085: - viewer - binary file viewer, obtained from `PetscViewerBinaryOpen()` or
1086: HDF5 file viewer, obtained from `PetscViewerHDF5Open()`
1088: Level: intermediate
1090: Notes:
1091: Defaults to the standard `VECSEQ` or `VECMPI`, if you want some other type of `Vec` call `VecSetFromOptions()`
1092: before calling this.
1094: The input file must contain the full global vector, as
1095: written by the routine `VecView()`.
1097: If the type or size of `vec` is not set before a call to `VecLoad()`, PETSc
1098: sets the type and the local and global sizes based on the vector it is reading in. If type and/or
1099: sizes are already set, then the same are used.
1101: If using the binary viewer and the blocksize of the vector is greater than one then you must provide a unique prefix to
1102: the vector with `PetscObjectSetOptionsPrefix`((`PetscObject`)vec,"uniqueprefix"); BEFORE calling `VecView()` on the
1103: vector to be stored and then set that same unique prefix on the vector that you pass to VecLoad(). The blocksize
1104: information is stored in an ASCII file with the same name as the binary file plus a ".info" appended to the
1105: filename. If you copy the binary file, make sure you copy the associated .info file with it.
1107: If using HDF5, you must assign the `Vec` the same name as was used in the Vec
1108: that was stored in the file using `PetscObjectSetName(). Otherwise you will
1109: get the error message: "Cannot H5DOpen2() with `Vec` name NAMEOFOBJECT".
1111: If the HDF5 file contains a two dimensional array the first dimension is treated as the block size
1112: in loading the vector. Hence, for example, using MATLAB notation h5create('vector.dat','/Test_Vec',[27 1]);
1113: will load a vector of size 27 and block size 27 thus resulting in all 27 entries being on the first process of
1114: vectors communicator and the rest of the processes having zero entries
1116: Notes for advanced users when using the binary viewer:
1117: Most users should not need to know the details of the binary storage
1118: format, since `VecLoad()` and `VecView()` completely hide these details.
1119: But for anyone who's interested, the standard binary vector storage
1120: format is
1121: .vb
1122: PetscInt VEC_FILE_CLASSID
1123: PetscInt number of rows
1124: PetscScalar *values of all entries
1125: .ve
1127: In addition, PETSc automatically uses byte swapping to work on all machines; the files
1128: are written ALWAYS using big-endian ordering. On small-endian machines the numbers
1129: are converted to the small-endian format when they are read in from the file.
1130: See PetscBinaryRead() and PetscBinaryWrite() to see how this may be done.
1132: .seealso: [](ch_vectors), `Vec`, `PetscViewerBinaryOpen()`, `VecView()`, `MatLoad()`
1133: @*/
1134: PetscErrorCode VecLoad(Vec vec, PetscViewer viewer)
1135: {
1136: PetscBool isbinary, ishdf5, isadios, isexodusii;
1137: PetscViewerFormat format;
1139: PetscFunctionBegin;
1142: PetscCheckSameComm(vec, 1, viewer, 2);
1143: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary));
1144: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERHDF5, &ishdf5));
1145: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERADIOS, &isadios));
1146: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWEREXODUSII, &isexodusii));
1147: PetscCheck(isbinary || ishdf5 || isadios || isexodusii, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid viewer; open viewer with PetscViewerBinaryOpen()");
1149: PetscCall(VecSetErrorIfLocked(vec, 1));
1150: if (!((PetscObject)vec)->type_name && !vec->ops->create) PetscCall(VecSetType(vec, VECSTANDARD));
1151: PetscCall(PetscLogEventBegin(VEC_Load, viewer, 0, 0, 0));
1152: PetscCall(PetscViewerGetFormat(viewer, &format));
1153: if (format == PETSC_VIEWER_NATIVE && vec->ops->loadnative) {
1154: PetscUseTypeMethod(vec, loadnative, viewer);
1155: } else {
1156: PetscUseTypeMethod(vec, load, viewer);
1157: }
1158: PetscCall(PetscLogEventEnd(VEC_Load, viewer, 0, 0, 0));
1159: PetscFunctionReturn(PETSC_SUCCESS);
1160: }
1162: /*@
1163: VecReciprocal - Replaces each component of a vector by its reciprocal.
1165: Logically Collective
1167: Input Parameter:
1168: . vec - the vector
1170: Output Parameter:
1171: . vec - the vector reciprocal
1173: Level: intermediate
1175: Note:
1176: Vector entries with value 0.0 are not changed
1178: .seealso: [](ch_vectors), `Vec`, `VecLog()`, `VecExp()`, `VecSqrtAbs()`
1179: @*/
1180: PetscErrorCode VecReciprocal(Vec vec)
1181: {
1182: PetscFunctionBegin;
1183: PetscCall(VecReciprocalAsync_Private(vec, NULL));
1184: PetscFunctionReturn(PETSC_SUCCESS);
1185: }
1187: /*@C
1188: VecSetOperation - Allows the user to override a particular vector operation.
1190: Logically Collective; No Fortran Support
1192: Input Parameters:
1193: + vec - The vector to modify
1194: . op - The name of the operation
1195: - f - The function that provides the operation.
1197: Notes:
1198: `f` may be `NULL` to remove the operation from `vec`. Depending on the operation this may be
1199: allowed, however some always expect a valid function. In these cases an error will be raised
1200: when calling the interface routine in question.
1202: See `VecOperation` for an up-to-date list of override-able operations. The operations listed
1203: there have the form `VECOP_<OPERATION>`, where `<OPERATION>` is the suffix (in all capital
1204: letters) of the public interface routine (e.g., `VecView()` -> `VECOP_VIEW`).
1206: Overriding a particular `Vec`'s operation has no affect on any other `Vec`s past, present,
1207: or future. The user should also note that overriding a method is "destructive"; the previous
1208: method is not retained in any way.
1210: Level: advanced
1212: Example Usage:
1213: .vb
1214: // some new VecView() implementation, must have the same signature as the function it seeks
1215: // to replace
1216: PetscErrorCode UserVecView(Vec x, PetscViewer viewer)
1217: {
1218: PetscFunctionBeginUser;
1219: // ...
1220: PetscFunctionReturn(PETSC_SUCCESS);
1221: }
1223: // Create a VECMPI which has a pre-defined VecView() implementation
1224: VecCreateMPI(comm, n, N, &x);
1225: // Calls the VECMPI implementation for VecView()
1226: VecView(x, viewer);
1228: VecSetOperation(x, VECOP_VIEW, (void (*)(void))UserVecView);
1229: // Now calls UserVecView()
1230: VecView(x, viewer);
1231: .ve
1233: .seealso: [](ch_vectors), `Vec`, `VecCreate()`, `MatShellSetOperation()`
1234: @*/
1235: PetscErrorCode VecSetOperation(Vec vec, VecOperation op, void (*f)(void))
1236: {
1237: PetscFunctionBegin;
1239: if (op == VECOP_VIEW && !vec->ops->viewnative) {
1240: vec->ops->viewnative = vec->ops->view;
1241: } else if (op == VECOP_LOAD && !vec->ops->loadnative) {
1242: vec->ops->loadnative = vec->ops->load;
1243: }
1244: ((void (**)(void))vec->ops)[(int)op] = f;
1245: PetscFunctionReturn(PETSC_SUCCESS);
1246: }
1248: /*@
1249: VecStashSetInitialSize - sets the sizes of the vec-stash, that is
1250: used during the assembly process to store values that belong to
1251: other processors.
1253: Not Collective, different processes can have different size stashes
1255: Input Parameters:
1256: + vec - the vector
1257: . size - the initial size of the stash.
1258: - bsize - the initial size of the block-stash(if used).
1260: Options Database Keys:
1261: + -vecstash_initial_size <size> or <size0,size1,...sizep-1> - set initial size
1262: - -vecstash_block_initial_size <bsize> or <bsize0,bsize1,...bsizep-1> - set initial block size
1264: Level: intermediate
1266: Notes:
1267: The block-stash is used for values set with `VecSetValuesBlocked()` while
1268: the stash is used for values set with `VecSetValues()`
1270: Run with the option -info and look for output of the form
1271: VecAssemblyBegin_MPIXXX:Stash has MM entries, uses nn mallocs.
1272: to determine the appropriate value, MM, to use for size and
1273: VecAssemblyBegin_MPIXXX:Block-Stash has BMM entries, uses nn mallocs.
1274: to determine the value, BMM to use for bsize
1276: PETSc attempts to smartly manage the stash size so there is little likelihood setting a
1277: a specific value here will affect performance
1279: .seealso: [](ch_vectors), `Vec`, `VecSetBlockSize()`, `VecSetValues()`, `VecSetValuesBlocked()`, `VecStashView()`
1280: @*/
1281: PetscErrorCode VecStashSetInitialSize(Vec vec, PetscInt size, PetscInt bsize)
1282: {
1283: PetscFunctionBegin;
1285: PetscCall(VecStashSetInitialSize_Private(&vec->stash, size));
1286: PetscCall(VecStashSetInitialSize_Private(&vec->bstash, bsize));
1287: PetscFunctionReturn(PETSC_SUCCESS);
1288: }
1290: /*@
1291: VecSetRandom - Sets all components of a vector to random numbers.
1293: Logically Collective
1295: Input Parameters:
1296: + x - the vector
1297: - rctx - the random number context, formed by `PetscRandomCreate()`, or use `NULL` and it will create one internally.
1299: Output Parameter:
1300: . x - the vector
1302: Example of Usage:
1303: .vb
1304: PetscRandomCreate(PETSC_COMM_WORLD,&rctx);
1305: VecSetRandom(x,rctx);
1306: PetscRandomDestroy(&rctx);
1307: .ve
1309: Level: intermediate
1311: .seealso: [](ch_vectors), `Vec`, `VecSet()`, `VecSetValues()`, `PetscRandomCreate()`, `PetscRandomDestroy()`
1312: @*/
1313: PetscErrorCode VecSetRandom(Vec x, PetscRandom rctx)
1314: {
1315: PetscRandom randObj = NULL;
1317: PetscFunctionBegin;
1321: VecCheckAssembled(x);
1322: PetscCall(VecSetErrorIfLocked(x, 1));
1324: if (!rctx) {
1325: PetscCall(PetscRandomCreate(PetscObjectComm((PetscObject)x), &randObj));
1326: PetscCall(PetscRandomSetType(randObj, x->defaultrandtype));
1327: PetscCall(PetscRandomSetFromOptions(randObj));
1328: rctx = randObj;
1329: }
1331: PetscCall(PetscLogEventBegin(VEC_SetRandom, x, rctx, 0, 0));
1332: PetscUseTypeMethod(x, setrandom, rctx);
1333: PetscCall(PetscLogEventEnd(VEC_SetRandom, x, rctx, 0, 0));
1335: PetscCall(PetscRandomDestroy(&randObj));
1336: PetscCall(PetscObjectStateIncrease((PetscObject)x));
1337: PetscFunctionReturn(PETSC_SUCCESS);
1338: }
1340: /*@
1341: VecZeroEntries - puts a `0.0` in each element of a vector
1343: Logically Collective
1345: Input Parameter:
1346: . vec - The vector
1348: Level: beginner
1350: Note:
1351: If the norm of the vector is known to be zero then this skips the unneeded zeroing process
1353: .seealso: [](ch_vectors), `Vec`, `VecCreate()`, `VecSetOptionsPrefix()`, `VecSet()`, `VecSetValues()`
1354: @*/
1355: PetscErrorCode VecZeroEntries(Vec vec)
1356: {
1357: PetscFunctionBegin;
1358: PetscCall(VecSet(vec, 0));
1359: PetscFunctionReturn(PETSC_SUCCESS);
1360: }
1362: /*
1363: VecSetTypeFromOptions_Private - Sets the type of vector from user options. Defaults to a PETSc sequential vector on one
1364: processor and a PETSc MPI vector on more than one processor.
1366: Collective
1368: Input Parameter:
1369: . vec - The vector
1371: Level: intermediate
1373: .seealso: [](ch_vectors), `Vec`, `VecSetFromOptions()`, `VecSetType()`
1374: */
1375: static PetscErrorCode VecSetTypeFromOptions_Private(Vec vec, PetscOptionItems *PetscOptionsObject)
1376: {
1377: PetscBool opt;
1378: VecType defaultType;
1379: char typeName[256];
1380: PetscMPIInt size;
1382: PetscFunctionBegin;
1383: if (((PetscObject)vec)->type_name) defaultType = ((PetscObject)vec)->type_name;
1384: else {
1385: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)vec), &size));
1386: if (size > 1) defaultType = VECMPI;
1387: else defaultType = VECSEQ;
1388: }
1390: PetscCall(VecRegisterAll());
1391: PetscCall(PetscOptionsFList("-vec_type", "Vector type", "VecSetType", VecList, defaultType, typeName, 256, &opt));
1392: if (opt) {
1393: PetscCall(VecSetType(vec, typeName));
1394: } else {
1395: PetscCall(VecSetType(vec, defaultType));
1396: }
1397: PetscFunctionReturn(PETSC_SUCCESS);
1398: }
1400: /*@
1401: VecSetFromOptions - Configures the vector from the options database.
1403: Collective
1405: Input Parameter:
1406: . vec - The vector
1408: Level: beginner
1410: Notes:
1411: To see all options, run your program with the -help option.
1413: Must be called after `VecCreate()` but before the vector is used.
1415: .seealso: [](ch_vectors), `Vec`, `VecCreate()`, `VecSetOptionsPrefix()`
1416: @*/
1417: PetscErrorCode VecSetFromOptions(Vec vec)
1418: {
1419: PetscBool flg;
1420: PetscInt bind_below = 0;
1422: PetscFunctionBegin;
1425: PetscObjectOptionsBegin((PetscObject)vec);
1426: /* Handle vector type options */
1427: PetscCall(VecSetTypeFromOptions_Private(vec, PetscOptionsObject));
1429: /* Handle specific vector options */
1430: PetscTryTypeMethod(vec, setfromoptions, PetscOptionsObject);
1432: /* Bind to CPU if below a user-specified size threshold.
1433: * This perhaps belongs in the options for the GPU Vec types, but VecBindToCPU() does nothing when called on non-GPU types,
1434: * and putting it here makes is more maintainable than duplicating this for all. */
1435: 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));
1436: if (flg && vec->map->n < bind_below) PetscCall(VecBindToCPU(vec, PETSC_TRUE));
1438: /* process any options handlers added with PetscObjectAddOptionsHandler() */
1439: PetscCall(PetscObjectProcessOptionsHandlers((PetscObject)vec, PetscOptionsObject));
1440: PetscOptionsEnd();
1441: PetscFunctionReturn(PETSC_SUCCESS);
1442: }
1444: /*@
1445: VecSetSizes - Sets the local and global sizes, and checks to determine compatibility of the sizes
1447: Collective
1449: Input Parameters:
1450: + v - the vector
1451: . n - the local size (or `PETSC_DECIDE` to have it set)
1452: - N - the global size (or `PETSC_DETERMINE` to have it set)
1454: Level: intermediate
1456: Notes:
1457: `N` cannot be `PETSC_DETERMINE` if `n` is `PETSC_DECIDE`
1459: If one processor calls this with `N` of `PETSC_DETERMINE` then all processors must, otherwise the program will hang.
1461: .seealso: [](ch_vectors), `Vec`, `VecGetSize()`, `PetscSplitOwnership()`
1462: @*/
1463: PetscErrorCode VecSetSizes(Vec v, PetscInt n, PetscInt N)
1464: {
1465: PetscFunctionBegin;
1467: if (N >= 0) {
1469: PetscCheck(n <= N, PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Local size %" PetscInt_FMT " cannot be larger than global size %" PetscInt_FMT, n, N);
1470: }
1471: 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,
1472: v->map->n, v->map->N);
1473: v->map->n = n;
1474: v->map->N = N;
1475: PetscTryTypeMethod(v, create);
1476: v->ops->create = NULL;
1477: PetscFunctionReturn(PETSC_SUCCESS);
1478: }
1480: /*@
1481: VecSetBlockSize - Sets the block size for future calls to `VecSetValuesBlocked()`
1482: and `VecSetValuesBlockedLocal()`.
1484: Logically Collective
1486: Input Parameters:
1487: + v - the vector
1488: - bs - the blocksize
1490: Level: advanced
1492: Note:
1493: All vectors obtained by `VecDuplicate()` inherit the same blocksize.
1495: Vectors obtained with `DMCreateGlobalVector()` and `DMCreateLocalVector()` generally already have a blocksize set based on the state of the `DM`
1497: .seealso: [](ch_vectors), `Vec`, `VecSetValuesBlocked()`, `VecSetLocalToGlobalMapping()`, `VecGetBlockSize()`
1498: @*/
1499: PetscErrorCode VecSetBlockSize(Vec v, PetscInt bs)
1500: {
1501: PetscFunctionBegin;
1504: PetscCall(PetscLayoutSetBlockSize(v->map, bs));
1505: v->bstash.bs = bs; /* use the same blocksize for the vec's block-stash */
1506: PetscFunctionReturn(PETSC_SUCCESS);
1507: }
1509: /*@
1510: VecGetBlockSize - Gets the blocksize for the vector, i.e. what is used for `VecSetValuesBlocked()`
1511: and `VecSetValuesBlockedLocal()`.
1513: Not Collective
1515: Input Parameter:
1516: . v - the vector
1518: Output Parameter:
1519: . bs - the blocksize
1521: Level: advanced
1523: Note:
1524: All vectors obtained by `VecDuplicate()` inherit the same blocksize.
1526: .seealso: [](ch_vectors), `Vec`, `VecSetValuesBlocked()`, `VecSetLocalToGlobalMapping()`, `VecSetBlockSize()`
1527: @*/
1528: PetscErrorCode VecGetBlockSize(Vec v, PetscInt *bs)
1529: {
1530: PetscFunctionBegin;
1532: PetscAssertPointer(bs, 2);
1533: PetscCall(PetscLayoutGetBlockSize(v->map, bs));
1534: PetscFunctionReturn(PETSC_SUCCESS);
1535: }
1537: /*@C
1538: VecSetOptionsPrefix - Sets the prefix used for searching for all
1539: `Vec` options in the database.
1541: Logically Collective
1543: Input Parameters:
1544: + v - the `Vec` context
1545: - prefix - the prefix to prepend to all option names
1547: Level: advanced
1549: Note:
1550: A hyphen (-) must NOT be given at the beginning of the prefix name.
1551: The first character of all runtime options is AUTOMATICALLY the hyphen.
1553: .seealso: [](ch_vectors), `Vec`, `VecSetFromOptions()`
1554: @*/
1555: PetscErrorCode VecSetOptionsPrefix(Vec v, const char prefix[])
1556: {
1557: PetscFunctionBegin;
1559: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)v, prefix));
1560: PetscFunctionReturn(PETSC_SUCCESS);
1561: }
1563: /*@C
1564: VecAppendOptionsPrefix - Appends to the prefix used for searching for all
1565: `Vec` options in the database.
1567: Logically Collective
1569: Input Parameters:
1570: + v - the `Vec` context
1571: - prefix - the prefix to prepend to all option names
1573: Level: advanced
1575: Note:
1576: A hyphen (-) must NOT be given at the beginning of the prefix name.
1577: The first character of all runtime options is AUTOMATICALLY the hyphen.
1579: .seealso: [](ch_vectors), `Vec`, `VecGetOptionsPrefix()`
1580: @*/
1581: PetscErrorCode VecAppendOptionsPrefix(Vec v, const char prefix[])
1582: {
1583: PetscFunctionBegin;
1585: PetscCall(PetscObjectAppendOptionsPrefix((PetscObject)v, prefix));
1586: PetscFunctionReturn(PETSC_SUCCESS);
1587: }
1589: /*@C
1590: VecGetOptionsPrefix - Sets the prefix used for searching for all
1591: Vec options in the database.
1593: Not Collective
1595: Input Parameter:
1596: . v - the `Vec` context
1598: Output Parameter:
1599: . prefix - pointer to the prefix string used
1601: Level: advanced
1603: Fortran Notes:
1604: The user must pass in a string `prefix` of
1605: sufficient length to hold the prefix.
1607: .seealso: [](ch_vectors), `Vec`, `VecAppendOptionsPrefix()`
1608: @*/
1609: PetscErrorCode VecGetOptionsPrefix(Vec v, const char *prefix[])
1610: {
1611: PetscFunctionBegin;
1613: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)v, prefix));
1614: PetscFunctionReturn(PETSC_SUCCESS);
1615: }
1617: /*@
1618: VecSetUp - Sets up the internal vector data structures for the later use.
1620: Collective
1622: Input Parameter:
1623: . v - the `Vec` context
1625: Level: advanced
1627: Notes:
1628: For basic use of the `Vec` classes the user need not explicitly call
1629: `VecSetUp()`, since these actions will happen automatically.
1631: .seealso: [](ch_vectors), `Vec`, `VecCreate()`, `VecDestroy()`
1632: @*/
1633: PetscErrorCode VecSetUp(Vec v)
1634: {
1635: PetscMPIInt size;
1637: PetscFunctionBegin;
1639: PetscCheck(v->map->n >= 0 || v->map->N >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Sizes not set");
1640: if (!((PetscObject)v)->type_name) {
1641: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)v), &size));
1642: if (size == 1) {
1643: PetscCall(VecSetType(v, VECSEQ));
1644: } else {
1645: PetscCall(VecSetType(v, VECMPI));
1646: }
1647: }
1648: PetscFunctionReturn(PETSC_SUCCESS);
1649: }
1651: /*
1652: These currently expose the PetscScalar/PetscReal in updating the
1653: cached norm. If we push those down into the implementation these
1654: will become independent of PetscScalar/PetscReal
1655: */
1657: PetscErrorCode VecCopyAsync_Private(Vec x, Vec y, PetscDeviceContext dctx)
1658: {
1659: PetscBool flgs[4];
1660: PetscReal norms[4] = {0.0, 0.0, 0.0, 0.0};
1662: PetscFunctionBegin;
1667: if (x == y) PetscFunctionReturn(PETSC_SUCCESS);
1668: VecCheckSameLocalSize(x, 1, y, 2);
1669: VecCheckAssembled(x);
1670: PetscCall(VecSetErrorIfLocked(y, 2));
1672: #if !defined(PETSC_USE_MIXED_PRECISION)
1673: for (PetscInt i = 0; i < 4; i++) PetscCall(PetscObjectComposedDataGetReal((PetscObject)x, NormIds[i], norms[i], flgs[i]));
1674: #endif
1676: PetscCall(PetscLogEventBegin(VEC_Copy, x, y, 0, 0));
1677: #if defined(PETSC_USE_MIXED_PRECISION)
1678: extern PetscErrorCode VecGetArray(Vec, double **);
1679: extern PetscErrorCode VecRestoreArray(Vec, double **);
1680: extern PetscErrorCode VecGetArray(Vec, float **);
1681: extern PetscErrorCode VecRestoreArray(Vec, float **);
1682: extern PetscErrorCode VecGetArrayRead(Vec, const double **);
1683: extern PetscErrorCode VecRestoreArrayRead(Vec, const double **);
1684: extern PetscErrorCode VecGetArrayRead(Vec, const float **);
1685: extern PetscErrorCode VecRestoreArrayRead(Vec, const float **);
1686: if ((((PetscObject)x)->precision == PETSC_PRECISION_SINGLE) && (((PetscObject)y)->precision == PETSC_PRECISION_DOUBLE)) {
1687: PetscInt i, n;
1688: const float *xx;
1689: double *yy;
1690: PetscCall(VecGetArrayRead(x, &xx));
1691: PetscCall(VecGetArray(y, &yy));
1692: PetscCall(VecGetLocalSize(x, &n));
1693: for (i = 0; i < n; i++) yy[i] = xx[i];
1694: PetscCall(VecRestoreArrayRead(x, &xx));
1695: PetscCall(VecRestoreArray(y, &yy));
1696: } else if ((((PetscObject)x)->precision == PETSC_PRECISION_DOUBLE) && (((PetscObject)y)->precision == PETSC_PRECISION_SINGLE)) {
1697: PetscInt i, n;
1698: float *yy;
1699: const double *xx;
1700: PetscCall(VecGetArrayRead(x, &xx));
1701: PetscCall(VecGetArray(y, &yy));
1702: PetscCall(VecGetLocalSize(x, &n));
1703: for (i = 0; i < n; i++) yy[i] = (float)xx[i];
1704: PetscCall(VecRestoreArrayRead(x, &xx));
1705: PetscCall(VecRestoreArray(y, &yy));
1706: } else PetscUseTypeMethod(x, copy, y);
1707: #else
1708: VecMethodDispatch(x, dctx, VecAsyncFnName(Copy), copy, (Vec, Vec, PetscDeviceContext), y);
1709: #endif
1711: PetscCall(PetscObjectStateIncrease((PetscObject)y));
1712: #if !defined(PETSC_USE_MIXED_PRECISION)
1713: for (PetscInt i = 0; i < 4; i++) {
1714: if (flgs[i]) PetscCall(PetscObjectComposedDataSetReal((PetscObject)y, NormIds[i], norms[i]));
1715: }
1716: #endif
1718: PetscCall(PetscLogEventEnd(VEC_Copy, x, y, 0, 0));
1719: PetscFunctionReturn(PETSC_SUCCESS);
1720: }
1722: /*@
1723: VecCopy - Copies a vector `y = x`
1725: Logically Collective
1727: Input Parameter:
1728: . x - the vector
1730: Output Parameter:
1731: . y - the copy
1733: Level: beginner
1735: Note:
1736: For default parallel PETSc vectors, both `x` and `y` must be distributed in
1737: the same manner; local copies are done.
1739: Developer Notes:
1740: `PetscCheckSameTypeAndComm`(x,1,y,2) is not used on these vectors because we allow one
1741: of the vectors to be sequential and one to be parallel so long as both have the same
1742: local sizes. This is used in some internal functions in PETSc.
1744: .seealso: [](ch_vectors), `Vec`, `VecDuplicate()`
1745: @*/
1746: PetscErrorCode VecCopy(Vec x, Vec y)
1747: {
1748: PetscFunctionBegin;
1749: PetscCall(VecCopyAsync_Private(x, y, NULL));
1750: PetscFunctionReturn(PETSC_SUCCESS);
1751: }
1753: PetscErrorCode VecSwapAsync_Private(Vec x, Vec y, PetscDeviceContext dctx)
1754: {
1755: PetscReal normxs[4], normys[4];
1756: PetscBool flgxs[4], flgys[4];
1758: PetscFunctionBegin;
1763: PetscCheckSameTypeAndComm(x, 1, y, 2);
1764: VecCheckSameSize(x, 1, y, 2);
1765: VecCheckAssembled(x);
1766: VecCheckAssembled(y);
1767: PetscCall(VecSetErrorIfLocked(x, 1));
1768: PetscCall(VecSetErrorIfLocked(y, 2));
1770: for (PetscInt i = 0; i < 4; i++) {
1771: PetscCall(PetscObjectComposedDataGetReal((PetscObject)x, NormIds[i], normxs[i], flgxs[i]));
1772: PetscCall(PetscObjectComposedDataGetReal((PetscObject)y, NormIds[i], normys[i], flgys[i]));
1773: }
1775: PetscCall(PetscLogEventBegin(VEC_Swap, x, y, 0, 0));
1776: VecMethodDispatch(x, dctx, VecAsyncFnName(Swap), swap, (Vec, Vec, PetscDeviceContext), y);
1777: PetscCall(PetscLogEventEnd(VEC_Swap, x, y, 0, 0));
1779: PetscCall(PetscObjectStateIncrease((PetscObject)x));
1780: PetscCall(PetscObjectStateIncrease((PetscObject)y));
1781: for (PetscInt i = 0; i < 4; i++) {
1782: if (flgxs[i]) PetscCall(PetscObjectComposedDataSetReal((PetscObject)y, NormIds[i], normxs[i]));
1783: if (flgys[i]) PetscCall(PetscObjectComposedDataSetReal((PetscObject)x, NormIds[i], normys[i]));
1784: }
1785: PetscFunctionReturn(PETSC_SUCCESS);
1786: }
1787: /*@
1788: VecSwap - Swaps the values between two vectors, `x` and `y`.
1790: Logically Collective
1792: Input Parameters:
1793: + x - the first vector
1794: - y - the second vector
1796: Level: advanced
1798: .seealso: [](ch_vectors), `Vec`, `VecSet()`
1799: @*/
1800: PetscErrorCode VecSwap(Vec x, Vec y)
1801: {
1802: PetscFunctionBegin;
1803: PetscCall(VecSwapAsync_Private(x, y, NULL));
1804: PetscFunctionReturn(PETSC_SUCCESS);
1805: }
1807: /*@C
1808: VecStashViewFromOptions - Processes command line options to determine if/how a `VecStash` object is to be viewed.
1810: Collective
1812: Input Parameters:
1813: + obj - the `Vec` containing a stash
1814: . bobj - optional other object that provides the prefix
1815: - optionname - option to activate viewing
1817: Level: intermediate
1819: Developer Notes:
1820: This cannot use `PetscObjectViewFromOptions()` because it takes a `Vec` as an argument but does not use `VecView()`
1822: .seealso: [](ch_vectors), `Vec`, `VecStashSetInitialSize()`
1823: @*/
1824: PetscErrorCode VecStashViewFromOptions(Vec obj, PetscObject bobj, const char optionname[])
1825: {
1826: PetscViewer viewer;
1827: PetscBool flg;
1828: PetscViewerFormat format;
1829: char *prefix;
1831: PetscFunctionBegin;
1832: prefix = bobj ? bobj->prefix : ((PetscObject)obj)->prefix;
1833: PetscCall(PetscOptionsGetViewer(PetscObjectComm((PetscObject)obj), ((PetscObject)obj)->options, prefix, optionname, &viewer, &format, &flg));
1834: if (flg) {
1835: PetscCall(PetscViewerPushFormat(viewer, format));
1836: PetscCall(VecStashView(obj, viewer));
1837: PetscCall(PetscViewerPopFormat(viewer));
1838: PetscCall(PetscViewerDestroy(&viewer));
1839: }
1840: PetscFunctionReturn(PETSC_SUCCESS);
1841: }
1843: /*@
1844: VecStashView - Prints the entries in the vector stash and block stash.
1846: Collective
1848: Input Parameters:
1849: + v - the vector
1850: - viewer - the viewer
1852: Level: advanced
1854: .seealso: [](ch_vectors), `Vec`, `VecSetBlockSize()`, `VecSetValues()`, `VecSetValuesBlocked()`
1855: @*/
1856: PetscErrorCode VecStashView(Vec v, PetscViewer viewer)
1857: {
1858: PetscMPIInt rank;
1859: PetscInt i, j;
1860: PetscBool match;
1861: VecStash *s;
1862: PetscScalar val;
1864: PetscFunctionBegin;
1867: PetscCheckSameComm(v, 1, viewer, 2);
1869: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &match));
1870: PetscCheck(match, PETSC_COMM_SELF, PETSC_ERR_SUP, "Stash viewer only works with ASCII viewer not %s", ((PetscObject)v)->type_name);
1871: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE));
1872: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)v), &rank));
1873: s = &v->bstash;
1875: /* print block stash */
1876: PetscCall(PetscViewerASCIIPushSynchronized(viewer));
1877: PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "[%d]Vector Block stash size %" PetscInt_FMT " block size %" PetscInt_FMT "\n", rank, s->n, s->bs));
1878: for (i = 0; i < s->n; i++) {
1879: PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "[%d] Element %" PetscInt_FMT " ", rank, s->idx[i]));
1880: for (j = 0; j < s->bs; j++) {
1881: val = s->array[i * s->bs + j];
1882: #if defined(PETSC_USE_COMPLEX)
1883: PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "(%18.16e %18.16e) ", (double)PetscRealPart(val), (double)PetscImaginaryPart(val)));
1884: #else
1885: PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "%18.16e ", (double)val));
1886: #endif
1887: }
1888: PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "\n"));
1889: }
1890: PetscCall(PetscViewerFlush(viewer));
1892: s = &v->stash;
1894: /* print basic stash */
1895: PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "[%d]Vector stash size %" PetscInt_FMT "\n", rank, s->n));
1896: for (i = 0; i < s->n; i++) {
1897: val = s->array[i];
1898: #if defined(PETSC_USE_COMPLEX)
1899: PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "[%d] Element %" PetscInt_FMT " (%18.16e %18.16e) ", rank, s->idx[i], (double)PetscRealPart(val), (double)PetscImaginaryPart(val)));
1900: #else
1901: PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "[%d] Element %" PetscInt_FMT " %18.16e\n", rank, s->idx[i], (double)val));
1902: #endif
1903: }
1904: PetscCall(PetscViewerFlush(viewer));
1905: PetscCall(PetscViewerASCIIPopSynchronized(viewer));
1906: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE));
1907: PetscFunctionReturn(PETSC_SUCCESS);
1908: }
1910: PetscErrorCode PetscOptionsGetVec(PetscOptions options, const char prefix[], const char key[], Vec v, PetscBool *set)
1911: {
1912: PetscInt i, N, rstart, rend;
1913: PetscScalar *xx;
1914: PetscReal *xreal;
1915: PetscBool iset;
1917: PetscFunctionBegin;
1918: PetscCall(VecGetOwnershipRange(v, &rstart, &rend));
1919: PetscCall(VecGetSize(v, &N));
1920: PetscCall(PetscCalloc1(N, &xreal));
1921: PetscCall(PetscOptionsGetRealArray(options, prefix, key, xreal, &N, &iset));
1922: if (iset) {
1923: PetscCall(VecGetArray(v, &xx));
1924: for (i = rstart; i < rend; i++) xx[i - rstart] = xreal[i];
1925: PetscCall(VecRestoreArray(v, &xx));
1926: }
1927: PetscCall(PetscFree(xreal));
1928: if (set) *set = iset;
1929: PetscFunctionReturn(PETSC_SUCCESS);
1930: }
1932: /*@
1933: VecGetLayout - get `PetscLayout` describing a vector layout
1935: Not Collective
1937: Input Parameter:
1938: . x - the vector
1940: Output Parameter:
1941: . map - the layout
1943: Level: developer
1945: Note:
1946: The layout determines what vector elements are contained on each MPI process
1948: .seealso: [](ch_vectors), `PetscLayout`, `Vec`, `VecGetSizes()`, `VecGetOwnershipRange()`, `VecGetOwnershipRanges()`
1949: @*/
1950: PetscErrorCode VecGetLayout(Vec x, PetscLayout *map)
1951: {
1952: PetscFunctionBegin;
1954: PetscAssertPointer(map, 2);
1955: *map = x->map;
1956: PetscFunctionReturn(PETSC_SUCCESS);
1957: }
1959: /*@
1960: VecSetLayout - set `PetscLayout` describing vector layout
1962: Not Collective
1964: Input Parameters:
1965: + x - the vector
1966: - map - the layout
1968: Level: developer
1970: Note:
1971: It is normally only valid to replace the layout with a layout known to be equivalent.
1973: .seealso: [](ch_vectors), `Vec`, `PetscLayout`, `VecGetLayout()`, `VecGetSizes()`, `VecGetOwnershipRange()`, `VecGetOwnershipRanges()`
1974: @*/
1975: PetscErrorCode VecSetLayout(Vec x, PetscLayout map)
1976: {
1977: PetscFunctionBegin;
1979: PetscCall(PetscLayoutReference(map, &x->map));
1980: PetscFunctionReturn(PETSC_SUCCESS);
1981: }
1983: PetscErrorCode VecSetInf(Vec xin)
1984: {
1985: // use of variables one and zero over just doing 1.0/0.0 is deliberate. MSVC complains that
1986: // we are dividing by zero in the latter case (ostensibly because dividing by 0 is UB, but
1987: // only for *integers* not floats).
1988: const PetscScalar one = 1.0, zero = 0.0;
1989: PetscScalar inf;
1991: PetscFunctionBegin;
1992: PetscCall(PetscFPTrapPush(PETSC_FP_TRAP_OFF));
1993: inf = one / zero;
1994: PetscCall(PetscFPTrapPop());
1995: if (xin->ops->set) {
1996: PetscUseTypeMethod(xin, set, inf);
1997: } else {
1998: PetscInt n;
1999: PetscScalar *xx;
2001: PetscCall(VecGetLocalSize(xin, &n));
2002: PetscCall(VecGetArrayWrite(xin, &xx));
2003: for (PetscInt i = 0; i < n; ++i) xx[i] = inf;
2004: PetscCall(VecRestoreArrayWrite(xin, &xx));
2005: }
2006: PetscFunctionReturn(PETSC_SUCCESS);
2007: }
2009: /*@
2010: VecBindToCPU - marks a vector to temporarily stay on the CPU and perform computations on the CPU
2012: Logically collective
2014: Input Parameters:
2015: + v - the vector
2016: - flg - bind to the CPU if value of `PETSC_TRUE`
2018: Level: intermediate
2020: .seealso: [](ch_vectors), `Vec`, `VecBoundToCPU()`
2021: @*/
2022: PetscErrorCode VecBindToCPU(Vec v, PetscBool flg)
2023: {
2024: PetscFunctionBegin;
2027: #if defined(PETSC_HAVE_DEVICE)
2028: if (v->boundtocpu == flg) PetscFunctionReturn(PETSC_SUCCESS);
2029: v->boundtocpu = flg;
2030: PetscTryTypeMethod(v, bindtocpu, flg);
2031: #endif
2032: PetscFunctionReturn(PETSC_SUCCESS);
2033: }
2035: /*@
2036: VecBoundToCPU - query if a vector is bound to the CPU
2038: Not collective
2040: Input Parameter:
2041: . v - the vector
2043: Output Parameter:
2044: . flg - the logical flag
2046: Level: intermediate
2048: .seealso: [](ch_vectors), `Vec`, `VecBindToCPU()`
2049: @*/
2050: PetscErrorCode VecBoundToCPU(Vec v, PetscBool *flg)
2051: {
2052: PetscFunctionBegin;
2054: PetscAssertPointer(flg, 2);
2055: #if defined(PETSC_HAVE_DEVICE)
2056: *flg = v->boundtocpu;
2057: #else
2058: *flg = PETSC_TRUE;
2059: #endif
2060: PetscFunctionReturn(PETSC_SUCCESS);
2061: }
2063: /*@
2064: VecSetBindingPropagates - Sets whether the state of being bound to the CPU for a GPU vector type propagates to child and some other associated objects
2066: Input Parameters:
2067: + v - the vector
2068: - flg - flag indicating whether the boundtocpu flag should be propagated
2070: Level: developer
2072: Notes:
2073: 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.
2074: The created vectors will also have their bindingpropagates flag set to true.
2076: Developer Notes:
2077: If a `DMDA` has the `-dm_bind_below option` set to true, then vectors created by `DMCreateGlobalVector()` will have `VecSetBindingPropagates()` called on them to
2078: set their bindingpropagates flag to true.
2080: .seealso: [](ch_vectors), `Vec`, `MatSetBindingPropagates()`, `VecGetBindingPropagates()`
2081: @*/
2082: PetscErrorCode VecSetBindingPropagates(Vec v, PetscBool flg)
2083: {
2084: PetscFunctionBegin;
2086: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA) || defined(PETSC_HAVE_HIP)
2087: v->bindingpropagates = flg;
2088: #endif
2089: PetscFunctionReturn(PETSC_SUCCESS);
2090: }
2092: /*@
2093: VecGetBindingPropagates - Gets whether the state of being bound to the CPU for a GPU vector type propagates to child and some other associated objects
2095: Input Parameter:
2096: . v - the vector
2098: Output Parameter:
2099: . flg - flag indicating whether the boundtocpu flag will be propagated
2101: Level: developer
2103: .seealso: [](ch_vectors), `Vec`, `VecSetBindingPropagates()`
2104: @*/
2105: PetscErrorCode VecGetBindingPropagates(Vec v, PetscBool *flg)
2106: {
2107: PetscFunctionBegin;
2109: PetscAssertPointer(flg, 2);
2110: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA) || defined(PETSC_HAVE_HIP)
2111: *flg = v->bindingpropagates;
2112: #else
2113: *flg = PETSC_FALSE;
2114: #endif
2115: PetscFunctionReturn(PETSC_SUCCESS);
2116: }
2118: /*@C
2119: VecSetPinnedMemoryMin - Set the minimum data size for which pinned memory will be used for host (CPU) allocations.
2121: Logically Collective
2123: Input Parameters:
2124: + v - the vector
2125: - mbytes - minimum data size in bytes
2127: Options Database Key:
2128: . -vec_pinned_memory_min <size> - minimum size (in bytes) for an allocation to use pinned memory on host.
2130: Level: developer
2132: Note:
2133: Specifying -1 ensures that pinned memory will never be used.
2135: .seealso: [](ch_vectors), `Vec`, `VecGetPinnedMemoryMin()`
2136: @*/
2137: PetscErrorCode VecSetPinnedMemoryMin(Vec v, size_t mbytes)
2138: {
2139: PetscFunctionBegin;
2141: #if PetscDefined(HAVE_DEVICE)
2142: v->minimum_bytes_pinned_memory = mbytes;
2143: #endif
2144: PetscFunctionReturn(PETSC_SUCCESS);
2145: }
2147: /*@C
2148: VecGetPinnedMemoryMin - Get the minimum data size for which pinned memory will be used for host (CPU) allocations.
2150: Logically Collective
2152: Input Parameter:
2153: . v - the vector
2155: Output Parameter:
2156: . mbytes - minimum data size in bytes
2158: Level: developer
2160: .seealso: [](ch_vectors), `Vec`, `VecSetPinnedMemoryMin()`
2161: @*/
2162: PetscErrorCode VecGetPinnedMemoryMin(Vec v, size_t *mbytes)
2163: {
2164: PetscFunctionBegin;
2166: PetscAssertPointer(mbytes, 2);
2167: #if PetscDefined(HAVE_DEVICE)
2168: *mbytes = v->minimum_bytes_pinned_memory;
2169: #endif
2170: PetscFunctionReturn(PETSC_SUCCESS);
2171: }
2173: /*@
2174: VecGetOffloadMask - Get the offload mask of a `Vec`
2176: Not Collective
2178: Input Parameter:
2179: . v - the vector
2181: Output Parameter:
2182: . mask - corresponding `PetscOffloadMask` enum value.
2184: Level: intermediate
2186: .seealso: [](ch_vectors), `Vec`, `VecCreateSeqCUDA()`, `VecCreateSeqViennaCL()`, `VecGetArray()`, `VecGetType()`
2187: @*/
2188: PetscErrorCode VecGetOffloadMask(Vec v, PetscOffloadMask *mask)
2189: {
2190: PetscFunctionBegin;
2192: PetscAssertPointer(mask, 2);
2193: *mask = v->offloadmask;
2194: PetscFunctionReturn(PETSC_SUCCESS);
2195: }
2197: #if !defined(PETSC_HAVE_VIENNACL)
2198: PETSC_EXTERN PetscErrorCode VecViennaCLGetCLContext(Vec v, PETSC_UINTPTR_T *ctx)
2199: {
2200: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "PETSc must be configured with --with-opencl to get a Vec's cl_context");
2201: }
2203: PETSC_EXTERN PetscErrorCode VecViennaCLGetCLQueue(Vec v, PETSC_UINTPTR_T *queue)
2204: {
2205: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "PETSc must be configured with --with-opencl to get a Vec's cl_command_queue");
2206: }
2208: PETSC_EXTERN PetscErrorCode VecViennaCLGetCLMem(Vec v, PETSC_UINTPTR_T *queue)
2209: {
2210: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "PETSc must be configured with --with-opencl to get a Vec's cl_mem");
2211: }
2213: PETSC_EXTERN PetscErrorCode VecViennaCLGetCLMemRead(Vec v, PETSC_UINTPTR_T *queue)
2214: {
2215: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "PETSc must be configured with --with-opencl to get a Vec's cl_mem");
2216: }
2218: PETSC_EXTERN PetscErrorCode VecViennaCLGetCLMemWrite(Vec v, PETSC_UINTPTR_T *queue)
2219: {
2220: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "PETSc must be configured with --with-opencl to get a Vec's cl_mem");
2221: }
2223: PETSC_EXTERN PetscErrorCode VecViennaCLRestoreCLMemWrite(Vec v)
2224: {
2225: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "PETSc must be configured with --with-opencl to restore a Vec's cl_mem");
2226: }
2227: #endif
2229: static PetscErrorCode VecErrorWeightedNorms_Basic(Vec U, Vec Y, Vec E, NormType wnormtype, PetscReal atol, Vec vatol, PetscReal rtol, Vec vrtol, PetscReal ignore_max, PetscReal *norm, PetscInt *norm_loc, PetscReal *norma, PetscInt *norma_loc, PetscReal *normr, PetscInt *normr_loc)
2230: {
2231: const PetscScalar *u, *y;
2232: const PetscScalar *atola = NULL, *rtola = NULL, *erra = NULL;
2233: PetscInt n, n_loc = 0, na_loc = 0, nr_loc = 0;
2234: PetscReal nrm = 0, nrma = 0, nrmr = 0, err_loc[6];
2236: PetscFunctionBegin;
2237: #define SkipSmallValue(a, b, tol) \
2238: if (PetscAbsScalar(a) < tol || PetscAbsScalar(b) < tol) continue
2240: PetscCall(VecGetLocalSize(U, &n));
2241: PetscCall(VecGetArrayRead(U, &u));
2242: PetscCall(VecGetArrayRead(Y, &y));
2243: if (E) PetscCall(VecGetArrayRead(E, &erra));
2244: if (vatol) PetscCall(VecGetArrayRead(vatol, &atola));
2245: if (vrtol) PetscCall(VecGetArrayRead(vrtol, &rtola));
2246: for (PetscInt i = 0; i < n; i++) {
2247: PetscReal err, tol, tola, tolr;
2249: SkipSmallValue(y[i], u[i], ignore_max);
2250: atol = atola ? PetscRealPart(atola[i]) : atol;
2251: rtol = rtola ? PetscRealPart(rtola[i]) : rtol;
2252: err = erra ? PetscAbsScalar(erra[i]) : PetscAbsScalar(y[i] - u[i]);
2253: tola = atol;
2254: tolr = rtol * PetscMax(PetscAbsScalar(u[i]), PetscAbsScalar(y[i]));
2255: tol = tola + tolr;
2256: if (tola > 0.) {
2257: if (wnormtype == NORM_INFINITY) nrma = PetscMax(nrma, err / tola);
2258: else nrma += PetscSqr(err / tola);
2259: na_loc++;
2260: }
2261: if (tolr > 0.) {
2262: if (wnormtype == NORM_INFINITY) nrmr = PetscMax(nrmr, err / tolr);
2263: else nrmr += PetscSqr(err / tolr);
2264: nr_loc++;
2265: }
2266: if (tol > 0.) {
2267: if (wnormtype == NORM_INFINITY) nrm = PetscMax(nrm, err / tol);
2268: else nrm += PetscSqr(err / tol);
2269: n_loc++;
2270: }
2271: }
2272: if (E) PetscCall(VecRestoreArrayRead(E, &erra));
2273: if (vatol) PetscCall(VecRestoreArrayRead(vatol, &atola));
2274: if (vrtol) PetscCall(VecRestoreArrayRead(vrtol, &rtola));
2275: PetscCall(VecRestoreArrayRead(U, &u));
2276: PetscCall(VecRestoreArrayRead(Y, &y));
2277: #undef SkipSmallValue
2279: err_loc[0] = nrm;
2280: err_loc[1] = nrma;
2281: err_loc[2] = nrmr;
2282: err_loc[3] = (PetscReal)n_loc;
2283: err_loc[4] = (PetscReal)na_loc;
2284: err_loc[5] = (PetscReal)nr_loc;
2285: if (wnormtype == NORM_2) {
2286: PetscCall(MPIU_Allreduce(MPI_IN_PLACE, err_loc, 6, MPIU_REAL, MPIU_SUM, PetscObjectComm((PetscObject)U)));
2287: } else {
2288: PetscCall(MPIU_Allreduce(MPI_IN_PLACE, err_loc, 3, MPIU_REAL, MPIU_MAX, PetscObjectComm((PetscObject)U)));
2289: PetscCall(MPIU_Allreduce(MPI_IN_PLACE, err_loc + 3, 3, MPIU_REAL, MPIU_SUM, PetscObjectComm((PetscObject)U)));
2290: }
2291: if (wnormtype == NORM_2) {
2292: *norm = PetscSqrtReal(err_loc[0]);
2293: *norma = PetscSqrtReal(err_loc[1]);
2294: *normr = PetscSqrtReal(err_loc[2]);
2295: } else {
2296: *norm = err_loc[0];
2297: *norma = err_loc[1];
2298: *normr = err_loc[2];
2299: }
2300: *norm_loc = (PetscInt)err_loc[3];
2301: *norma_loc = (PetscInt)err_loc[4];
2302: *normr_loc = (PetscInt)err_loc[5];
2303: PetscFunctionReturn(PETSC_SUCCESS);
2304: }
2306: /*@
2307: VecErrorWeightedNorms - compute a weighted norm of the difference between two vectors
2309: Collective
2311: Input Parameters:
2312: + U - first vector to be compared
2313: . Y - second vector to be compared
2314: . E - optional third vector representing the error (if not provided, the error is ||U-Y||)
2315: . wnormtype - norm type
2316: . atol - scalar for absolute tolerance
2317: . vatol - vector representing per-entry absolute tolerances (can be ``NULL``)
2318: . rtol - scalar for relative tolerance
2319: . vrtol - vector representing per-entry relative tolerances (can be ``NULL``)
2320: - ignore_max - ignore values smaller then this value in absolute terms.
2322: Output Parameters:
2323: + norm - weighted norm
2324: . norm_loc - number of vector locations used for the weighted norm
2325: . norma - weighted norm based on the absolute tolerance
2326: . norma_loc - number of vector locations used for the absolute weighted norm
2327: . normr - weighted norm based on the relative tolerance
2328: - normr_loc - number of vector locations used for the relative weighted norm
2330: Level: developer
2332: Notes:
2333: This is primarily used for computing weighted local truncation errors in ``TS``.
2335: .seealso: [](ch_vectors), `Vec`, `NormType`, ``TSErrorWeightedNorm()``, ``TSErrorWeightedENorm()``
2336: @*/
2337: PetscErrorCode VecErrorWeightedNorms(Vec U, Vec Y, Vec E, NormType wnormtype, PetscReal atol, Vec vatol, PetscReal rtol, Vec vrtol, PetscReal ignore_max, PetscReal *norm, PetscInt *norm_loc, PetscReal *norma, PetscInt *norma_loc, PetscReal *normr, PetscInt *normr_loc)
2338: {
2339: PetscFunctionBegin;
2344: if (E) {
2347: }
2350: if (vatol) {
2353: }
2355: if (vrtol) {
2358: }
2360: PetscAssertPointer(norm, 10);
2361: PetscAssertPointer(norm_loc, 11);
2362: PetscAssertPointer(norma, 12);
2363: PetscAssertPointer(norma_loc, 13);
2364: PetscAssertPointer(normr, 14);
2365: PetscAssertPointer(normr_loc, 15);
2366: PetscCheck(wnormtype == NORM_2 || wnormtype == NORM_INFINITY, PetscObjectComm((PetscObject)U), PETSC_ERR_SUP, "No support for norm type %s", NormTypes[wnormtype]);
2368: /* There are potentially 5 vectors involved, some of them may happen to be of different type or bound to cpu.
2369: Here we check that they all implement the same operation and call it if so.
2370: Otherwise, we call the _Basic implementation that always works (provided VecGetArrayRead is implemented). */
2371: PetscBool sameop = (PetscBool)(U->ops->errorwnorm && U->ops->errorwnorm == Y->ops->errorwnorm);
2372: if (sameop && E) sameop = (PetscBool)(U->ops->errorwnorm == E->ops->errorwnorm);
2373: if (sameop && vatol) sameop = (PetscBool)(U->ops->errorwnorm == vatol->ops->errorwnorm);
2374: if (sameop && vrtol) sameop = (PetscBool)(U->ops->errorwnorm == vrtol->ops->errorwnorm);
2375: if (sameop) PetscUseTypeMethod(U, errorwnorm, Y, E, wnormtype, atol, vatol, rtol, vrtol, ignore_max, norm, norm_loc, norma, norma_loc, normr, normr_loc);
2376: else PetscCall(VecErrorWeightedNorms_Basic(U, Y, E, wnormtype, atol, vatol, rtol, vrtol, ignore_max, norm, norm_loc, norma, norma_loc, normr, normr_loc));
2377: PetscFunctionReturn(PETSC_SUCCESS);
2378: }