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>
7: /* Logging support */
8: PetscClassId VEC_CLASSID;
9: PetscLogEvent VEC_View, VEC_Max, VEC_Min, VEC_Dot, VEC_MDot, VEC_TDot;
10: PetscLogEvent VEC_Norm, VEC_Normalize, VEC_Scale, VEC_Copy, VEC_Set, VEC_AXPY, VEC_AYPX, VEC_WAXPY;
11: PetscLogEvent VEC_MTDot, VEC_MAXPY, VEC_Swap, VEC_AssemblyBegin, VEC_ScatterBegin, VEC_ScatterEnd;
12: PetscLogEvent VEC_AssemblyEnd, VEC_PointwiseMult, VEC_SetValues, VEC_Load;
13: PetscLogEvent VEC_SetRandom, VEC_ReduceArithmetic, VEC_ReduceCommunication,VEC_ReduceBegin,VEC_ReduceEnd,VEC_Ops;
14: PetscLogEvent VEC_DotNorm2, VEC_AXPBYPCZ;
15: PetscLogEvent VEC_ViennaCLCopyFromGPU, VEC_ViennaCLCopyToGPU;
16: PetscLogEvent VEC_CUDACopyFromGPU, VEC_CUDACopyToGPU;
17: PetscLogEvent VEC_CUDACopyFromGPUSome, VEC_CUDACopyToGPUSome;
18: PetscLogEvent VEC_HIPCopyFromGPU, VEC_HIPCopyToGPU;
19: PetscLogEvent VEC_HIPCopyFromGPUSome, VEC_HIPCopyToGPUSome;
21: /*@
22: VecStashGetInfo - Gets how many values are currently in the vector stash, i.e. need
23: to be communicated to other processors during the VecAssemblyBegin/End() process
25: Not collective
27: Input Parameter:
28: . vec - the vector
30: Output Parameters:
31: + nstash - the size of the stash
32: . reallocs - the number of additional mallocs incurred.
33: . bnstash - the size of the block stash
34: - breallocs - the number of additional mallocs incurred.in the block stash
36: Level: advanced
38: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), Vec, VecStashSetInitialSize(), VecStashView()
40: @*/
41: PetscErrorCode VecStashGetInfo(Vec vec,PetscInt *nstash,PetscInt *reallocs,PetscInt *bnstash,PetscInt *breallocs)
42: {
46: VecStashGetInfo_Private(&vec->stash,nstash,reallocs);
47: VecStashGetInfo_Private(&vec->bstash,bnstash,breallocs);
48: return(0);
49: }
51: /*@
52: VecSetLocalToGlobalMapping - Sets a local numbering to global numbering used
53: by the routine VecSetValuesLocal() to allow users to insert vector entries
54: using a local (per-processor) numbering.
56: Logically Collective on Vec
58: Input Parameters:
59: + x - vector
60: - mapping - mapping created with ISLocalToGlobalMappingCreate() or ISLocalToGlobalMappingCreateIS()
62: Notes:
63: All vectors obtained with VecDuplicate() from this vector inherit the same mapping.
65: Level: intermediate
67: seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetValuesLocal(),
68: VecSetLocalToGlobalMapping(), VecSetValuesBlockedLocal()
69: @*/
70: PetscErrorCode VecSetLocalToGlobalMapping(Vec x,ISLocalToGlobalMapping mapping)
71: {
78: if (x->ops->setlocaltoglobalmapping) {
79: (*x->ops->setlocaltoglobalmapping)(x,mapping);
80: } else {
81: PetscLayoutSetISLocalToGlobalMapping(x->map,mapping);
82: }
83: return(0);
84: }
86: /*@
87: VecGetLocalToGlobalMapping - Gets the local-to-global numbering set by VecSetLocalToGlobalMapping()
89: Not Collective
91: Input Parameter:
92: . X - the vector
94: Output Parameter:
95: . mapping - the mapping
97: Level: advanced
99: .seealso: VecSetValuesLocal()
100: @*/
101: PetscErrorCode VecGetLocalToGlobalMapping(Vec X,ISLocalToGlobalMapping *mapping)
102: {
107: *mapping = X->map->mapping;
108: return(0);
109: }
111: /*@
112: VecAssemblyBegin - Begins assembling the vector. This routine should
113: be called after completing all calls to VecSetValues().
115: Collective on Vec
117: Input Parameter:
118: . vec - the vector
120: Level: beginner
122: .seealso: VecAssemblyEnd(), VecSetValues()
123: @*/
124: PetscErrorCode VecAssemblyBegin(Vec vec)
125: {
131: VecStashViewFromOptions(vec,NULL,"-vec_view_stash");
132: PetscLogEventBegin(VEC_AssemblyBegin,vec,0,0,0);
133: if (vec->ops->assemblybegin) {
134: (*vec->ops->assemblybegin)(vec);
135: }
136: PetscLogEventEnd(VEC_AssemblyBegin,vec,0,0,0);
137: PetscObjectStateIncrease((PetscObject)vec);
138: return(0);
139: }
141: /*@
142: VecAssemblyEnd - Completes assembling the vector. This routine should
143: be called after VecAssemblyBegin().
145: Collective on Vec
147: Input Parameter:
148: . vec - the vector
150: Options Database Keys:
151: + -vec_view - Prints vector in ASCII format
152: . -vec_view ::ascii_matlab - Prints vector in ASCII MATLAB format to stdout
153: . -vec_view matlab:filename - Prints vector in MATLAB format to matlaboutput.mat
154: . -vec_view draw - Activates vector viewing using drawing tools
155: . -display <name> - Sets display name (default is host)
156: . -draw_pause <sec> - Sets number of seconds to pause after display
157: - -vec_view socket - Activates vector viewing using a socket
159: Level: beginner
161: .seealso: VecAssemblyBegin(), VecSetValues()
162: @*/
163: PetscErrorCode VecAssemblyEnd(Vec vec)
164: {
169: PetscLogEventBegin(VEC_AssemblyEnd,vec,0,0,0);
171: if (vec->ops->assemblyend) {
172: (*vec->ops->assemblyend)(vec);
173: }
174: PetscLogEventEnd(VEC_AssemblyEnd,vec,0,0,0);
175: VecViewFromOptions(vec,NULL,"-vec_view");
176: return(0);
177: }
179: /*@
180: VecPointwiseMax - Computes the componentwise maximum w_i = max(x_i, y_i).
182: Logically Collective on Vec
184: Input Parameters:
185: . x, y - the vectors
187: Output Parameter:
188: . w - the result
190: Level: advanced
192: Notes:
193: any subset of the x, y, and w may be the same vector.
194: For complex numbers compares only the real part
196: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
197: @*/
198: PetscErrorCode VecPointwiseMax(Vec w,Vec x,Vec y)
199: {
211: VecCheckSameSize(w,1,x,2);
212: VecCheckSameSize(w,1,y,3);
213: VecSetErrorIfLocked(w,1);
214: (*w->ops->pointwisemax)(w,x,y);
215: PetscObjectStateIncrease((PetscObject)w);
216: return(0);
217: }
219: /*@
220: VecPointwiseMin - Computes the componentwise minimum w_i = min(x_i, y_i).
222: Logically Collective on Vec
224: Input Parameters:
225: . x, y - the vectors
227: Output Parameter:
228: . w - the result
230: Level: advanced
232: Notes:
233: any subset of the x, y, and w may be the same vector.
234: For complex numbers compares only the real part
236: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
237: @*/
238: PetscErrorCode VecPointwiseMin(Vec w,Vec x,Vec y)
239: {
251: VecCheckSameSize(w,1,x,2);
252: VecCheckSameSize(w,1,y,3);
253: VecSetErrorIfLocked(w,1);
254: (*w->ops->pointwisemin)(w,x,y);
255: PetscObjectStateIncrease((PetscObject)w);
256: return(0);
257: }
259: /*@
260: VecPointwiseMaxAbs - Computes the componentwise maximum of the absolute values w_i = max(abs(x_i), abs(y_i)).
262: Logically Collective on Vec
264: Input Parameters:
265: . x, y - the vectors
267: Output Parameter:
268: . w - the result
270: Level: advanced
272: Notes:
273: any subset of the x, y, and w may be the same vector.
275: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMin(), VecPointwiseMax(), VecMaxPointwiseDivide()
276: @*/
277: PetscErrorCode VecPointwiseMaxAbs(Vec w,Vec x,Vec y)
278: {
290: VecCheckSameSize(w,1,x,2);
291: VecCheckSameSize(w,1,y,3);
292: VecSetErrorIfLocked(w,1);
293: (*w->ops->pointwisemaxabs)(w,x,y);
294: PetscObjectStateIncrease((PetscObject)w);
295: return(0);
296: }
298: /*@
299: VecPointwiseDivide - Computes the componentwise division w = x/y.
301: Logically Collective on Vec
303: Input Parameters:
304: . x, y - the vectors
306: Output Parameter:
307: . w - the result
309: Level: advanced
311: Notes:
312: any subset of the x, y, and w may be the same vector.
314: .seealso: VecPointwiseMult(), VecPointwiseMax(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
315: @*/
316: PetscErrorCode VecPointwiseDivide(Vec w,Vec x,Vec y)
317: {
329: VecCheckSameSize(w,1,x,2);
330: VecCheckSameSize(w,1,y,3);
331: VecSetErrorIfLocked(w,1);
332: (*w->ops->pointwisedivide)(w,x,y);
333: PetscObjectStateIncrease((PetscObject)w);
334: return(0);
335: }
337: /*@
338: VecDuplicate - Creates a new vector of the same type as an existing vector.
340: Collective on Vec
342: Input Parameters:
343: . v - a vector to mimic
345: Output Parameter:
346: . newv - location to put new vector
348: Notes:
349: VecDuplicate() DOES NOT COPY the vector entries, but rather allocates storage
350: for the new vector. Use VecCopy() to copy a vector.
352: Use VecDestroy() to free the space. Use VecDuplicateVecs() to get several
353: vectors.
355: Level: beginner
357: .seealso: VecDestroy(), VecDuplicateVecs(), VecCreate(), VecCopy()
358: @*/
359: PetscErrorCode VecDuplicate(Vec v,Vec *newv)
360: {
367: (*v->ops->duplicate)(v,newv);
368: PetscObjectStateIncrease((PetscObject)*newv);
369: return(0);
370: }
372: /*@C
373: VecDestroy - Destroys a vector.
375: Collective on Vec
377: Input Parameters:
378: . v - the vector
380: Level: beginner
382: .seealso: VecDuplicate(), VecDestroyVecs()
383: @*/
384: PetscErrorCode VecDestroy(Vec *v)
385: {
389: if (!*v) return(0);
391: if (--((PetscObject)(*v))->refct > 0) {*v = NULL; return(0);}
393: PetscObjectSAWsViewOff((PetscObject)*v);
394: /* destroy the internal part */
395: if ((*v)->ops->destroy) {
396: (*(*v)->ops->destroy)(*v);
397: }
398: PetscFree((*v)->defaultrandtype);
399: /* destroy the external/common part */
400: PetscLayoutDestroy(&(*v)->map);
401: PetscHeaderDestroy(v);
402: return(0);
403: }
405: /*@C
406: VecDuplicateVecs - Creates several vectors of the same type as an existing vector.
408: Collective on Vec
410: Input Parameters:
411: + m - the number of vectors to obtain
412: - v - a vector to mimic
414: Output Parameter:
415: . V - location to put pointer to array of vectors
417: Notes:
418: Use VecDestroyVecs() to free the space. Use VecDuplicate() to form a single
419: vector.
421: Fortran Note:
422: The Fortran interface is slightly different from that given below, it
423: requires one to pass in V a Vec (integer) array of size at least m.
424: See the Fortran chapter of the users manual and petsc/src/vec/vec/examples for details.
426: Level: intermediate
428: .seealso: VecDestroyVecs(), VecDuplicate(), VecCreate(), VecDuplicateVecsF90()
429: @*/
430: PetscErrorCode VecDuplicateVecs(Vec v,PetscInt m,Vec *V[])
431: {
438: (*v->ops->duplicatevecs)(v,m,V);
439: return(0);
440: }
442: /*@C
443: VecDestroyVecs - Frees a block of vectors obtained with VecDuplicateVecs().
445: Collective on Vec
447: Input Parameters:
448: + vv - pointer to pointer to array of vector pointers, if NULL no vectors are destroyed
449: - m - the number of vectors previously obtained, if zero no vectors are destroyed
451: Fortran Note:
452: The Fortran interface is slightly different from that given below.
453: See the Fortran chapter of the users manual
455: Level: intermediate
457: .seealso: VecDuplicateVecs(), VecDestroyVecsf90()
458: @*/
459: PetscErrorCode VecDestroyVecs(PetscInt m,Vec *vv[])
460: {
465: if (m < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Trying to destroy negative number of vectors %D",m);
466: if (!m || !*vv) {*vv = NULL; return(0);}
469: (*(**vv)->ops->destroyvecs)(m,*vv);
470: *vv = NULL;
471: return(0);
472: }
474: /*@C
475: VecViewFromOptions - View from Options
477: Collective on Vec
479: Input Parameters:
480: + A - the vector
481: . obj - Optional object
482: - name - command line option
484: Level: intermediate
485: .seealso: Vec, VecView, PetscObjectViewFromOptions(), VecCreate()
486: @*/
487: PetscErrorCode VecViewFromOptions(Vec A,PetscObject obj,const char name[])
488: {
493: PetscObjectViewFromOptions((PetscObject)A,obj,name);
494: return(0);
495: }
497: /*@C
498: VecView - Views a vector object.
500: Collective on Vec
502: Input Parameters:
503: + vec - the vector
504: - viewer - an optional visualization context
506: Notes:
507: The available visualization contexts include
508: + PETSC_VIEWER_STDOUT_SELF - for sequential vectors
509: . PETSC_VIEWER_STDOUT_WORLD - for parallel vectors created on PETSC_COMM_WORLD
510: - PETSC_VIEWER_STDOUT_(comm) - for parallel vectors created on MPI communicator comm
512: You can change the format the vector is printed using the
513: option PetscViewerPushFormat().
515: The user can open alternative viewers with
516: + PetscViewerASCIIOpen() - Outputs vector to a specified file
517: . PetscViewerBinaryOpen() - Outputs vector in binary to a
518: specified file; corresponding input uses VecLoad()
519: . PetscViewerDrawOpen() - Outputs vector to an X window display
520: . PetscViewerSocketOpen() - Outputs vector to Socket viewer
521: - PetscViewerHDF5Open() - Outputs vector to HDF5 file viewer
523: The user can call PetscViewerPushFormat() to specify the output
524: format of ASCII printed objects (when using PETSC_VIEWER_STDOUT_SELF,
525: PETSC_VIEWER_STDOUT_WORLD and PetscViewerASCIIOpen). Available formats include
526: + PETSC_VIEWER_DEFAULT - default, prints vector contents
527: . PETSC_VIEWER_ASCII_MATLAB - prints vector contents in MATLAB format
528: . PETSC_VIEWER_ASCII_INDEX - prints vector contents, including indices of vector elements
529: - PETSC_VIEWER_ASCII_COMMON - prints vector contents, using a
530: format common among all vector types
532: Notes:
533: You can pass any number of vector objects, or other PETSc objects to the same viewer.
535: In the debugger you can do "call VecView(v,0)" to display the vector. (The same holds for any PETSc object viewer).
537: Notes for binary viewer:
538: If you pass multiple vectors to a binary viewer you can read them back in in the same order
539: with VecLoad().
541: If the blocksize of the vector is greater than one then you must provide a unique prefix to
542: the vector with PetscObjectSetOptionsPrefix((PetscObject)vec,"uniqueprefix"); BEFORE calling VecView() on the
543: vector to be stored and then set that same unique prefix on the vector that you pass to VecLoad(). The blocksize
544: information is stored in an ASCII file with the same name as the binary file plus a ".info" appended to the
545: filename. If you copy the binary file, make sure you copy the associated .info file with it.
547: See the manual page for VecLoad() on the exact format the binary viewer stores
548: the values in the file.
550: Notes for HDF5 Viewer:
551: The name of the Vec (given with PetscObjectSetName() is the name that is used
552: for the object in the HDF5 file. If you wish to store the same Vec into multiple
553: datasets in the same file (typically with different values), you must change its
554: name each time before calling the VecView(). To load the same vector,
555: the name of the Vec object passed to VecLoad() must be the same.
557: If the block size of the vector is greater than 1 then it is used as the first dimension in the HDF5 array.
558: If the function PetscViewerHDF5SetBaseDimension2()is called then even if the block size is one it will
559: be used as the first dimension in the HDF5 array (that is the HDF5 array will always be two dimensional)
560: See also PetscViewerHDF5SetTimestep() which adds an additional complication to reading and writing Vecs
561: with the HDF5 viewer.
563: Level: beginner
565: .seealso: PetscViewerASCIIOpen(), PetscViewerDrawOpen(), PetscDrawLGCreate(),
566: PetscViewerSocketOpen(), PetscViewerBinaryOpen(), VecLoad(), PetscViewerCreate(),
567: PetscRealView(), PetscScalarView(), PetscIntView(), PetscViewerHDF5SetTimestep()
568: @*/
569: PetscErrorCode VecView(Vec vec,PetscViewer viewer)
570: {
571: PetscErrorCode ierr;
572: PetscBool iascii;
573: PetscViewerFormat format;
574: PetscMPIInt size;
579: if (!viewer) {PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)vec),&viewer);}
581: PetscViewerGetFormat(viewer,&format);
582: MPI_Comm_size(PetscObjectComm((PetscObject)vec),&size);
583: if (size == 1 && format == PETSC_VIEWER_LOAD_BALANCE) return(0);
585: if (vec->stash.n || vec->bstash.n) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call VecAssemblyBegin/End() before viewing this vector");
587: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
588: if (iascii) {
589: PetscInt rows,bs;
591: PetscObjectPrintClassNamePrefixType((PetscObject)vec,viewer);
592: if (format == PETSC_VIEWER_ASCII_INFO || format == PETSC_VIEWER_ASCII_INFO_DETAIL) {
593: PetscViewerASCIIPushTab(viewer);
594: VecGetSize(vec,&rows);
595: VecGetBlockSize(vec,&bs);
596: if (bs != 1) {
597: PetscViewerASCIIPrintf(viewer,"length=%D, bs=%D\n",rows,bs);
598: } else {
599: PetscViewerASCIIPrintf(viewer,"length=%D\n",rows);
600: }
601: PetscViewerASCIIPopTab(viewer);
602: }
603: }
604: VecLockReadPush(vec);
605: PetscLogEventBegin(VEC_View,vec,viewer,0,0);
606: if ((format == PETSC_VIEWER_NATIVE || format == PETSC_VIEWER_LOAD_BALANCE) && vec->ops->viewnative) {
607: (*vec->ops->viewnative)(vec,viewer);
608: } else {
609: (*vec->ops->view)(vec,viewer);
610: }
611: VecLockReadPop(vec);
612: PetscLogEventEnd(VEC_View,vec,viewer,0,0);
613: return(0);
614: }
616: #if defined(PETSC_USE_DEBUG)
617: #include <../src/sys/totalview/tv_data_display.h>
618: PETSC_UNUSED static int TV_display_type(const struct _p_Vec *v)
619: {
620: const PetscScalar *values;
621: char type[32];
622: PetscErrorCode ierr;
624: TV_add_row("Local rows", "int", &v->map->n);
625: TV_add_row("Global rows", "int", &v->map->N);
626: TV_add_row("Typename", TV_ascii_string_type, ((PetscObject)v)->type_name);
627: VecGetArrayRead((Vec)v,&values);
628: PetscSNPrintf(type,32,"double[%d]",v->map->n);
629: TV_add_row("values",type, values);
630: VecRestoreArrayRead((Vec)v,&values);
631: return TV_format_OK;
632: }
633: #endif
635: /*@
636: VecGetSize - Returns the global number of elements of the vector.
638: Not Collective
640: Input Parameter:
641: . x - the vector
643: Output Parameters:
644: . size - the global length of the vector
646: Level: beginner
648: .seealso: VecGetLocalSize()
649: @*/
650: PetscErrorCode VecGetSize(Vec x,PetscInt *size)
651: {
658: (*x->ops->getsize)(x,size);
659: return(0);
660: }
662: /*@
663: VecGetLocalSize - Returns the number of elements of the vector stored
664: in local memory.
666: Not Collective
668: Input Parameter:
669: . x - the vector
671: Output Parameter:
672: . size - the length of the local piece of the vector
674: Level: beginner
676: .seealso: VecGetSize()
677: @*/
678: PetscErrorCode VecGetLocalSize(Vec x,PetscInt *size)
679: {
686: (*x->ops->getlocalsize)(x,size);
687: return(0);
688: }
690: /*@C
691: VecGetOwnershipRange - Returns the range of indices owned by
692: this processor, assuming that the vectors are laid out with the
693: first n1 elements on the first processor, next n2 elements on the
694: second, etc. For certain parallel layouts this range may not be
695: well defined.
697: Not Collective
699: Input Parameter:
700: . x - the vector
702: Output Parameters:
703: + low - the first local element, pass in NULL if not interested
704: - high - one more than the last local element, pass in NULL if not interested
706: Note:
707: The high argument is one more than the last element stored locally.
709: Fortran: PETSC_NULL_INTEGER should be used instead of NULL
711: Level: beginner
713: .seealso: MatGetOwnershipRange(), MatGetOwnershipRanges(), VecGetOwnershipRanges()
714: @*/
715: PetscErrorCode VecGetOwnershipRange(Vec x,PetscInt *low,PetscInt *high)
716: {
722: if (low) *low = x->map->rstart;
723: if (high) *high = x->map->rend;
724: return(0);
725: }
727: /*@C
728: VecGetOwnershipRanges - Returns the range of indices owned by EACH processor,
729: assuming that the vectors are laid out with the
730: first n1 elements on the first processor, next n2 elements on the
731: second, etc. For certain parallel layouts this range may not be
732: well defined.
734: Not Collective
736: Input Parameter:
737: . x - the vector
739: Output Parameters:
740: . range - array of length size+1 with the start and end+1 for each process
742: Note:
743: The high argument is one more than the last element stored locally.
745: Fortran: You must PASS in an array of length size+1
747: If the ranges are used after all vectors that share the ranges has been destroyed then the program will crash accessing ranges[].
749: Level: beginner
751: .seealso: MatGetOwnershipRange(), MatGetOwnershipRanges(), VecGetOwnershipRange()
752: @*/
753: PetscErrorCode VecGetOwnershipRanges(Vec x,const PetscInt *ranges[])
754: {
760: PetscLayoutGetRanges(x->map,ranges);
761: return(0);
762: }
764: /*@
765: VecSetOption - Sets an option for controling a vector's behavior.
767: Collective on Vec
769: Input Parameters:
770: + x - the vector
771: . op - the option
772: - flag - turn the option on or off
774: Supported Options:
775: + VEC_IGNORE_OFF_PROC_ENTRIES, which causes VecSetValues() to ignore
776: entries destined to be stored on a separate processor. This can be used
777: to eliminate the global reduction in the VecAssemblyXXXX() if you know
778: that you have only used VecSetValues() to set local elements
779: . VEC_IGNORE_NEGATIVE_INDICES, which means you can pass negative indices
780: in ix in calls to VecSetValues() or VecGetValues(). These rows are simply
781: ignored.
782: - VEC_SUBSET_OFF_PROC_ENTRIES, which causes VecAssemblyBegin() to assume that the off-process
783: entries will always be a subset (possibly equal) of the off-process entries set on the
784: first assembly which had a true VEC_SUBSET_OFF_PROC_ENTRIES and the vector has not
785: changed this flag afterwards. If this assembly is not such first assembly, then this
786: assembly can reuse the communication pattern setup in that first assembly, thus avoiding
787: a global reduction. Subsequent assemblies setting off-process values should use the same
788: InsertMode as the first assembly.
790: Developer Note:
791: The InsertMode restriction could be removed by packing the stash messages out of place.
793: Level: intermediate
795: @*/
796: PetscErrorCode VecSetOption(Vec x,VecOption op,PetscBool flag)
797: {
803: if (x->ops->setoption) {
804: (*x->ops->setoption)(x,op,flag);
805: }
806: return(0);
807: }
809: /* Default routines for obtaining and releasing; */
810: /* may be used by any implementation */
811: PetscErrorCode VecDuplicateVecs_Default(Vec w,PetscInt m,Vec *V[])
812: {
814: PetscInt i;
819: if (m <= 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"m must be > 0: m = %D",m);
820: PetscMalloc1(m,V);
821: for (i=0; i<m; i++) {VecDuplicate(w,*V+i);}
822: return(0);
823: }
825: PetscErrorCode VecDestroyVecs_Default(PetscInt m,Vec v[])
826: {
828: PetscInt i;
832: for (i=0; i<m; i++) {VecDestroy(&v[i]);}
833: PetscFree(v);
834: return(0);
835: }
837: /*@
838: VecResetArray - Resets a vector to use its default memory. Call this
839: after the use of VecPlaceArray().
841: Not Collective
843: Input Parameters:
844: . vec - the vector
846: Level: developer
848: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecPlaceArray()
850: @*/
851: PetscErrorCode VecResetArray(Vec vec)
852: {
858: if (vec->ops->resetarray) {
859: (*vec->ops->resetarray)(vec);
860: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Cannot reset array in this type of vector");
861: PetscObjectStateIncrease((PetscObject)vec);
862: return(0);
863: }
865: /*@C
866: VecLoad - Loads a vector that has been stored in binary or HDF5 format
867: with VecView().
869: Collective on PetscViewer
871: Input Parameters:
872: + vec - the newly loaded vector, this needs to have been created with VecCreate() or
873: some related function before a call to VecLoad().
874: - viewer - binary file viewer, obtained from PetscViewerBinaryOpen() or
875: HDF5 file viewer, obtained from PetscViewerHDF5Open()
877: Level: intermediate
879: Notes:
880: Defaults to the standard Seq or MPI Vec, if you want some other type of Vec call VecSetFromOptions()
881: before calling this.
883: The input file must contain the full global vector, as
884: written by the routine VecView().
886: If the type or size of vec is not set before a call to VecLoad, PETSc
887: sets the type and the local and global sizes. If type and/or
888: sizes are already set, then the same are used.
890: If using the binary viewer and the blocksize of the vector is greater than one then you must provide a unique prefix to
891: the vector with PetscObjectSetOptionsPrefix((PetscObject)vec,"uniqueprefix"); BEFORE calling VecView() on the
892: vector to be stored and then set that same unique prefix on the vector that you pass to VecLoad(). The blocksize
893: information is stored in an ASCII file with the same name as the binary file plus a ".info" appended to the
894: filename. If you copy the binary file, make sure you copy the associated .info file with it.
896: If using HDF5, you must assign the Vec the same name as was used in the Vec
897: that was stored in the file using PetscObjectSetName(). Otherwise you will
898: get the error message: "Cannot H5DOpen2() with Vec name NAMEOFOBJECT".
900: If the HDF5 file contains a two dimensional array the first dimension is treated as the block size
901: in loading the vector. Hence, for example, using Matlab notation h5create('vector.dat','/Test_Vec',[27 1]);
902: will load a vector of size 27 and block size 27 thus resulting in all 27 entries being on the first process of
903: vectors communicator and the rest of the processes having zero entries
905: Notes for advanced users when using the binary viewer:
906: Most users should not need to know the details of the binary storage
907: format, since VecLoad() and VecView() completely hide these details.
908: But for anyone who's interested, the standard binary vector storage
909: format is
910: .vb
911: PetscInt VEC_FILE_CLASSID
912: PetscInt number of rows
913: PetscScalar *values of all entries
914: .ve
916: In addition, PETSc automatically uses byte swapping to work on all machines; the files
917: are written ALWAYS using big-endian ordering. On small-endian machines the numbers
918: are converted to the small-endian format when they are read in from the file.
919: See PetscBinaryRead() and PetscBinaryWrite() to see how this may be done.
921: .seealso: PetscViewerBinaryOpen(), VecView(), MatLoad(), VecLoad()
922: @*/
923: PetscErrorCode VecLoad(Vec vec, PetscViewer viewer)
924: {
925: PetscErrorCode ierr;
926: PetscBool isbinary,ishdf5,isadios,isexodusii;
927: PetscViewerFormat format;
933: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
934: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5);
935: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERADIOS,&isadios);
936: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWEREXODUSII,&isexodusii);
937: if (!isbinary && !ishdf5 && !isadios && !isexodusii) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Invalid viewer; open viewer with PetscViewerBinaryOpen()");
939: VecSetErrorIfLocked(vec,1);
940: if (!((PetscObject)vec)->type_name && !vec->ops->create) {
941: VecSetType(vec, VECSTANDARD);
942: }
943: PetscLogEventBegin(VEC_Load,viewer,0,0,0);
944: PetscViewerGetFormat(viewer,&format);
945: if (format == PETSC_VIEWER_NATIVE && vec->ops->loadnative) {
946: (*vec->ops->loadnative)(vec,viewer);
947: } else {
948: (*vec->ops->load)(vec,viewer);
949: }
950: PetscLogEventEnd(VEC_Load,viewer,0,0,0);
951: return(0);
952: }
954: /*@
955: VecReciprocal - Replaces each component of a vector by its reciprocal.
957: Logically Collective on Vec
959: Input Parameter:
960: . vec - the vector
962: Output Parameter:
963: . vec - the vector reciprocal
965: Level: intermediate
967: .seealso: VecLog(), VecExp(), VecSqrtAbs()
969: @*/
970: PetscErrorCode VecReciprocal(Vec vec)
971: {
977: if (vec->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
978: if (!vec->ops->reciprocal) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Vector does not support reciprocal operation");
979: VecSetErrorIfLocked(vec,1);
980: (*vec->ops->reciprocal)(vec);
981: PetscObjectStateIncrease((PetscObject)vec);
982: return(0);
983: }
985: /*@C
986: VecSetOperation - Allows user to set a vector operation.
988: Logically Collective on Vec
990: Input Parameters:
991: + vec - the vector
992: . op - the name of the operation
993: - f - the function that provides the operation.
995: Level: advanced
997: Usage:
998: $ PetscErrorCode userview(Vec,PetscViewer);
999: $ VecCreateMPI(comm,m,M,&x);
1000: $ VecSetOperation(x,VECOP_VIEW,(void(*)(void))userview);
1002: Notes:
1003: See the file include/petscvec.h for a complete list of matrix
1004: operations, which all have the form VECOP_<OPERATION>, where
1005: <OPERATION> is the name (in all capital letters) of the
1006: user interface routine (e.g., VecView() -> VECOP_VIEW).
1008: This function is not currently available from Fortran.
1010: .seealso: VecCreate(), MatShellSetOperation()
1011: @*/
1012: PetscErrorCode VecSetOperation(Vec vec,VecOperation op, void (*f)(void))
1013: {
1016: if (op == VECOP_VIEW && !vec->ops->viewnative) {
1017: vec->ops->viewnative = vec->ops->view;
1018: } else if (op == VECOP_LOAD && !vec->ops->loadnative) {
1019: vec->ops->loadnative = vec->ops->load;
1020: }
1021: (((void(**)(void))vec->ops)[(int)op]) = f;
1022: return(0);
1023: }
1025: /*@
1026: VecStashSetInitialSize - sets the sizes of the vec-stash, that is
1027: used during the assembly process to store values that belong to
1028: other processors.
1030: Not Collective, different processes can have different size stashes
1032: Input Parameters:
1033: + vec - the vector
1034: . size - the initial size of the stash.
1035: - bsize - the initial size of the block-stash(if used).
1037: Options Database Keys:
1038: + -vecstash_initial_size <size> or <size0,size1,...sizep-1>
1039: - -vecstash_block_initial_size <bsize> or <bsize0,bsize1,...bsizep-1>
1041: Level: intermediate
1043: Notes:
1044: The block-stash is used for values set with VecSetValuesBlocked() while
1045: the stash is used for values set with VecSetValues()
1047: Run with the option -info and look for output of the form
1048: VecAssemblyBegin_MPIXXX:Stash has MM entries, uses nn mallocs.
1049: to determine the appropriate value, MM, to use for size and
1050: VecAssemblyBegin_MPIXXX:Block-Stash has BMM entries, uses nn mallocs.
1051: to determine the value, BMM to use for bsize
1053: .seealso: VecSetBlockSize(), VecSetValues(), VecSetValuesBlocked(), VecStashView()
1055: @*/
1056: PetscErrorCode VecStashSetInitialSize(Vec vec,PetscInt size,PetscInt bsize)
1057: {
1062: VecStashSetInitialSize_Private(&vec->stash,size);
1063: VecStashSetInitialSize_Private(&vec->bstash,bsize);
1064: return(0);
1065: }
1067: /*@
1068: VecConjugate - Conjugates a vector.
1070: Logically Collective on Vec
1072: Input Parameters:
1073: . x - the vector
1075: Level: intermediate
1077: @*/
1078: PetscErrorCode VecConjugate(Vec x)
1079: {
1080: #if defined(PETSC_USE_COMPLEX)
1086: if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1087: VecSetErrorIfLocked(x,1);
1088: (*x->ops->conjugate)(x);
1089: /* we need to copy norms here */
1090: PetscObjectStateIncrease((PetscObject)x);
1091: return(0);
1092: #else
1093: return(0);
1094: #endif
1095: }
1097: /*@
1098: VecPointwiseMult - Computes the componentwise multiplication w = x*y.
1100: Logically Collective on Vec
1102: Input Parameters:
1103: . x, y - the vectors
1105: Output Parameter:
1106: . w - the result
1108: Level: advanced
1110: Notes:
1111: any subset of the x, y, and w may be the same vector.
1113: .seealso: VecPointwiseDivide(), VecPointwiseMax(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
1114: @*/
1115: PetscErrorCode VecPointwiseMult(Vec w,Vec x,Vec y)
1116: {
1128: VecCheckSameSize(w,1,x,2);
1129: VecCheckSameSize(w,2,y,3);
1130: VecSetErrorIfLocked(w,1);
1131: PetscLogEventBegin(VEC_PointwiseMult,x,y,w,0);
1132: (*w->ops->pointwisemult)(w,x,y);
1133: PetscLogEventEnd(VEC_PointwiseMult,x,y,w,0);
1134: PetscObjectStateIncrease((PetscObject)w);
1135: return(0);
1136: }
1138: /*@
1139: VecSetRandom - Sets all components of a vector to random numbers.
1141: Logically Collective on Vec
1143: Input Parameters:
1144: + x - the vector
1145: - rctx - the random number context, formed by PetscRandomCreate(), or NULL and
1146: it will create one internally.
1148: Output Parameter:
1149: . x - the vector
1151: Example of Usage:
1152: .vb
1153: PetscRandomCreate(PETSC_COMM_WORLD,&rctx);
1154: VecSetRandom(x,rctx);
1155: PetscRandomDestroy(&rctx);
1156: .ve
1158: Level: intermediate
1160: .seealso: VecSet(), VecSetValues(), PetscRandomCreate(), PetscRandomDestroy()
1161: @*/
1162: PetscErrorCode VecSetRandom(Vec x,PetscRandom rctx)
1163: {
1165: PetscRandom randObj = NULL;
1171: if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1172: VecSetErrorIfLocked(x,1);
1174: if (!rctx) {
1175: PetscRandomCreate(PetscObjectComm((PetscObject)x),&randObj);
1176: PetscRandomSetType(randObj,x->defaultrandtype);
1177: PetscRandomSetFromOptions(randObj);
1178: rctx = randObj;
1179: }
1181: PetscLogEventBegin(VEC_SetRandom,x,rctx,0,0);
1182: (*x->ops->setrandom)(x,rctx);
1183: PetscLogEventEnd(VEC_SetRandom,x,rctx,0,0);
1185: PetscRandomDestroy(&randObj);
1186: PetscObjectStateIncrease((PetscObject)x);
1187: return(0);
1188: }
1190: /*@
1191: VecZeroEntries - puts a 0.0 in each element of a vector
1193: Logically Collective on Vec
1195: Input Parameter:
1196: . vec - The vector
1198: Level: beginner
1200: .seealso: VecCreate(), VecSetOptionsPrefix(), VecSet(), VecSetValues()
1201: @*/
1202: PetscErrorCode VecZeroEntries(Vec vec)
1203: {
1207: VecSet(vec,0);
1208: return(0);
1209: }
1211: /*
1212: VecSetTypeFromOptions_Private - Sets the type of vector from user options. Defaults to a PETSc sequential vector on one
1213: processor and a PETSc MPI vector on more than one processor.
1215: Collective on Vec
1217: Input Parameter:
1218: . vec - The vector
1220: Level: intermediate
1222: .seealso: VecSetFromOptions(), VecSetType()
1223: */
1224: static PetscErrorCode VecSetTypeFromOptions_Private(PetscOptionItems *PetscOptionsObject,Vec vec)
1225: {
1226: PetscBool opt;
1227: VecType defaultType;
1228: char typeName[256];
1229: PetscMPIInt size;
1233: if (((PetscObject)vec)->type_name) defaultType = ((PetscObject)vec)->type_name;
1234: else {
1235: MPI_Comm_size(PetscObjectComm((PetscObject)vec), &size);
1236: if (size > 1) defaultType = VECMPI;
1237: else defaultType = VECSEQ;
1238: }
1240: VecRegisterAll();
1241: PetscOptionsFList("-vec_type","Vector type","VecSetType",VecList,defaultType,typeName,256,&opt);
1242: if (opt) {
1243: VecSetType(vec, typeName);
1244: } else {
1245: VecSetType(vec, defaultType);
1246: }
1247: return(0);
1248: }
1250: /*@
1251: VecSetFromOptions - Configures the vector from the options database.
1253: Collective on Vec
1255: Input Parameter:
1256: . vec - The vector
1258: Notes:
1259: To see all options, run your program with the -help option, or consult the users manual.
1260: Must be called after VecCreate() but before the vector is used.
1262: Level: beginner
1264: .seealso: VecCreate(), VecSetOptionsPrefix()
1265: @*/
1266: PetscErrorCode VecSetFromOptions(Vec vec)
1267: {
1273: PetscObjectOptionsBegin((PetscObject)vec);
1274: /* Handle vector type options */
1275: VecSetTypeFromOptions_Private(PetscOptionsObject,vec);
1277: /* Handle specific vector options */
1278: if (vec->ops->setfromoptions) {
1279: (*vec->ops->setfromoptions)(PetscOptionsObject,vec);
1280: }
1282: /* process any options handlers added with PetscObjectAddOptionsHandler() */
1283: PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject)vec);
1284: PetscOptionsEnd();
1285: return(0);
1286: }
1288: /*@
1289: VecSetSizes - Sets the local and global sizes, and checks to determine compatibility
1291: Collective on Vec
1293: Input Parameters:
1294: + v - the vector
1295: . n - the local size (or PETSC_DECIDE to have it set)
1296: - N - the global size (or PETSC_DECIDE)
1298: Notes:
1299: n and N cannot be both PETSC_DECIDE
1300: If one processor calls this with N of PETSC_DECIDE then all processors must, otherwise the program will hang.
1302: Level: intermediate
1304: .seealso: VecGetSize(), PetscSplitOwnership()
1305: @*/
1306: PetscErrorCode VecSetSizes(Vec v, PetscInt n, PetscInt N)
1307: {
1313: if (N >= 0 && n > N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Local size %D cannot be larger than global size %D",n,N);
1314: if ((v->map->n >= 0 || v->map->N >= 0) && (v->map->n != n || v->map->N != N)) SETERRQ4(PETSC_COMM_SELF,PETSC_ERR_SUP,"Cannot change/reset vector sizes to %D local %D global after previously setting them to %D local %D global",n,N,v->map->n,v->map->N);
1315: v->map->n = n;
1316: v->map->N = N;
1317: if (v->ops->create) {
1318: (*v->ops->create)(v);
1319: v->ops->create = NULL;
1320: }
1321: return(0);
1322: }
1324: /*@
1325: VecSetBlockSize - Sets the blocksize for future calls to VecSetValuesBlocked()
1326: and VecSetValuesBlockedLocal().
1328: Logically Collective on Vec
1330: Input Parameters:
1331: + v - the vector
1332: - bs - the blocksize
1334: Notes:
1335: All vectors obtained by VecDuplicate() inherit the same blocksize.
1337: Level: advanced
1339: .seealso: VecSetValuesBlocked(), VecSetLocalToGlobalMapping(), VecGetBlockSize()
1341: @*/
1342: PetscErrorCode VecSetBlockSize(Vec v,PetscInt bs)
1343: {
1349: PetscLayoutSetBlockSize(v->map,bs);
1350: v->bstash.bs = bs; /* use the same blocksize for the vec's block-stash */
1351: return(0);
1352: }
1354: /*@
1355: VecGetBlockSize - Gets the blocksize for the vector, i.e. what is used for VecSetValuesBlocked()
1356: and VecSetValuesBlockedLocal().
1358: Not Collective
1360: Input Parameter:
1361: . v - the vector
1363: Output Parameter:
1364: . bs - the blocksize
1366: Notes:
1367: All vectors obtained by VecDuplicate() inherit the same blocksize.
1369: Level: advanced
1371: .seealso: VecSetValuesBlocked(), VecSetLocalToGlobalMapping(), VecSetBlockSize()
1373: @*/
1374: PetscErrorCode VecGetBlockSize(Vec v,PetscInt *bs)
1375: {
1381: PetscLayoutGetBlockSize(v->map,bs);
1382: return(0);
1383: }
1385: /*@C
1386: VecSetOptionsPrefix - Sets the prefix used for searching for all
1387: Vec options in the database.
1389: Logically Collective on Vec
1391: Input Parameters:
1392: + v - the Vec context
1393: - prefix - the prefix to prepend to all option names
1395: Notes:
1396: A hyphen (-) must NOT be given at the beginning of the prefix name.
1397: The first character of all runtime options is AUTOMATICALLY the hyphen.
1399: Level: advanced
1401: .seealso: VecSetFromOptions()
1402: @*/
1403: PetscErrorCode VecSetOptionsPrefix(Vec v,const char prefix[])
1404: {
1409: PetscObjectSetOptionsPrefix((PetscObject)v,prefix);
1410: return(0);
1411: }
1413: /*@C
1414: VecAppendOptionsPrefix - Appends to the prefix used for searching for all
1415: Vec options in the database.
1417: Logically Collective on Vec
1419: Input Parameters:
1420: + v - the Vec context
1421: - prefix - the prefix to prepend to all option names
1423: Notes:
1424: A hyphen (-) must NOT be given at the beginning of the prefix name.
1425: The first character of all runtime options is AUTOMATICALLY the hyphen.
1427: Level: advanced
1429: .seealso: VecGetOptionsPrefix()
1430: @*/
1431: PetscErrorCode VecAppendOptionsPrefix(Vec v,const char prefix[])
1432: {
1437: PetscObjectAppendOptionsPrefix((PetscObject)v,prefix);
1438: return(0);
1439: }
1441: /*@C
1442: VecGetOptionsPrefix - Sets the prefix used for searching for all
1443: Vec options in the database.
1445: Not Collective
1447: Input Parameter:
1448: . v - the Vec context
1450: Output Parameter:
1451: . prefix - pointer to the prefix string used
1453: Notes:
1454: On the fortran side, the user should pass in a string 'prefix' of
1455: sufficient length to hold the prefix.
1457: Level: advanced
1459: .seealso: VecAppendOptionsPrefix()
1460: @*/
1461: PetscErrorCode VecGetOptionsPrefix(Vec v,const char *prefix[])
1462: {
1467: PetscObjectGetOptionsPrefix((PetscObject)v,prefix);
1468: return(0);
1469: }
1471: /*@
1472: VecSetUp - Sets up the internal vector data structures for the later use.
1474: Collective on Vec
1476: Input Parameters:
1477: . v - the Vec context
1479: Notes:
1480: For basic use of the Vec classes the user need not explicitly call
1481: VecSetUp(), since these actions will happen automatically.
1483: Level: advanced
1485: .seealso: VecCreate(), VecDestroy()
1486: @*/
1487: PetscErrorCode VecSetUp(Vec v)
1488: {
1489: PetscMPIInt size;
1494: if (v->map->n < 0 && v->map->N < 0) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Sizes not set");
1495: if (!((PetscObject)v)->type_name) {
1496: MPI_Comm_size(PetscObjectComm((PetscObject)v), &size);
1497: if (size == 1) {
1498: VecSetType(v, VECSEQ);
1499: } else {
1500: VecSetType(v, VECMPI);
1501: }
1502: }
1503: return(0);
1504: }
1506: /*
1507: These currently expose the PetscScalar/PetscReal in updating the
1508: cached norm. If we push those down into the implementation these
1509: will become independent of PetscScalar/PetscReal
1510: */
1512: /*@
1513: VecCopy - Copies a vector. y <- x
1515: Logically Collective on Vec
1517: Input Parameter:
1518: . x - the vector
1520: Output Parameter:
1521: . y - the copy
1523: Notes:
1524: For default parallel PETSc vectors, both x and y must be distributed in
1525: the same manner; local copies are done.
1527: Developer Notes:
1529: of the vectors to be sequential and one to be parallel so long as both have the same
1530: local sizes. This is used in some internal functions in PETSc.
1532: Level: beginner
1534: .seealso: VecDuplicate()
1535: @*/
1536: PetscErrorCode VecCopy(Vec x,Vec y)
1537: {
1538: PetscBool flgs[4];
1539: PetscReal norms[4] = {0.0,0.0,0.0,0.0};
1541: PetscInt i;
1548: if (x == y) return(0);
1549: VecCheckSameLocalSize(x,1,y,2);
1550: if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1551: VecSetErrorIfLocked(y,2);
1553: #if !defined(PETSC_USE_MIXED_PRECISION)
1554: for (i=0; i<4; i++) {
1555: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[i],norms[i],flgs[i]);
1556: }
1557: #endif
1559: PetscLogEventBegin(VEC_Copy,x,y,0,0);
1560: #if defined(PETSC_USE_MIXED_PRECISION)
1561: extern PetscErrorCode VecGetArray(Vec,double**);
1562: extern PetscErrorCode VecRestoreArray(Vec,double**);
1563: extern PetscErrorCode VecGetArray(Vec,float**);
1564: extern PetscErrorCode VecRestoreArray(Vec,float**);
1565: extern PetscErrorCode VecGetArrayRead(Vec,const double**);
1566: extern PetscErrorCode VecRestoreArrayRead(Vec,const double**);
1567: extern PetscErrorCode VecGetArrayRead(Vec,const float**);
1568: extern PetscErrorCode VecRestoreArrayRead(Vec,const float**);
1569: if ((((PetscObject)x)->precision == PETSC_PRECISION_SINGLE) && (((PetscObject)y)->precision == PETSC_PRECISION_DOUBLE)) {
1570: PetscInt i,n;
1571: const float *xx;
1572: double *yy;
1573: VecGetArrayRead(x,&xx);
1574: VecGetArray(y,&yy);
1575: VecGetLocalSize(x,&n);
1576: for (i=0; i<n; i++) yy[i] = xx[i];
1577: VecRestoreArrayRead(x,&xx);
1578: VecRestoreArray(y,&yy);
1579: } else if ((((PetscObject)x)->precision == PETSC_PRECISION_DOUBLE) && (((PetscObject)y)->precision == PETSC_PRECISION_SINGLE)) {
1580: PetscInt i,n;
1581: float *yy;
1582: const double *xx;
1583: VecGetArrayRead(x,&xx);
1584: VecGetArray(y,&yy);
1585: VecGetLocalSize(x,&n);
1586: for (i=0; i<n; i++) yy[i] = (float) xx[i];
1587: VecRestoreArrayRead(x,&xx);
1588: VecRestoreArray(y,&yy);
1589: } else {
1590: (*x->ops->copy)(x,y);
1591: }
1592: #else
1593: (*x->ops->copy)(x,y);
1594: #endif
1596: PetscObjectStateIncrease((PetscObject)y);
1597: #if !defined(PETSC_USE_MIXED_PRECISION)
1598: for (i=0; i<4; i++) {
1599: if (flgs[i]) {
1600: PetscObjectComposedDataSetReal((PetscObject)y,NormIds[i],norms[i]);
1601: }
1602: }
1603: #endif
1605: PetscLogEventEnd(VEC_Copy,x,y,0,0);
1606: return(0);
1607: }
1609: /*@
1610: VecSwap - Swaps the vectors x and y.
1612: Logically Collective on Vec
1614: Input Parameters:
1615: . x, y - the vectors
1617: Level: advanced
1619: @*/
1620: PetscErrorCode VecSwap(Vec x,Vec y)
1621: {
1622: PetscReal normxs[4]={0.0,0.0,0.0,0.0},normys[4]={0.0,0.0,0.0,0.0};
1623: PetscBool flgxs[4],flgys[4];
1625: PetscInt i;
1633: VecCheckSameSize(x,1,y,2);
1634: if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1635: if (y->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1636: VecSetErrorIfLocked(x,1);
1637: VecSetErrorIfLocked(y,2);
1639: PetscLogEventBegin(VEC_Swap,x,y,0,0);
1640: for (i=0; i<4; i++) {
1641: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[i],normxs[i],flgxs[i]);
1642: PetscObjectComposedDataGetReal((PetscObject)y,NormIds[i],normys[i],flgys[i]);
1643: }
1644: (*x->ops->swap)(x,y);
1645: PetscObjectStateIncrease((PetscObject)x);
1646: PetscObjectStateIncrease((PetscObject)y);
1647: for (i=0; i<4; i++) {
1648: if (flgxs[i]) {
1649: PetscObjectComposedDataSetReal((PetscObject)y,NormIds[i],normxs[i]);
1650: }
1651: if (flgys[i]) {
1652: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[i],normys[i]);
1653: }
1654: }
1655: PetscLogEventEnd(VEC_Swap,x,y,0,0);
1656: return(0);
1657: }
1659: /*
1660: VecStashViewFromOptions - Processes command line options to determine if/how an VecStash object is to be viewed.
1662: Collective on VecStash
1664: Input Parameters:
1665: + obj - the VecStash object
1666: . bobj - optional other object that provides the prefix
1667: - optionname - option to activate viewing
1669: Level: intermediate
1671: Developer Note: This cannot use PetscObjectViewFromOptions() because it takes a Vec as an argument but does not use VecView
1673: */
1674: PetscErrorCode VecStashViewFromOptions(Vec obj,PetscObject bobj,const char optionname[])
1675: {
1676: PetscErrorCode ierr;
1677: PetscViewer viewer;
1678: PetscBool flg;
1679: PetscViewerFormat format;
1680: char *prefix;
1683: prefix = bobj ? bobj->prefix : ((PetscObject)obj)->prefix;
1684: PetscOptionsGetViewer(PetscObjectComm((PetscObject)obj),((PetscObject)obj)->options,prefix,optionname,&viewer,&format,&flg);
1685: if (flg) {
1686: PetscViewerPushFormat(viewer,format);
1687: VecStashView(obj,viewer);
1688: PetscViewerPopFormat(viewer);
1689: PetscViewerDestroy(&viewer);
1690: }
1691: return(0);
1692: }
1694: /*@
1695: VecStashView - Prints the entries in the vector stash and block stash.
1697: Collective on Vec
1699: Input Parameters:
1700: + v - the vector
1701: - viewer - the viewer
1703: Level: advanced
1705: .seealso: VecSetBlockSize(), VecSetValues(), VecSetValuesBlocked()
1707: @*/
1708: PetscErrorCode VecStashView(Vec v,PetscViewer viewer)
1709: {
1711: PetscMPIInt rank;
1712: PetscInt i,j;
1713: PetscBool match;
1714: VecStash *s;
1715: PetscScalar val;
1722: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&match);
1723: if (!match) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_SUP,"Stash viewer only works with ASCII viewer not %s\n",((PetscObject)v)->type_name);
1724: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
1725: MPI_Comm_rank(PetscObjectComm((PetscObject)v),&rank);
1726: s = &v->bstash;
1728: /* print block stash */
1729: PetscViewerASCIIPushSynchronized(viewer);
1730: PetscViewerASCIISynchronizedPrintf(viewer,"[%d]Vector Block stash size %D block size %D\n",rank,s->n,s->bs);
1731: for (i=0; i<s->n; i++) {
1732: PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Element %D ",rank,s->idx[i]);
1733: for (j=0; j<s->bs; j++) {
1734: val = s->array[i*s->bs+j];
1735: #if defined(PETSC_USE_COMPLEX)
1736: PetscViewerASCIISynchronizedPrintf(viewer,"(%18.16e %18.16e) ",PetscRealPart(val),PetscImaginaryPart(val));
1737: #else
1738: PetscViewerASCIISynchronizedPrintf(viewer,"%18.16e ",val);
1739: #endif
1740: }
1741: PetscViewerASCIISynchronizedPrintf(viewer,"\n");
1742: }
1743: PetscViewerFlush(viewer);
1745: s = &v->stash;
1747: /* print basic stash */
1748: PetscViewerASCIISynchronizedPrintf(viewer,"[%d]Vector stash size %D\n",rank,s->n);
1749: for (i=0; i<s->n; i++) {
1750: val = s->array[i];
1751: #if defined(PETSC_USE_COMPLEX)
1752: PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Element %D (%18.16e %18.16e) ",rank,s->idx[i],PetscRealPart(val),PetscImaginaryPart(val));
1753: #else
1754: PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Element %D %18.16e\n",rank,s->idx[i],val);
1755: #endif
1756: }
1757: PetscViewerFlush(viewer);
1758: PetscViewerASCIIPopSynchronized(viewer);
1759: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
1760: return(0);
1761: }
1763: PetscErrorCode PetscOptionsGetVec(PetscOptions options,const char prefix[],const char key[],Vec v,PetscBool *set)
1764: {
1765: PetscInt i,N,rstart,rend;
1767: PetscScalar *xx;
1768: PetscReal *xreal;
1769: PetscBool iset;
1772: VecGetOwnershipRange(v,&rstart,&rend);
1773: VecGetSize(v,&N);
1774: PetscCalloc1(N,&xreal);
1775: PetscOptionsGetRealArray(options,prefix,key,xreal,&N,&iset);
1776: if (iset) {
1777: VecGetArray(v,&xx);
1778: for (i=rstart; i<rend; i++) xx[i-rstart] = xreal[i];
1779: VecRestoreArray(v,&xx);
1780: }
1781: PetscFree(xreal);
1782: if (set) *set = iset;
1783: return(0);
1784: }
1786: /*@
1787: VecGetLayout - get PetscLayout describing vector layout
1789: Not Collective
1791: Input Parameter:
1792: . x - the vector
1794: Output Parameter:
1795: . map - the layout
1797: Level: developer
1799: .seealso: VecGetSizes(), VecGetOwnershipRange(), VecGetOwnershipRanges()
1800: @*/
1801: PetscErrorCode VecGetLayout(Vec x,PetscLayout *map)
1802: {
1806: *map = x->map;
1807: return(0);
1808: }
1810: /*@
1811: VecSetLayout - set PetscLayout describing vector layout
1813: Not Collective
1815: Input Parameters:
1816: + x - the vector
1817: - map - the layout
1819: Notes:
1820: It is normally only valid to replace the layout with a layout known to be equivalent.
1822: Level: developer
1824: .seealso: VecGetLayout(), VecGetSizes(), VecGetOwnershipRange(), VecGetOwnershipRanges()
1825: @*/
1826: PetscErrorCode VecSetLayout(Vec x,PetscLayout map)
1827: {
1832: PetscLayoutReference(map,&x->map);
1833: return(0);
1834: }
1836: PetscErrorCode VecSetInf(Vec xin)
1837: {
1838: PetscInt i,n = xin->map->n;
1839: PetscScalar *xx;
1840: PetscScalar zero=0.0,one=1.0,inf=one/zero;
1844: if (xin->ops->set) { /* can be called by a subset of processes, do not use collective routines */
1845: (*xin->ops->set)(xin,inf);
1846: } else {
1847: VecGetArrayWrite(xin,&xx);
1848: for (i=0; i<n; i++) xx[i] = inf;
1849: VecRestoreArrayWrite(xin,&xx);
1850: }
1851: return(0);
1852: }
1854: /*@
1855: VecBindToCPU - marks a vector to temporarily stay on the CPU and perform computations on the CPU
1857: Logically collective on Vec
1859: Input Parameters:
1860: + v - the vector
1861: - flg - bind to the CPU if value of PETSC_TRUE
1863: Level: intermediate
1864: @*/
1865: PetscErrorCode VecBindToCPU(Vec v,PetscBool flg)
1866: {
1870: #if defined(PETSC_HAVE_DEVICE)
1871: if (v->boundtocpu == flg) return(0);
1872: v->boundtocpu = flg;
1873: if (v->ops->bindtocpu) {
1875: (*v->ops->bindtocpu)(v,flg);
1876: }
1877: #endif
1878: return(0);
1879: }
1881: /*@
1882: VecBoundToCPU - query if a vector is bound to the CPU
1884: Not collective
1886: Input Parameter:
1887: . v - the vector
1889: Output Parameter:
1890: . flg - the logical flag
1892: Level: intermediate
1894: .seealso: VecBindToCPU()
1895: @*/
1896: PetscErrorCode VecBoundToCPU(Vec v,PetscBool *flg)
1897: {
1901: #if defined(PETSC_HAVE_DEVICE)
1902: *flg = v->boundtocpu;
1903: #else
1904: *flg = PETSC_TRUE;
1905: #endif
1906: return(0);
1907: }
1909: /*@C
1910: VecSetPinnedMemoryMin - Set the minimum data size for which pinned memory will be used for host (CPU) allocations.
1912: Logically Collective on Vec
1914: Input Parameters:
1915: + v - the vector
1916: - mbytes - minimum data size in bytes
1918: Options Database Keys:
1920: . -vec_pinned_memory_min <size> - minimum size (in bytes) for an allocation to use pinned memory on host.
1921: Note that this takes a PetscScalar, to accommodate large values;
1922: specifying -1 ensures that pinned memory will never be used.
1924: Level: developer
1926: .seealso: VecGetPinnedMemoryMin()
1927: @*/
1928: PetscErrorCode VecSetPinnedMemoryMin(Vec v,size_t mbytes)
1929: {
1930: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA) || defined(PETSC_HAVE_HIP)
1932: v->minimum_bytes_pinned_memory = mbytes;
1933: return(0);
1934: #else
1935: return 0;
1936: #endif
1937: }
1939: /*@C
1940: VecGetPinnedMemoryMin - Get the minimum data size for which pinned memory will be used for host (CPU) allocations.
1942: Logically Collective on Vec
1944: Input Parameters:
1945: . v - the vector
1947: Output Parameters:
1948: . mbytes - minimum data size in bytes
1950: Level: developer
1952: .seealso: VecSetPinnedMemoryMin()
1953: @*/
1954: PetscErrorCode VecGetPinnedMemoryMin(Vec v,size_t *mbytes)
1955: {
1956: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA) || defined(PETSC_HAVE_HIP)
1958: *mbytes = v->minimum_bytes_pinned_memory;
1959: return(0);
1960: #else
1961: return 0;
1962: #endif
1963: }
1965: /*@
1966: VecGetOffloadMask - Get the offload mask of a Vec.
1968: Not Collective
1970: Input Parameters:
1971: . v - the vector
1973: Output Parameters:
1974: . mask - corresponding PetscOffloadMask enum value.
1976: Level: intermediate
1978: .seealso: VecCreateSeqCUDA(), VecCreateSeqViennaCL(), VecGetArray(), VecGetType()
1979: @*/
1980: PetscErrorCode VecGetOffloadMask(Vec v,PetscOffloadMask* mask)
1981: {
1983: *mask = v->offloadmask;
1984: return(0);
1985: }
1987: #if !defined(PETSC_HAVE_VIENNACL)
1988: PETSC_EXTERN PetscErrorCode VecViennaCLGetCLContext(Vec v,PETSC_UINTPTR_T* ctx)
1989: {
1990: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"PETSc must be configured with --with-opencl to get a Vec's cl_context");
1991: }
1993: PETSC_EXTERN PetscErrorCode VecViennaCLGetCLQueue(Vec v,PETSC_UINTPTR_T* queue)
1994: {
1995: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"PETSc must be configured with --with-opencl to get a Vec's cl_command_queue");
1996: }
1998: PETSC_EXTERN PetscErrorCode VecViennaCLGetCLMem(Vec v,PETSC_UINTPTR_T* queue)
1999: {
2000: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"PETSc must be configured with --with-opencl to get a Vec's cl_mem");
2001: }
2003: PETSC_EXTERN PetscErrorCode VecViennaCLGetCLMemRead(Vec v,PETSC_UINTPTR_T* queue)
2004: {
2005: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"PETSc must be configured with --with-opencl to get a Vec's cl_mem");
2006: }
2008: PETSC_EXTERN PetscErrorCode VecViennaCLGetCLMemWrite(Vec v,PETSC_UINTPTR_T* queue)
2009: {
2010: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"PETSc must be configured with --with-opencl to get a Vec's cl_mem");
2011: }
2013: PETSC_EXTERN PetscErrorCode VecViennaCLRestoreCLMemWrite(Vec v)
2014: {
2015: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"PETSc must be configured with --with-opencl to restore a Vec's cl_mem");
2016: }
2017: #endif