Actual source code: vector.c
petsc-3.4.5 2014-06-29
2: /*
3: Provides the interface functions for vector operations that do NOT have PetscScalar/PetscReal in the signature
4: These are the vector functions the user calls.
5: */
6: #include <petsc-private/vecimpl.h> /*I "petscvec.h" I*/
8: /* Logging support */
9: PetscClassId VEC_CLASSID;
10: PetscLogEvent VEC_View, VEC_Max, VEC_Min, VEC_DotBarrier, VEC_Dot, VEC_MDotBarrier, 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_NormBarrier, VEC_MAXPY, VEC_Swap, VEC_AssemblyBegin, VEC_ScatterBegin, VEC_ScatterEnd;
13: PetscLogEvent VEC_AssemblyEnd, VEC_PointwiseMult, VEC_SetValues, VEC_Load, VEC_ScatterBarrier;
14: PetscLogEvent VEC_SetRandom, VEC_ReduceArithmetic, VEC_ReduceBarrier, VEC_ReduceCommunication,VEC_ReduceBegin,VEC_ReduceEnd,VEC_Ops;
15: PetscLogEvent VEC_DotNormBarrier, VEC_DotNorm, VEC_AXPBYPCZ, VEC_CUSPCopyFromGPU, VEC_CUSPCopyToGPU;
16: PetscLogEvent VEC_CUSPCopyFromGPUSome, VEC_CUSPCopyToGPUSome;
18: extern PetscErrorCode VecStashGetInfo_Private(VecStash*,PetscInt*,PetscInt*);
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: }
53: /*@
54: VecSetLocalToGlobalMapping - Sets a local numbering to global numbering used
55: by the routine VecSetValuesLocal() to allow users to insert vector entries
56: using a local (per-processor) numbering.
58: Logically Collective on Vec
60: Input Parameters:
61: + x - vector
62: - mapping - mapping created with ISLocalToGlobalMappingCreate() or ISLocalToGlobalMappingCreateIS()
64: Notes:
65: All vectors obtained with VecDuplicate() from this vector inherit the same mapping.
67: Level: intermediate
69: Concepts: vector^setting values with local numbering
71: seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetValuesLocal(),
72: VecSetLocalToGlobalMappingBlock(), VecSetValuesBlockedLocal()
73: @*/
74: PetscErrorCode VecSetLocalToGlobalMapping(Vec x,ISLocalToGlobalMapping mapping)
75: {
82: if (x->ops->setlocaltoglobalmapping) {
83: (*x->ops->setlocaltoglobalmapping)(x,mapping);
84: } else {
85: PetscLayoutSetISLocalToGlobalMapping(x->map,mapping);
86: }
87: return(0);
88: }
92: /*@
93: VecSetLocalToGlobalMappingBlock - Sets a local numbering to global numbering used
94: by the routine VecSetValuesBlockedLocal() to allow users to insert vector entries
95: using a local (per-processor) numbering.
97: Logically Collective on Vec
99: Input Parameters:
100: + x - vector
101: - mapping - mapping created with ISLocalToGlobalMappingCreate() or ISLocalToGlobalMappingCreateIS()
103: Notes:
104: All vectors obtained with VecDuplicate() from this vector inherit the same mapping.
106: Level: intermediate
108: Concepts: vector^setting values blocked with local numbering
110: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetValuesLocal(),
111: VecSetLocalToGlobalMapping(), VecSetValuesBlockedLocal()
112: @*/
113: PetscErrorCode VecSetLocalToGlobalMappingBlock(Vec x,ISLocalToGlobalMapping mapping)
114: {
121: PetscLayoutSetISLocalToGlobalMappingBlock(x->map,mapping);
122: return(0);
123: }
127: /*@
128: VecGetLocalToGlobalMapping - Gets the local-to-global numbering set by VecSetLocalToGlobalMapping()
130: Not Collective
132: Input Parameter:
133: . X - the vector
135: Output Parameter:
136: . mapping - the mapping
138: Level: advanced
140: Concepts: vectors^local to global mapping
141: Concepts: local to global mapping^for vectors
143: .seealso: VecSetValuesLocal(), VecGetLocalToGlobalMappingBlock()
144: @*/
145: PetscErrorCode VecGetLocalToGlobalMapping(Vec X,ISLocalToGlobalMapping *mapping)
146: {
151: *mapping = X->map->mapping;
152: return(0);
153: }
157: /*@
158: VecGetLocalToGlobalMappingBlock - Gets the local-to-global numbering set by VecSetLocalToGlobalMappingBlock()
160: Not Collective
162: Input Parameters:
163: . X - the vector
165: Output Parameters:
166: . mapping - the mapping
168: Level: advanced
170: Concepts: vectors^local to global mapping blocked
171: Concepts: local to global mapping^for vectors, blocked
173: .seealso: VecSetValuesBlockedLocal(), VecGetLocalToGlobalMapping()
174: @*/
175: PetscErrorCode VecGetLocalToGlobalMappingBlock(Vec X,ISLocalToGlobalMapping *mapping)
176: {
181: *mapping = X->map->bmapping;
182: return(0);
183: }
187: /*@
188: VecAssemblyBegin - Begins assembling the vector. This routine should
189: be called after completing all calls to VecSetValues().
191: Collective on Vec
193: Input Parameter:
194: . vec - the vector
196: Level: beginner
198: Concepts: assembly^vectors
200: .seealso: VecAssemblyEnd(), VecSetValues()
201: @*/
202: PetscErrorCode VecAssemblyBegin(Vec vec)
203: {
205: PetscBool flg = PETSC_FALSE;
211: PetscOptionsGetBool(((PetscObject)vec)->prefix,"-vec_view_stash",&flg,NULL);
212: if (flg) {
213: PetscViewer viewer;
214: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)vec),&viewer);
215: VecStashView(vec,viewer);
216: }
218: PetscLogEventBegin(VEC_AssemblyBegin,vec,0,0,0);
219: if (vec->ops->assemblybegin) {
220: (*vec->ops->assemblybegin)(vec);
221: }
222: PetscLogEventEnd(VEC_AssemblyBegin,vec,0,0,0);
223: PetscObjectStateIncrease((PetscObject)vec);
224: return(0);
225: }
229: /*
230: VecViewFromOptions - Processes command line options to determine if/how a vector is to be viewed. Called from higher level packages.
232: Collective on Vec
234: Input Parameters:
235: + vec - the vector
236: . prefix - prefix to use for viewing, or NULL to use prefix of 'vec'
237: - optionname - option to activate viewing
239: Level: intermediate
241: .keywords: Vec, view, options, database
242: .seealso: MatViewFromOptions()
243: */
244: PetscErrorCode VecViewFromOptions(Vec vec,const char prefix[],const char optionname[])
245: {
246: PetscErrorCode ierr;
247: PetscBool flg;
248: PetscViewer viewer;
249: PetscViewerFormat format;
252: if (prefix) {
253: PetscOptionsGetViewer(PetscObjectComm((PetscObject)vec),prefix,optionname,&viewer,&format,&flg);
254: } else {
255: PetscOptionsGetViewer(PetscObjectComm((PetscObject)vec),((PetscObject)vec)->prefix,optionname,&viewer,&format,&flg);
256: }
257: if (flg) {
258: PetscViewerPushFormat(viewer,format);
259: VecView(vec,viewer);
260: PetscViewerPopFormat(viewer);
261: PetscViewerDestroy(&viewer);
262: }
263: return(0);
264: }
268: /*@
269: VecAssemblyEnd - Completes assembling the vector. This routine should
270: be called after VecAssemblyBegin().
272: Collective on Vec
274: Input Parameter:
275: . vec - the vector
277: Options Database Keys:
278: + -vec_view - Prints vector in ASCII format
279: . -vec_view ::ascii_matlab - Prints vector in ASCII MATLAB format to stdout
280: . -vec_view matlab:filename - Prints vector in MATLAB format to matlaboutput.mat
281: . -vec_view draw - Activates vector viewing using drawing tools
282: . -display <name> - Sets display name (default is host)
283: . -draw_pause <sec> - Sets number of seconds to pause after display
284: - -vec_view socket - Activates vector viewing using a socket
286: Level: beginner
288: .seealso: VecAssemblyBegin(), VecSetValues()
289: @*/
290: PetscErrorCode VecAssemblyEnd(Vec vec)
291: {
296: PetscLogEventBegin(VEC_AssemblyEnd,vec,0,0,0);
298: if (vec->ops->assemblyend) {
299: (*vec->ops->assemblyend)(vec);
300: }
301: PetscLogEventEnd(VEC_AssemblyEnd,vec,0,0,0);
302: if (vec->viewonassembly) {
303: PetscViewerPushFormat(vec->viewonassembly,vec->viewformatonassembly);
304: VecView(vec,vec->viewonassembly);
305: PetscViewerPopFormat(vec->viewonassembly);
306: }
307: return(0);
308: }
312: /*@
313: VecPointwiseMax - Computes the componentwise maximum w_i = max(x_i, y_i).
315: Logically Collective on Vec
317: Input Parameters:
318: . x, y - the vectors
320: Output Parameter:
321: . w - the result
323: Level: advanced
325: Notes: any subset of the x, y, and w may be the same vector.
326: For complex numbers compares only the real part
328: Concepts: vector^pointwise multiply
330: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
331: @*/
332: PetscErrorCode VecPointwiseMax(Vec w,Vec x,Vec y)
333: {
345: if (x->map->N != y->map->N || x->map->N != w->map->N) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths");
346: if (x->map->n != y->map->n || x->map->n != w->map->n) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");
348: (*w->ops->pointwisemax)(w,x,y);
349: PetscObjectStateIncrease((PetscObject)w);
350: return(0);
351: }
356: /*@
357: VecPointwiseMin - Computes the componentwise minimum w_i = min(x_i, y_i).
359: Logically Collective on Vec
361: Input Parameters:
362: . x, y - the vectors
364: Output Parameter:
365: . w - the result
367: Level: advanced
369: Notes: any subset of the x, y, and w may be the same vector.
370: For complex numbers compares only the real part
372: Concepts: vector^pointwise multiply
374: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
375: @*/
376: PetscErrorCode VecPointwiseMin(Vec w,Vec x,Vec y)
377: {
389: if (x->map->N != y->map->N || x->map->N != w->map->N) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths");
390: if (x->map->n != y->map->n || x->map->n != w->map->n) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");
392: (*w->ops->pointwisemin)(w,x,y);
393: PetscObjectStateIncrease((PetscObject)w);
394: return(0);
395: }
399: /*@
400: VecPointwiseMaxAbs - Computes the componentwise maximum of the absolute values w_i = max(abs(x_i), abs(y_i)).
402: Logically Collective on Vec
404: Input Parameters:
405: . x, y - the vectors
407: Output Parameter:
408: . w - the result
410: Level: advanced
412: Notes: any subset of the x, y, and w may be the same vector.
414: Concepts: vector^pointwise multiply
416: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMin(), VecPointwiseMax(), VecMaxPointwiseDivide()
417: @*/
418: PetscErrorCode VecPointwiseMaxAbs(Vec w,Vec x,Vec y)
419: {
431: if (x->map->N != y->map->N || x->map->N != w->map->N) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths");
432: if (x->map->n != y->map->n || x->map->n != w->map->n) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");
434: (*w->ops->pointwisemaxabs)(w,x,y);
435: PetscObjectStateIncrease((PetscObject)w);
436: return(0);
437: }
441: /*@
442: VecPointwiseDivide - Computes the componentwise division w = x/y.
444: Logically Collective on Vec
446: Input Parameters:
447: . x, y - the vectors
449: Output Parameter:
450: . w - the result
452: Level: advanced
454: Notes: any subset of the x, y, and w may be the same vector.
456: Concepts: vector^pointwise divide
458: .seealso: VecPointwiseMult(), VecPointwiseMax(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
459: @*/
460: PetscErrorCode VecPointwiseDivide(Vec w,Vec x,Vec y)
461: {
473: if (x->map->N != y->map->N || x->map->N != w->map->N) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths");
474: if (x->map->n != y->map->n || x->map->n != w->map->n) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");
476: (*w->ops->pointwisedivide)(w,x,y);
477: PetscObjectStateIncrease((PetscObject)w);
478: return(0);
479: }
484: /*@
485: VecDuplicate - Creates a new vector of the same type as an existing vector.
487: Collective on Vec
489: Input Parameters:
490: . v - a vector to mimic
492: Output Parameter:
493: . newv - location to put new vector
495: Notes:
496: VecDuplicate() DOES NOT COPY the vector entries, but rather allocates storage
497: for the new vector. Use VecCopy() to copy a vector.
499: Use VecDestroy() to free the space. Use VecDuplicateVecs() to get several
500: vectors.
502: Level: beginner
504: .seealso: VecDestroy(), VecDuplicateVecs(), VecCreate(), VecCopy()
505: @*/
506: PetscErrorCode VecDuplicate(Vec v,Vec *newv)
507: {
514: (*v->ops->duplicate)(v,newv);
515: PetscObjectStateIncrease((PetscObject)*newv);
516: return(0);
517: }
521: /*@
522: VecDestroy - Destroys a vector.
524: Collective on Vec
526: Input Parameters:
527: . v - the vector
529: Level: beginner
531: .seealso: VecDuplicate(), VecDestroyVecs()
532: @*/
533: PetscErrorCode VecDestroy(Vec *v)
534: {
538: if (!*v) return(0);
540: if (--((PetscObject)(*v))->refct > 0) {*v = 0; return(0);}
542: PetscObjectAMSViewOff((PetscObject)*v);
543: PetscViewerDestroy(&(*v)->viewonassembly);
544: /* destroy the internal part */
545: if ((*v)->ops->destroy) {
546: (*(*v)->ops->destroy)(*v);
547: }
548: /* destroy the external/common part */
549: PetscLayoutDestroy(&(*v)->map);
550: PetscHeaderDestroy(v);
551: return(0);
552: }
556: /*@C
557: VecDuplicateVecs - Creates several vectors of the same type as an existing vector.
559: Collective on Vec
561: Input Parameters:
562: + m - the number of vectors to obtain
563: - v - a vector to mimic
565: Output Parameter:
566: . V - location to put pointer to array of vectors
568: Notes:
569: Use VecDestroyVecs() to free the space. Use VecDuplicate() to form a single
570: vector.
572: Fortran Note:
573: The Fortran interface is slightly different from that given below, it
574: requires one to pass in V a Vec (integer) array of size at least m.
575: See the Fortran chapter of the users manual and petsc/src/vec/vec/examples for details.
577: Level: intermediate
579: .seealso: VecDestroyVecs(), VecDuplicate(), VecCreate(), VecDuplicateVecsF90()
580: @*/
581: PetscErrorCode VecDuplicateVecs(Vec v,PetscInt m,Vec *V[])
582: {
589: (*v->ops->duplicatevecs)(v, m,V);
590: return(0);
591: }
595: /*@C
596: VecDestroyVecs - Frees a block of vectors obtained with VecDuplicateVecs().
598: Collective on Vec
600: Input Parameters:
601: + vv - pointer to pointer to array of vector pointers
602: - m - the number of vectors previously obtained
604: Fortran Note:
605: The Fortran interface is slightly different from that given below.
606: See the Fortran chapter of the users manual
608: Level: intermediate
610: .seealso: VecDuplicateVecs(), VecDestroyVecsf90()
611: @*/
612: PetscErrorCode VecDestroyVecs(PetscInt m,Vec *vv[])
613: {
618: if (!*vv) return(0);
621: if (m < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Trying to destroy negative number of vectors %D",m);
622: (*(**vv)->ops->destroyvecs)(m,*vv);
623: *vv = 0;
624: return(0);
625: }
629: /*@C
630: VecView - Views a vector object.
632: Collective on Vec
634: Input Parameters:
635: + vec - the vector
636: - viewer - an optional visualization context
638: Notes:
639: The available visualization contexts include
640: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
641: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
642: output where only the first processor opens
643: the file. All other processors send their
644: data to the first processor to print.
646: You can change the format the vector is printed using the
647: option PetscViewerSetFormat().
649: The user can open alternative visualization contexts with
650: + PetscViewerASCIIOpen() - Outputs vector to a specified file
651: . PetscViewerBinaryOpen() - Outputs vector in binary to a
652: specified file; corresponding input uses VecLoad()
653: . PetscViewerDrawOpen() - Outputs vector to an X window display
654: - PetscViewerSocketOpen() - Outputs vector to Socket viewer
656: The user can call PetscViewerSetFormat() to specify the output
657: format of ASCII printed objects (when using PETSC_VIEWER_STDOUT_SELF,
658: PETSC_VIEWER_STDOUT_WORLD and PetscViewerASCIIOpen). Available formats include
659: + PETSC_VIEWER_DEFAULT - default, prints vector contents
660: . PETSC_VIEWER_ASCII_MATLAB - prints vector contents in MATLAB format
661: . PETSC_VIEWER_ASCII_INDEX - prints vector contents, including indices of vector elements
662: - PETSC_VIEWER_ASCII_COMMON - prints vector contents, using a
663: format common among all vector types
665: Notes for HDF5 Viewer: the name of the Vec (given with PetscObjectSetName() is the name that is used
666: for the object in the HDF5 file. If you wish to store the same vector to the HDF5 viewer (with different values,
667: obviously) several times, you must change its name each time before calling the VecView(). The name you use
668: here should equal the name that you use in the Vec object that you use with VecLoad().
670: See the manual page for VecLoad() on the exact format the binary viewer stores
671: the values in the file.
673: Level: beginner
675: Concepts: vector^printing
676: Concepts: vector^saving to disk
678: .seealso: PetscViewerASCIIOpen(), PetscViewerDrawOpen(), PetscDrawLGCreate(),
679: PetscViewerSocketOpen(), PetscViewerBinaryOpen(), VecLoad(), PetscViewerCreate(),
680: PetscRealView(), PetscScalarView(), PetscIntView()
681: @*/
682: PetscErrorCode VecView(Vec vec,PetscViewer viewer)
683: {
684: PetscErrorCode ierr;
685: PetscBool iascii;
686: PetscViewerFormat format;
691: if (!viewer) {
692: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)vec),&viewer);
693: }
696: if (vec->stash.n || vec->bstash.n) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call VecAssemblyBegin/End() before viewing this vector");
698: PetscLogEventBegin(VEC_View,vec,viewer,0,0);
699: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
700: if (iascii) {
701: PetscInt rows,bs;
703: PetscViewerGetFormat(viewer,&format);
704: if (format == PETSC_VIEWER_ASCII_INFO || format == PETSC_VIEWER_ASCII_INFO_DETAIL) {
705: PetscObjectPrintClassNamePrefixType((PetscObject)vec,viewer,"Vector Object");
706: PetscViewerASCIIPushTab(viewer);
707: VecGetSize(vec,&rows);
708: VecGetBlockSize(vec,&bs);
709: if (bs != 1) {
710: PetscViewerASCIIPrintf(viewer,"length=%D, bs=%D\n",rows,bs);
711: } else {
712: PetscViewerASCIIPrintf(viewer,"length=%D\n",rows);
713: }
714: PetscViewerASCIIPopTab(viewer);
715: }
716: }
717: (*vec->ops->view)(vec,viewer);
718: PetscLogEventEnd(VEC_View,vec,viewer,0,0);
719: return(0);
720: }
722: #if defined(PETSC_USE_DEBUG)
723: #include <../src/sys/totalview/tv_data_display.h>
724: PETSC_UNUSED static int TV_display_type(const struct _p_Vec *v)
725: {
726: const PetscScalar *values;
727: char type[32];
728: PetscErrorCode ierr;
731: TV_add_row("Local rows", "int", &v->map->n);
732: TV_add_row("Global rows", "int", &v->map->N);
733: TV_add_row("Typename", TV_ascii_string_type, ((PetscObject)v)->type_name);
734: VecGetArrayRead((Vec)v,&values);
735: PetscSNPrintf(type,32,"double[%d]",v->map->n);
736: TV_add_row("values",type, values);
737: VecRestoreArrayRead((Vec)v,&values);
738: return TV_format_OK;
739: }
740: #endif
744: /*@
745: VecGetSize - Returns the global number of elements of the vector.
747: Not Collective
749: Input Parameter:
750: . x - the vector
752: Output Parameters:
753: . size - the global length of the vector
755: Level: beginner
757: Concepts: vector^local size
759: .seealso: VecGetLocalSize()
760: @*/
761: PetscErrorCode VecGetSize(Vec x,PetscInt *size)
762: {
769: (*x->ops->getsize)(x,size);
770: return(0);
771: }
775: /*@
776: VecGetLocalSize - Returns the number of elements of the vector stored
777: in local memory. This routine may be implementation dependent, so use
778: with care.
780: Not Collective
782: Input Parameter:
783: . x - the vector
785: Output Parameter:
786: . size - the length of the local piece of the vector
788: Level: beginner
790: Concepts: vector^size
792: .seealso: VecGetSize()
793: @*/
794: PetscErrorCode VecGetLocalSize(Vec x,PetscInt *size)
795: {
802: (*x->ops->getlocalsize)(x,size);
803: return(0);
804: }
808: /*@C
809: VecGetOwnershipRange - Returns the range of indices owned by
810: this processor, assuming that the vectors are laid out with the
811: first n1 elements on the first processor, next n2 elements on the
812: second, etc. For certain parallel layouts this range may not be
813: well defined.
815: Not Collective
817: Input Parameter:
818: . x - the vector
820: Output Parameters:
821: + low - the first local element, pass in NULL if not interested
822: - high - one more than the last local element, pass in NULL if not interested
824: Note:
825: The high argument is one more than the last element stored locally.
827: Fortran: NULL_INTEGER should be used instead of NULL
829: Level: beginner
831: Concepts: ownership^of vectors
832: Concepts: vector^ownership of elements
834: .seealso: MatGetOwnershipRange(), MatGetOwnershipRanges(), VecGetOwnershipRanges()
835: @*/
836: PetscErrorCode VecGetOwnershipRange(Vec x,PetscInt *low,PetscInt *high)
837: {
843: if (low) *low = x->map->rstart;
844: if (high) *high = x->map->rend;
845: return(0);
846: }
850: /*@C
851: VecGetOwnershipRanges - Returns the range of indices owned by EACH processor,
852: assuming that the vectors are laid out with the
853: first n1 elements on the first processor, next n2 elements on the
854: second, etc. For certain parallel layouts this range may not be
855: well defined.
857: Not Collective
859: Input Parameter:
860: . x - the vector
862: Output Parameters:
863: . range - array of length size+1 with the start and end+1 for each process
865: Note:
866: The high argument is one more than the last element stored locally.
868: Fortran: You must PASS in an array of length size+1
870: Level: beginner
872: Concepts: ownership^of vectors
873: Concepts: vector^ownership of elements
875: .seealso: MatGetOwnershipRange(), MatGetOwnershipRanges(), VecGetOwnershipRange()
876: @*/
877: PetscErrorCode VecGetOwnershipRanges(Vec x,const PetscInt *ranges[])
878: {
884: PetscLayoutGetRanges(x->map,ranges);
885: return(0);
886: }
890: /*@
891: VecSetOption - Sets an option for controling a vector's behavior.
893: Collective on Vec
895: Input Parameter:
896: + x - the vector
897: . op - the option
898: - flag - turn the option on or off
900: Supported Options:
901: + VEC_IGNORE_OFF_PROC_ENTRIES, which causes VecSetValues() to ignore
902: entries destined to be stored on a separate processor. This can be used
903: to eliminate the global reduction in the VecAssemblyXXXX() if you know
904: that you have only used VecSetValues() to set local elements
905: . VEC_IGNORE_NEGATIVE_INDICES, which means you can pass negative indices
906: in ix in calls to VecSetValues() or VecGetValues(). These rows are simply
907: ignored.
909: Level: intermediate
911: @*/
912: PetscErrorCode VecSetOption(Vec x,VecOption op,PetscBool flag)
913: {
919: if (x->ops->setoption) {
920: (*x->ops->setoption)(x,op,flag);
921: }
922: return(0);
923: }
927: /* Default routines for obtaining and releasing; */
928: /* may be used by any implementation */
929: PetscErrorCode VecDuplicateVecs_Default(Vec w,PetscInt m,Vec *V[])
930: {
932: PetscInt i;
937: if (m <= 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"m must be > 0: m = %D",m);
938: PetscMalloc(m*sizeof(Vec*),V);
939: for (i=0; i<m; i++) {VecDuplicate(w,*V+i);}
940: return(0);
941: }
945: PetscErrorCode VecDestroyVecs_Default(PetscInt m,Vec v[])
946: {
948: PetscInt i;
952: if (m <= 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"m must be > 0: m = %D",m);
953: for (i=0; i<m; i++) {VecDestroy(&v[i]);}
954: PetscFree(v);
955: return(0);
956: }
960: /*@
961: VecResetArray - Resets a vector to use its default memory. Call this
962: after the use of VecPlaceArray().
964: Not Collective
966: Input Parameters:
967: . vec - the vector
969: Level: developer
971: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecPlaceArray()
973: @*/
974: PetscErrorCode VecResetArray(Vec vec)
975: {
981: if (vec->ops->resetarray) {
982: (*vec->ops->resetarray)(vec);
983: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Cannot reset array in this type of vector");
984: PetscObjectStateIncrease((PetscObject)vec);
985: return(0);
986: }
990: /*@C
991: VecLoad - Loads a vector that has been stored in binary or HDF5 format
992: with VecView().
994: Collective on PetscViewer
996: Input Parameters:
997: + newvec - the newly loaded vector, this needs to have been created with VecCreate() or
998: some related function before a call to VecLoad().
999: - viewer - binary file viewer, obtained from PetscViewerBinaryOpen() or
1000: HDF5 file viewer, obtained from PetscViewerHDF5Open()
1002: Level: intermediate
1004: Notes:
1005: Defaults to the standard Seq or MPI Vec, if you want some other type of Vec call VecSetFromOptions()
1006: before calling this.
1008: The input file must contain the full global vector, as
1009: written by the routine VecView().
1011: If the type or size of newvec is not set before a call to VecLoad, PETSc
1012: sets the type and the local and global sizes.If type and/or
1013: sizes are already set, then the same are used.
1015: IF using HDF5, you must assign the Vec the same name as was used in the Vec
1016: that was stored in the file using PetscObjectSetName(). Otherwise you will
1017: get the error message: "Cannot H5DOpen2() with Vec name NAMEOFOBJECT"
1019: Notes for advanced users:
1020: Most users should not need to know the details of the binary storage
1021: format, since VecLoad() and VecView() completely hide these details.
1022: But for anyone who's interested, the standard binary matrix storage
1023: format is
1024: .vb
1025: int VEC_FILE_CLASSID
1026: int number of rows
1027: PetscScalar *values of all entries
1028: .ve
1030: In addition, PETSc automatically does the byte swapping for
1031: machines that store the bytes reversed, e.g. DEC alpha, freebsd,
1032: linux, Windows and the paragon; thus if you write your own binary
1033: read/write routines you have to swap the bytes; see PetscBinaryRead()
1034: and PetscBinaryWrite() to see how this may be done.
1036: Concepts: vector^loading from file
1038: .seealso: PetscViewerBinaryOpen(), VecView(), MatLoad(), VecLoad()
1039: @*/
1040: PetscErrorCode VecLoad(Vec newvec, PetscViewer viewer)
1041: {
1043: PetscBool isbinary,ishdf5;
1048: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
1049: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5);
1050: if (!isbinary && !ishdf5) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Invalid viewer; open viewer with PetscViewerBinaryOpen()");
1052: PetscLogEventBegin(VEC_Load,viewer,0,0,0);
1053: if (!((PetscObject)newvec)->type_name && !newvec->ops->create) {
1054: VecSetType(newvec, VECSTANDARD);
1055: }
1056: (*newvec->ops->load)(newvec,viewer);
1057: PetscLogEventEnd(VEC_Load,viewer,0,0,0);
1058: return(0);
1059: }
1064: /*@
1065: VecReciprocal - Replaces each component of a vector by its reciprocal.
1067: Logically Collective on Vec
1069: Input Parameter:
1070: . vec - the vector
1072: Output Parameter:
1073: . vec - the vector reciprocal
1075: Level: intermediate
1077: Concepts: vector^reciprocal
1079: .seealso: VecLog(), VecExp(), VecSqrtAbs()
1081: @*/
1082: PetscErrorCode VecReciprocal(Vec vec)
1083: {
1089: if (vec->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1090: if (!vec->ops->reciprocal) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Vector does not support reciprocal operation");
1091: (*vec->ops->reciprocal)(vec);
1092: PetscObjectStateIncrease((PetscObject)vec);
1093: return(0);
1094: }
1098: PetscErrorCode VecSetOperation(Vec vec,VecOperation op, void (*f)(void))
1099: {
1102: (((void(**)(void))vec->ops)[(int)op]) = f;
1103: return(0);
1104: }
1109: /*@
1110: VecStashSetInitialSize - sets the sizes of the vec-stash, that is
1111: used during the assembly process to store values that belong to
1112: other processors.
1114: Not Collective, different processes can have different size stashes
1116: Input Parameters:
1117: + vec - the vector
1118: . size - the initial size of the stash.
1119: - bsize - the initial size of the block-stash(if used).
1121: Options Database Keys:
1122: + -vecstash_initial_size <size> or <size0,size1,...sizep-1>
1123: - -vecstash_block_initial_size <bsize> or <bsize0,bsize1,...bsizep-1>
1125: Level: intermediate
1127: Notes:
1128: The block-stash is used for values set with VecSetValuesBlocked() while
1129: the stash is used for values set with VecSetValues()
1131: Run with the option -info and look for output of the form
1132: VecAssemblyBegin_MPIXXX:Stash has MM entries, uses nn mallocs.
1133: to determine the appropriate value, MM, to use for size and
1134: VecAssemblyBegin_MPIXXX:Block-Stash has BMM entries, uses nn mallocs.
1135: to determine the value, BMM to use for bsize
1137: Concepts: vector^stash
1138: Concepts: stash^vector
1140: .seealso: VecSetBlockSize(), VecSetValues(), VecSetValuesBlocked(), VecStashView()
1142: @*/
1143: PetscErrorCode VecStashSetInitialSize(Vec vec,PetscInt size,PetscInt bsize)
1144: {
1149: VecStashSetInitialSize_Private(&vec->stash,size);
1150: VecStashSetInitialSize_Private(&vec->bstash,bsize);
1151: return(0);
1152: }
1156: /*@
1157: VecConjugate - Conjugates a vector.
1159: Logically Collective on Vec
1161: Input Parameters:
1162: . x - the vector
1164: Level: intermediate
1166: Concepts: vector^conjugate
1168: @*/
1169: PetscErrorCode VecConjugate(Vec x)
1170: {
1171: #if defined(PETSC_USE_COMPLEX)
1177: if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1178: (*x->ops->conjugate)(x);
1179: /* we need to copy norms here */
1180: PetscObjectStateIncrease((PetscObject)x);
1181: return(0);
1182: #else
1183: return(0);
1184: #endif
1185: }
1189: /*@
1190: VecPointwiseMult - Computes the componentwise multiplication w = x*y.
1192: Logically Collective on Vec
1194: Input Parameters:
1195: . x, y - the vectors
1197: Output Parameter:
1198: . w - the result
1200: Level: advanced
1202: Notes: any subset of the x, y, and w may be the same vector.
1204: Concepts: vector^pointwise multiply
1206: .seealso: VecPointwiseDivide(), VecPointwiseMax(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
1207: @*/
1208: PetscErrorCode VecPointwiseMult(Vec w, Vec x,Vec y)
1209: {
1221: if (x->map->n != y->map->n || x->map->n != w->map->n) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");
1223: PetscLogEventBegin(VEC_PointwiseMult,x,y,w,0);
1224: (*w->ops->pointwisemult)(w,x,y);
1225: PetscLogEventEnd(VEC_PointwiseMult,x,y,w,0);
1226: PetscObjectStateIncrease((PetscObject)w);
1227: return(0);
1228: }
1232: /*@
1233: VecSetRandom - Sets all components of a vector to random numbers.
1235: Logically Collective on Vec
1237: Input Parameters:
1238: + x - the vector
1239: - rctx - the random number context, formed by PetscRandomCreate(), or NULL and
1240: it will create one internally.
1242: Output Parameter:
1243: . x - the vector
1245: Example of Usage:
1246: .vb
1247: PetscRandomCreate(PETSC_COMM_WORLD,&rctx);
1248: VecSetRandom(x,rctx);
1249: PetscRandomDestroy(rctx);
1250: .ve
1252: Level: intermediate
1254: Concepts: vector^setting to random
1255: Concepts: random^vector
1257: .seealso: VecSet(), VecSetValues(), PetscRandomCreate(), PetscRandomDestroy()
1258: @*/
1259: PetscErrorCode VecSetRandom(Vec x,PetscRandom rctx)
1260: {
1262: PetscRandom randObj = NULL;
1268: if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1270: if (!rctx) {
1271: MPI_Comm comm;
1272: PetscObjectGetComm((PetscObject)x,&comm);
1273: PetscRandomCreate(comm,&randObj);
1274: PetscRandomSetFromOptions(randObj);
1275: rctx = randObj;
1276: }
1278: PetscLogEventBegin(VEC_SetRandom,x,rctx,0,0);
1279: (*x->ops->setrandom)(x,rctx);
1280: PetscLogEventEnd(VEC_SetRandom,x,rctx,0,0);
1282: PetscRandomDestroy(&randObj);
1283: PetscObjectStateIncrease((PetscObject)x);
1284: return(0);
1285: }
1289: /*@
1290: VecZeroEntries - puts a 0.0 in each element of a vector
1292: Logically Collective on Vec
1294: Input Parameter:
1295: . vec - The vector
1297: Level: beginner
1299: Developer Note: This routine does not need to exist since the exact functionality is obtained with
1300: VecSet(vec,0); I guess someone added it to mirror the functionality of MatZeroEntries() but Mat is nothing
1301: like a Vec (one is an operator and one is an element of a vector space, yeah yeah dual blah blah blah) so
1302: this routine should not exist.
1304: .keywords: Vec, set, options, database
1305: .seealso: VecCreate(), VecSetOptionsPrefix(), VecSet(), VecSetValues()
1306: @*/
1307: PetscErrorCode VecZeroEntries(Vec vec)
1308: {
1312: VecSet(vec,0);
1313: return(0);
1314: }
1318: /*
1319: VecSetTypeFromOptions_Private - Sets the type of vector from user options. Defaults to a PETSc sequential vector on one
1320: processor and a PETSc MPI vector on more than one processor.
1322: Collective on Vec
1324: Input Parameter:
1325: . vec - The vector
1327: Level: intermediate
1329: .keywords: Vec, set, options, database, type
1330: .seealso: VecSetFromOptions(), VecSetType()
1331: */
1332: static PetscErrorCode VecSetTypeFromOptions_Private(Vec vec)
1333: {
1334: PetscBool opt;
1335: VecType defaultType;
1336: char typeName[256];
1337: PetscMPIInt size;
1341: if (((PetscObject)vec)->type_name) defaultType = ((PetscObject)vec)->type_name;
1342: else {
1343: MPI_Comm_size(PetscObjectComm((PetscObject)vec), &size);
1344: if (size > 1) defaultType = VECMPI;
1345: else defaultType = VECSEQ;
1346: }
1348: if (!VecRegisterAllCalled) {VecRegisterAll();}
1349: PetscOptionsList("-vec_type","Vector type","VecSetType",VecList,defaultType,typeName,256,&opt);
1350: if (opt) {
1351: VecSetType(vec, typeName);
1352: } else {
1353: VecSetType(vec, defaultType);
1354: }
1355: return(0);
1356: }
1360: /*@
1361: VecSetFromOptions - Configures the vector from the options database.
1363: Collective on Vec
1365: Input Parameter:
1366: . vec - The vector
1368: Notes: To see all options, run your program with the -help option, or consult the users manual.
1369: Must be called after VecCreate() but before the vector is used.
1371: Level: beginner
1373: Concepts: vectors^setting options
1374: Concepts: vectors^setting type
1376: .keywords: Vec, set, options, database
1377: .seealso: VecCreate(), VecSetOptionsPrefix()
1378: @*/
1379: PetscErrorCode VecSetFromOptions(Vec vec)
1380: {
1386: PetscObjectOptionsBegin((PetscObject)vec);
1387: /* Handle vector type options */
1388: VecSetTypeFromOptions_Private(vec);
1389: PetscViewerDestroy(&vec->viewonassembly);
1390: PetscOptionsViewer("-vec_view","Display vector with the viewer on VecAssemblyEnd()","VecView",&vec->viewonassembly,&vec->viewformatonassembly,NULL);
1392: /* Handle specific vector options */
1393: if (vec->ops->setfromoptions) {
1394: (*vec->ops->setfromoptions)(vec);
1395: }
1397: /* process any options handlers added with PetscObjectAddOptionsHandler() */
1398: PetscObjectProcessOptionsHandlers((PetscObject)vec);
1399: PetscOptionsEnd();
1400: return(0);
1401: }
1405: /*@
1406: VecSetSizes - Sets the local and global sizes, and checks to determine compatibility
1408: Collective on Vec
1410: Input Parameters:
1411: + v - the vector
1412: . n - the local size (or PETSC_DECIDE to have it set)
1413: - N - the global size (or PETSC_DECIDE)
1415: Notes:
1416: n and N cannot be both PETSC_DECIDE
1417: If one processor calls this with N of PETSC_DECIDE then all processors must, otherwise the program will hang.
1419: Level: intermediate
1421: .seealso: VecGetSize(), PetscSplitOwnership()
1422: @*/
1423: PetscErrorCode VecSetSizes(Vec v, PetscInt n, PetscInt N)
1424: {
1430: 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);
1431: 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);
1432: v->map->n = n;
1433: v->map->N = N;
1434: if (v->ops->create) {
1435: (*v->ops->create)(v);
1436: v->ops->create = 0;
1437: }
1438: return(0);
1439: }
1443: /*@
1444: VecSetBlockSize - Sets the blocksize for future calls to VecSetValuesBlocked()
1445: and VecSetValuesBlockedLocal().
1447: Logically Collective on Vec
1449: Input Parameter:
1450: + v - the vector
1451: - bs - the blocksize
1453: Notes:
1454: All vectors obtained by VecDuplicate() inherit the same blocksize.
1456: Level: advanced
1458: .seealso: VecSetValuesBlocked(), VecSetLocalToGlobalMappingBlock(), VecGetBlockSize()
1460: Concepts: block size^vectors
1461: @*/
1462: PetscErrorCode VecSetBlockSize(Vec v,PetscInt bs)
1463: {
1468: if (bs == v->map->bs) return(0);
1470: PetscLayoutSetBlockSize(v->map,bs);
1471: v->bstash.bs = bs; /* use the same blocksize for the vec's block-stash */
1472: return(0);
1473: }
1477: /*@
1478: VecGetBlockSize - Gets the blocksize for the vector, i.e. what is used for VecSetValuesBlocked()
1479: and VecSetValuesBlockedLocal().
1481: Not Collective
1483: Input Parameter:
1484: . v - the vector
1486: Output Parameter:
1487: . bs - the blocksize
1489: Notes:
1490: All vectors obtained by VecDuplicate() inherit the same blocksize.
1492: Level: advanced
1494: .seealso: VecSetValuesBlocked(), VecSetLocalToGlobalMappingBlock(), VecSetBlockSize()
1496: Concepts: vector^block size
1497: Concepts: block^vector
1499: @*/
1500: PetscErrorCode VecGetBlockSize(Vec v,PetscInt *bs)
1501: {
1507: PetscLayoutGetBlockSize(v->map,bs);
1508: return(0);
1509: }
1513: /*@C
1514: VecSetOptionsPrefix - Sets the prefix used for searching for all
1515: Vec options in the database.
1517: Logically Collective on Vec
1519: Input Parameter:
1520: + v - the Vec context
1521: - prefix - the prefix to prepend to all option names
1523: Notes:
1524: A hyphen (-) must NOT be given at the beginning of the prefix name.
1525: The first character of all runtime options is AUTOMATICALLY the hyphen.
1527: Level: advanced
1529: .keywords: Vec, set, options, prefix, database
1531: .seealso: VecSetFromOptions()
1532: @*/
1533: PetscErrorCode VecSetOptionsPrefix(Vec v,const char prefix[])
1534: {
1539: PetscObjectSetOptionsPrefix((PetscObject)v,prefix);
1540: return(0);
1541: }
1545: /*@C
1546: VecAppendOptionsPrefix - Appends to the prefix used for searching for all
1547: Vec options in the database.
1549: Logically Collective on Vec
1551: Input Parameters:
1552: + v - the Vec context
1553: - prefix - the prefix to prepend to all option names
1555: Notes:
1556: A hyphen (-) must NOT be given at the beginning of the prefix name.
1557: The first character of all runtime options is AUTOMATICALLY the hyphen.
1559: Level: advanced
1561: .keywords: Vec, append, options, prefix, database
1563: .seealso: VecGetOptionsPrefix()
1564: @*/
1565: PetscErrorCode VecAppendOptionsPrefix(Vec v,const char prefix[])
1566: {
1571: PetscObjectAppendOptionsPrefix((PetscObject)v,prefix);
1572: return(0);
1573: }
1577: /*@C
1578: VecGetOptionsPrefix - Sets the prefix used for searching for all
1579: Vec options in the database.
1581: Not Collective
1583: Input Parameter:
1584: . v - the Vec context
1586: Output Parameter:
1587: . prefix - pointer to the prefix string used
1589: Notes: On the fortran side, the user should pass in a string 'prefix' of
1590: sufficient length to hold the prefix.
1592: Level: advanced
1594: .keywords: Vec, get, options, prefix, database
1596: .seealso: VecAppendOptionsPrefix()
1597: @*/
1598: PetscErrorCode VecGetOptionsPrefix(Vec v,const char *prefix[])
1599: {
1604: PetscObjectGetOptionsPrefix((PetscObject)v,prefix);
1605: return(0);
1606: }
1610: /*@
1611: VecSetUp - Sets up the internal vector data structures for the later use.
1613: Collective on Vec
1615: Input Parameters:
1616: . v - the Vec context
1618: Notes:
1619: For basic use of the Vec classes the user need not explicitly call
1620: VecSetUp(), since these actions will happen automatically.
1622: Level: advanced
1624: .keywords: Vec, setup
1626: .seealso: VecCreate(), VecDestroy()
1627: @*/
1628: PetscErrorCode VecSetUp(Vec v)
1629: {
1630: PetscMPIInt size;
1635: if (!((PetscObject)v)->type_name) {
1636: MPI_Comm_size(PetscObjectComm((PetscObject)v), &size);
1637: if (size == 1) {
1638: VecSetType(v, VECSEQ);
1639: } else {
1640: VecSetType(v, VECMPI);
1641: }
1642: }
1643: return(0);
1644: }
1646: /*
1647: These currently expose the PetscScalar/PetscReal in updating the
1648: cached norm. If we push those down into the implementation these
1649: will become independent of PetscScalar/PetscReal
1650: */
1654: /*@
1655: VecCopy - Copies a vector. y <- x
1657: Logically Collective on Vec
1659: Input Parameter:
1660: . x - the vector
1662: Output Parameter:
1663: . y - the copy
1665: Notes:
1666: For default parallel PETSc vectors, both x and y must be distributed in
1667: the same manner; local copies are done.
1669: Level: beginner
1671: .seealso: VecDuplicate()
1672: @*/
1673: PetscErrorCode VecCopy(Vec x,Vec y)
1674: {
1675: PetscBool flgs[4];
1676: PetscReal norms[4] = {0.0,0.0,0.0,0.0};
1678: PetscInt i;
1685: if (x == y) return(0);
1686: if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1687: if (x->map->n != y->map->n) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths %d != %d", x->map->n, y->map->n);
1689: #if !defined(PETSC_USE_MIXED_PRECISION)
1690: for (i=0; i<4; i++) {
1691: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[i],norms[i],flgs[i]);
1692: }
1693: #endif
1695: PetscLogEventBegin(VEC_Copy,x,y,0,0);
1696: #if defined(PETSC_USE_MIXED_PRECISION)
1697: extern PetscErrorCode VecGetArray(Vec,double**);
1698: extern PetscErrorCode VecRestoreArray(Vec,double**);
1699: extern PetscErrorCode VecGetArray(Vec,float**);
1700: extern PetscErrorCode VecRestoreArray(Vec,float**);
1701: extern PetscErrorCode VecGetArrayRead(Vec,const double**);
1702: extern PetscErrorCode VecRestoreArrayRead(Vec,const double**);
1703: extern PetscErrorCode VecGetArrayRead(Vec,const float**);
1704: extern PetscErrorCode VecRestoreArrayRead(Vec,const float**);
1705: if ((((PetscObject)x)->precision == PETSC_PRECISION_SINGLE) && (((PetscObject)y)->precision == PETSC_PRECISION_DOUBLE)) {
1706: PetscInt i,n;
1707: const float *xx;
1708: double *yy;
1709: VecGetArrayRead(x,&xx);
1710: VecGetArray(y,&yy);
1711: VecGetLocalSize(x,&n);
1712: for (i=0; i<n; i++) yy[i] = xx[i];
1713: VecRestoreArrayRead(x,&xx);
1714: VecRestoreArray(y,&yy);
1715: } else if ((((PetscObject)x)->precision == PETSC_PRECISION_DOUBLE) && (((PetscObject)y)->precision == PETSC_PRECISION_SINGLE)) {
1716: PetscInt i,n;
1717: float *yy;
1718: const double *xx;
1719: VecGetArrayRead(x,&xx);
1720: VecGetArray(y,&yy);
1721: VecGetLocalSize(x,&n);
1722: for (i=0; i<n; i++) yy[i] = (float) xx[i];
1723: VecRestoreArrayRead(x,&xx);
1724: VecRestoreArray(y,&yy);
1725: } else {
1726: (*x->ops->copy)(x,y);
1727: }
1728: #else
1729: (*x->ops->copy)(x,y);
1730: #endif
1732: PetscObjectStateIncrease((PetscObject)y);
1733: #if !defined(PETSC_USE_MIXED_PRECISION)
1734: for (i=0; i<4; i++) {
1735: if (flgs[i]) {
1736: PetscObjectComposedDataSetReal((PetscObject)y,NormIds[i],norms[i]);
1737: }
1738: }
1739: #endif
1741: PetscLogEventEnd(VEC_Copy,x,y,0,0);
1742: return(0);
1743: }
1747: /*@
1748: VecSwap - Swaps the vectors x and y.
1750: Logically Collective on Vec
1752: Input Parameters:
1753: . x, y - the vectors
1755: Level: advanced
1757: Concepts: vector^swapping values
1759: @*/
1760: PetscErrorCode VecSwap(Vec x,Vec y)
1761: {
1762: PetscReal normxs[4]={0.0,0.0,0.0,0.0},normys[4]={0.0,0.0,0.0,0.0};
1763: PetscBool flgxs[4],flgys[4];
1765: PetscInt i;
1773: if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1774: if (y->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1775: if (x->map->N != y->map->N) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths");
1776: if (x->map->n != y->map->n) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");
1778: PetscLogEventBegin(VEC_Swap,x,y,0,0);
1779: for (i=0; i<4; i++) {
1780: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[i],normxs[i],flgxs[i]);
1781: PetscObjectComposedDataGetReal((PetscObject)y,NormIds[i],normys[i],flgys[i]);
1782: }
1783: (*x->ops->swap)(x,y);
1784: PetscObjectStateIncrease((PetscObject)x);
1785: PetscObjectStateIncrease((PetscObject)y);
1786: for (i=0; i<4; i++) {
1787: if (flgxs[i]) {
1788: PetscObjectComposedDataSetReal((PetscObject)y,NormIds[i],normxs[i]);
1789: }
1790: if (flgys[i]) {
1791: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[i],normys[i]);
1792: }
1793: }
1794: PetscLogEventEnd(VEC_Swap,x,y,0,0);
1795: return(0);
1796: }
1800: /*@
1801: VecStashView - Prints the entries in the vector stash and block stash.
1803: Collective on Vec
1805: Input Parameters:
1806: + v - the vector
1807: - viewer - the viewer
1809: Level: advanced
1811: Concepts: vector^stash
1812: Concepts: stash^vector
1814: .seealso: VecSetBlockSize(), VecSetValues(), VecSetValuesBlocked()
1816: @*/
1817: PetscErrorCode VecStashView(Vec v,PetscViewer viewer)
1818: {
1820: PetscMPIInt rank;
1821: PetscInt i,j;
1822: PetscBool match;
1823: VecStash *s;
1824: PetscScalar val;
1831: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&match);
1832: if (!match) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_SUP,"Stash viewer only works with ASCII viewer not %s\n",((PetscObject)v)->type_name);
1833: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
1834: MPI_Comm_rank(PetscObjectComm((PetscObject)v),&rank);
1835: s = &v->bstash;
1837: /* print block stash */
1838: PetscViewerASCIISynchronizedAllow(viewer,PETSC_TRUE);
1839: PetscViewerASCIISynchronizedPrintf(viewer,"[%d]Vector Block stash size %D block size %D\n",rank,s->n,s->bs);
1840: for (i=0; i<s->n; i++) {
1841: PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Element %D ",rank,s->idx[i]);
1842: for (j=0; j<s->bs; j++) {
1843: val = s->array[i*s->bs+j];
1844: #if defined(PETSC_USE_COMPLEX)
1845: PetscViewerASCIISynchronizedPrintf(viewer,"(%18.16e %18.16e) ",PetscRealPart(val),PetscImaginaryPart(val));
1846: #else
1847: PetscViewerASCIISynchronizedPrintf(viewer,"%18.16e ",val);
1848: #endif
1849: }
1850: PetscViewerASCIISynchronizedPrintf(viewer,"\n");
1851: }
1852: PetscViewerFlush(viewer);
1854: s = &v->stash;
1856: /* print basic stash */
1857: PetscViewerASCIISynchronizedPrintf(viewer,"[%d]Vector stash size %D\n",rank,s->n);
1858: for (i=0; i<s->n; i++) {
1859: val = s->array[i];
1860: #if defined(PETSC_USE_COMPLEX)
1861: PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Element %D (%18.16e %18.16e) ",rank,s->idx[i],PetscRealPart(val),PetscImaginaryPart(val));
1862: #else
1863: PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Element %D %18.16e\n",rank,s->idx[i],val);
1864: #endif
1865: }
1866: PetscViewerFlush(viewer);
1867: PetscViewerASCIISynchronizedAllow(viewer,PETSC_FALSE);
1869: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
1870: return(0);
1871: }
1875: PetscErrorCode PetscOptionsVec(const char key[],const char text[],const char man[],Vec v,PetscBool *set)
1876: {
1877: PetscInt i,N,rstart,rend;
1879: PetscScalar *xx;
1880: PetscReal *xreal;
1881: PetscBool iset;
1884: VecGetOwnershipRange(v,&rstart,&rend);
1885: VecGetSize(v,&N);
1886: PetscMalloc(N*sizeof(PetscReal),&xreal);
1887: PetscMemzero(xreal,N*sizeof(PetscReal));
1888: PetscOptionsRealArray(key,text,man,xreal,&N,&iset);
1889: if (iset) {
1890: VecGetArray(v,&xx);
1891: for (i=rstart; i<rend; i++) xx[i-rstart] = xreal[i];
1892: VecRestoreArray(v,&xx);
1893: }
1894: PetscFree(xreal);
1895: if (set) *set = iset;
1896: return(0);
1897: }
1901: /*@
1902: VecGetLayout - get PetscLayout describing vector layout
1904: Not Collective
1906: Input Arguments:
1907: . x - the vector
1909: Output Arguments:
1910: . map - the layout
1912: Level: developer
1914: .seealso: VecGetSizes(), VecGetOwnershipRange(), VecGetOwnershipRanges()
1915: @*/
1916: PetscErrorCode VecGetLayout(Vec x,PetscLayout *map)
1917: {
1921: *map = x->map;
1922: return(0);
1923: }
1927: /*@
1928: VecSetLayout - set PetscLayout describing vector layout
1930: Not Collective
1932: Input Arguments:
1933: + x - the vector
1934: - map - the layout
1936: Notes:
1937: It is normally only valid to replace the layout with a layout known to be equivalent.
1939: Level: developer
1941: .seealso: VecGetLayout(), VecGetSizes(), VecGetOwnershipRange(), VecGetOwnershipRanges()
1942: @*/
1943: PetscErrorCode VecSetLayout(Vec x,PetscLayout map)
1944: {
1949: PetscLayoutReference(map,&x->map);
1950: return(0);
1951: }