Actual source code: vector.c
petsc-3.5.4 2015-05-23
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;
17: PetscLogEvent VEC_ViennaCLCopyFromGPU, VEC_ViennaCLCopyToGPU;
19: extern PetscErrorCode VecStashGetInfo_Private(VecStash*,PetscInt*,PetscInt*);
22: /*@
23: VecStashGetInfo - Gets how many values are currently in the vector stash, i.e. need
24: to be communicated to other processors during the VecAssemblyBegin/End() process
26: Not collective
28: Input Parameter:
29: . vec - the vector
31: Output Parameters:
32: + nstash - the size of the stash
33: . reallocs - the number of additional mallocs incurred.
34: . bnstash - the size of the block stash
35: - breallocs - the number of additional mallocs incurred.in the block stash
37: Level: advanced
39: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), Vec, VecStashSetInitialSize(), VecStashView()
41: @*/
42: PetscErrorCode VecStashGetInfo(Vec vec,PetscInt *nstash,PetscInt *reallocs,PetscInt *bnstash,PetscInt *breallocs)
43: {
47: VecStashGetInfo_Private(&vec->stash,nstash,reallocs);
48: VecStashGetInfo_Private(&vec->bstash,bnstash,breallocs);
49: return(0);
50: }
54: /*@
55: VecSetLocalToGlobalMapping - Sets a local numbering to global numbering used
56: by the routine VecSetValuesLocal() to allow users to insert vector entries
57: using a local (per-processor) numbering.
59: Logically Collective on Vec
61: Input Parameters:
62: + x - vector
63: - mapping - mapping created with ISLocalToGlobalMappingCreate() or ISLocalToGlobalMappingCreateIS()
65: Notes:
66: All vectors obtained with VecDuplicate() from this vector inherit the same mapping.
68: Level: intermediate
70: Concepts: vector^setting values with local numbering
72: seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetValuesLocal(),
73: VecSetLocalToGlobalMapping(), VecSetValuesBlockedLocal()
74: @*/
75: PetscErrorCode VecSetLocalToGlobalMapping(Vec x,ISLocalToGlobalMapping mapping)
76: {
83: if (x->ops->setlocaltoglobalmapping) {
84: (*x->ops->setlocaltoglobalmapping)(x,mapping);
85: } else {
86: PetscLayoutSetISLocalToGlobalMapping(x->map,mapping);
87: }
88: return(0);
89: }
93: /*@
94: VecGetLocalToGlobalMapping - Gets the local-to-global numbering set by VecSetLocalToGlobalMapping()
96: Not Collective
98: Input Parameter:
99: . X - the vector
101: Output Parameter:
102: . mapping - the mapping
104: Level: advanced
106: Concepts: vectors^local to global mapping
107: Concepts: local to global mapping^for vectors
109: .seealso: VecSetValuesLocal()
110: @*/
111: PetscErrorCode VecGetLocalToGlobalMapping(Vec X,ISLocalToGlobalMapping *mapping)
112: {
117: *mapping = X->map->mapping;
118: return(0);
119: }
123: /*@
124: VecAssemblyBegin - Begins assembling the vector. This routine should
125: be called after completing all calls to VecSetValues().
127: Collective on Vec
129: Input Parameter:
130: . vec - the vector
132: Level: beginner
134: Concepts: assembly^vectors
136: .seealso: VecAssemblyEnd(), VecSetValues()
137: @*/
138: PetscErrorCode VecAssemblyBegin(Vec vec)
139: {
145: VecStashViewFromOptions(vec,NULL,"-vec_view_stash");
146: PetscLogEventBegin(VEC_AssemblyBegin,vec,0,0,0);
147: if (vec->ops->assemblybegin) {
148: (*vec->ops->assemblybegin)(vec);
149: }
150: PetscLogEventEnd(VEC_AssemblyBegin,vec,0,0,0);
151: PetscObjectStateIncrease((PetscObject)vec);
152: return(0);
153: }
157: /*@
158: VecAssemblyEnd - Completes assembling the vector. This routine should
159: be called after VecAssemblyBegin().
161: Collective on Vec
163: Input Parameter:
164: . vec - the vector
166: Options Database Keys:
167: + -vec_view - Prints vector in ASCII format
168: . -vec_view ::ascii_matlab - Prints vector in ASCII MATLAB format to stdout
169: . -vec_view matlab:filename - Prints vector in MATLAB format to matlaboutput.mat
170: . -vec_view draw - Activates vector viewing using drawing tools
171: . -display <name> - Sets display name (default is host)
172: . -draw_pause <sec> - Sets number of seconds to pause after display
173: - -vec_view socket - Activates vector viewing using a socket
175: Level: beginner
177: .seealso: VecAssemblyBegin(), VecSetValues()
178: @*/
179: PetscErrorCode VecAssemblyEnd(Vec vec)
180: {
185: PetscLogEventBegin(VEC_AssemblyEnd,vec,0,0,0);
187: if (vec->ops->assemblyend) {
188: (*vec->ops->assemblyend)(vec);
189: }
190: PetscLogEventEnd(VEC_AssemblyEnd,vec,0,0,0);
191: VecViewFromOptions(vec,NULL,"-vec_view");
192: return(0);
193: }
197: /*@
198: VecPointwiseMax - Computes the componentwise maximum w_i = max(x_i, y_i).
200: Logically Collective on Vec
202: Input Parameters:
203: . x, y - the vectors
205: Output Parameter:
206: . w - the result
208: Level: advanced
210: Notes: any subset of the x, y, and w may be the same vector.
211: For complex numbers compares only the real part
213: Concepts: vector^pointwise multiply
215: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
216: @*/
217: PetscErrorCode VecPointwiseMax(Vec w,Vec x,Vec y)
218: {
230: 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");
231: 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");
233: (*w->ops->pointwisemax)(w,x,y);
234: PetscObjectStateIncrease((PetscObject)w);
235: return(0);
236: }
241: /*@
242: VecPointwiseMin - Computes the componentwise minimum w_i = min(x_i, y_i).
244: Logically Collective on Vec
246: Input Parameters:
247: . x, y - the vectors
249: Output Parameter:
250: . w - the result
252: Level: advanced
254: Notes: any subset of the x, y, and w may be the same vector.
255: For complex numbers compares only the real part
257: Concepts: vector^pointwise multiply
259: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
260: @*/
261: PetscErrorCode VecPointwiseMin(Vec w,Vec x,Vec y)
262: {
274: 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");
275: 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");
277: (*w->ops->pointwisemin)(w,x,y);
278: PetscObjectStateIncrease((PetscObject)w);
279: return(0);
280: }
284: /*@
285: VecPointwiseMaxAbs - Computes the componentwise maximum of the absolute values w_i = max(abs(x_i), abs(y_i)).
287: Logically Collective on Vec
289: Input Parameters:
290: . x, y - the vectors
292: Output Parameter:
293: . w - the result
295: Level: advanced
297: Notes: any subset of the x, y, and w may be the same vector.
299: Concepts: vector^pointwise multiply
301: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMin(), VecPointwiseMax(), VecMaxPointwiseDivide()
302: @*/
303: PetscErrorCode VecPointwiseMaxAbs(Vec w,Vec x,Vec y)
304: {
316: 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");
317: 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");
319: (*w->ops->pointwisemaxabs)(w,x,y);
320: PetscObjectStateIncrease((PetscObject)w);
321: return(0);
322: }
326: /*@
327: VecPointwiseDivide - Computes the componentwise division w = x/y.
329: Logically Collective on Vec
331: Input Parameters:
332: . x, y - the vectors
334: Output Parameter:
335: . w - the result
337: Level: advanced
339: Notes: any subset of the x, y, and w may be the same vector.
341: Concepts: vector^pointwise divide
343: .seealso: VecPointwiseMult(), VecPointwiseMax(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
344: @*/
345: PetscErrorCode VecPointwiseDivide(Vec w,Vec x,Vec y)
346: {
358: 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");
359: 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");
361: (*w->ops->pointwisedivide)(w,x,y);
362: PetscObjectStateIncrease((PetscObject)w);
363: return(0);
364: }
369: /*@
370: VecDuplicate - Creates a new vector of the same type as an existing vector.
372: Collective on Vec
374: Input Parameters:
375: . v - a vector to mimic
377: Output Parameter:
378: . newv - location to put new vector
380: Notes:
381: VecDuplicate() DOES NOT COPY the vector entries, but rather allocates storage
382: for the new vector. Use VecCopy() to copy a vector.
384: Use VecDestroy() to free the space. Use VecDuplicateVecs() to get several
385: vectors.
387: Level: beginner
389: .seealso: VecDestroy(), VecDuplicateVecs(), VecCreate(), VecCopy()
390: @*/
391: PetscErrorCode VecDuplicate(Vec v,Vec *newv)
392: {
399: (*v->ops->duplicate)(v,newv);
400: PetscObjectStateIncrease((PetscObject)*newv);
401: return(0);
402: }
406: /*@
407: VecDestroy - Destroys a vector.
409: Collective on Vec
411: Input Parameters:
412: . v - the vector
414: Level: beginner
416: .seealso: VecDuplicate(), VecDestroyVecs()
417: @*/
418: PetscErrorCode VecDestroy(Vec *v)
419: {
423: if (!*v) return(0);
425: if (--((PetscObject)(*v))->refct > 0) {*v = 0; return(0);}
427: PetscObjectSAWsViewOff((PetscObject)*v);
428: /* destroy the internal part */
429: if ((*v)->ops->destroy) {
430: (*(*v)->ops->destroy)(*v);
431: }
432: /* destroy the external/common part */
433: PetscLayoutDestroy(&(*v)->map);
434: PetscHeaderDestroy(v);
435: return(0);
436: }
440: /*@C
441: VecDuplicateVecs - Creates several vectors of the same type as an existing vector.
443: Collective on Vec
445: Input Parameters:
446: + m - the number of vectors to obtain
447: - v - a vector to mimic
449: Output Parameter:
450: . V - location to put pointer to array of vectors
452: Notes:
453: Use VecDestroyVecs() to free the space. Use VecDuplicate() to form a single
454: vector.
456: Fortran Note:
457: The Fortran interface is slightly different from that given below, it
458: requires one to pass in V a Vec (integer) array of size at least m.
459: See the Fortran chapter of the users manual and petsc/src/vec/vec/examples for details.
461: Level: intermediate
463: .seealso: VecDestroyVecs(), VecDuplicate(), VecCreate(), VecDuplicateVecsF90()
464: @*/
465: PetscErrorCode VecDuplicateVecs(Vec v,PetscInt m,Vec *V[])
466: {
473: (*v->ops->duplicatevecs)(v, m,V);
474: return(0);
475: }
479: /*@C
480: VecDestroyVecs - Frees a block of vectors obtained with VecDuplicateVecs().
482: Collective on Vec
484: Input Parameters:
485: + vv - pointer to pointer to array of vector pointers
486: - m - the number of vectors previously obtained
488: Fortran Note:
489: The Fortran interface is slightly different from that given below.
490: See the Fortran chapter of the users manual
492: Level: intermediate
494: .seealso: VecDuplicateVecs(), VecDestroyVecsf90()
495: @*/
496: PetscErrorCode VecDestroyVecs(PetscInt m,Vec *vv[])
497: {
502: if (!*vv) return(0);
505: if (m < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Trying to destroy negative number of vectors %D",m);
506: (*(**vv)->ops->destroyvecs)(m,*vv);
507: *vv = 0;
508: return(0);
509: }
513: /*@C
514: VecView - Views a vector object.
516: Collective on Vec
518: Input Parameters:
519: + vec - the vector
520: - viewer - an optional visualization context
522: Notes:
523: The available visualization contexts include
524: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
525: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
526: output where only the first processor opens
527: the file. All other processors send their
528: data to the first processor to print.
530: You can change the format the vector is printed using the
531: option PetscViewerSetFormat().
533: The user can open alternative visualization contexts with
534: + PetscViewerASCIIOpen() - Outputs vector to a specified file
535: . PetscViewerBinaryOpen() - Outputs vector in binary to a
536: specified file; corresponding input uses VecLoad()
537: . PetscViewerDrawOpen() - Outputs vector to an X window display
538: - PetscViewerSocketOpen() - Outputs vector to Socket viewer
540: The user can call PetscViewerSetFormat() to specify the output
541: format of ASCII printed objects (when using PETSC_VIEWER_STDOUT_SELF,
542: PETSC_VIEWER_STDOUT_WORLD and PetscViewerASCIIOpen). Available formats include
543: + PETSC_VIEWER_DEFAULT - default, prints vector contents
544: . PETSC_VIEWER_ASCII_MATLAB - prints vector contents in MATLAB format
545: . PETSC_VIEWER_ASCII_INDEX - prints vector contents, including indices of vector elements
546: - PETSC_VIEWER_ASCII_COMMON - prints vector contents, using a
547: format common among all vector types
549: Notes for HDF5 Viewer: the name of the Vec (given with PetscObjectSetName() is the name that is used
550: for the object in the HDF5 file. If you wish to store the same vector to the HDF5 viewer (with different values,
551: obviously) several times, you must change its name each time before calling the VecView(). The name you use
552: here should equal the name that you use in the Vec object that you use with VecLoad().
554: See the manual page for VecLoad() on the exact format the binary viewer stores
555: the values in the file.
557: Level: beginner
559: Concepts: vector^printing
560: Concepts: vector^saving to disk
562: .seealso: PetscViewerASCIIOpen(), PetscViewerDrawOpen(), PetscDrawLGCreate(),
563: PetscViewerSocketOpen(), PetscViewerBinaryOpen(), VecLoad(), PetscViewerCreate(),
564: PetscRealView(), PetscScalarView(), PetscIntView()
565: @*/
566: PetscErrorCode VecView(Vec vec,PetscViewer viewer)
567: {
568: PetscErrorCode ierr;
569: PetscBool iascii;
570: PetscViewerFormat format;
575: if (!viewer) {
576: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)vec),&viewer);
577: }
580: if (vec->stash.n || vec->bstash.n) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call VecAssemblyBegin/End() before viewing this vector");
582: PetscLogEventBegin(VEC_View,vec,viewer,0,0);
583: PetscViewerGetFormat(viewer,&format);
584: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
585: if (iascii) {
586: PetscInt rows,bs;
588: PetscObjectPrintClassNamePrefixType((PetscObject)vec,viewer);
589: if (format == PETSC_VIEWER_ASCII_INFO || format == PETSC_VIEWER_ASCII_INFO_DETAIL) {
590: PetscViewerASCIIPushTab(viewer);
591: VecGetSize(vec,&rows);
592: VecGetBlockSize(vec,&bs);
593: if (bs != 1) {
594: PetscViewerASCIIPrintf(viewer,"length=%D, bs=%D\n",rows,bs);
595: } else {
596: PetscViewerASCIIPrintf(viewer,"length=%D\n",rows);
597: }
598: PetscViewerASCIIPopTab(viewer);
599: }
600: }
601: if (format == PETSC_VIEWER_NATIVE && vec->ops->viewnative) {
602: (*vec->ops->viewnative)(vec,viewer);
603: } else {
604: (*vec->ops->view)(vec,viewer);
605: }
606: PetscLogEventEnd(VEC_View,vec,viewer,0,0);
607: return(0);
608: }
610: #if defined(PETSC_USE_DEBUG)
611: #include <../src/sys/totalview/tv_data_display.h>
612: PETSC_UNUSED static int TV_display_type(const struct _p_Vec *v)
613: {
614: const PetscScalar *values;
615: char type[32];
616: PetscErrorCode ierr;
619: TV_add_row("Local rows", "int", &v->map->n);
620: TV_add_row("Global rows", "int", &v->map->N);
621: TV_add_row("Typename", TV_ascii_string_type, ((PetscObject)v)->type_name);
622: VecGetArrayRead((Vec)v,&values);
623: PetscSNPrintf(type,32,"double[%d]",v->map->n);
624: TV_add_row("values",type, values);
625: VecRestoreArrayRead((Vec)v,&values);
626: return TV_format_OK;
627: }
628: #endif
632: /*@
633: VecGetSize - Returns the global number of elements of the vector.
635: Not Collective
637: Input Parameter:
638: . x - the vector
640: Output Parameters:
641: . size - the global length of the vector
643: Level: beginner
645: Concepts: vector^local size
647: .seealso: VecGetLocalSize()
648: @*/
649: PetscErrorCode VecGetSize(Vec x,PetscInt *size)
650: {
657: (*x->ops->getsize)(x,size);
658: return(0);
659: }
663: /*@
664: VecGetLocalSize - Returns the number of elements of the vector stored
665: in local memory. This routine may be implementation dependent, so use
666: with care.
668: Not Collective
670: Input Parameter:
671: . x - the vector
673: Output Parameter:
674: . size - the length of the local piece of the vector
676: Level: beginner
678: Concepts: vector^size
680: .seealso: VecGetSize()
681: @*/
682: PetscErrorCode VecGetLocalSize(Vec x,PetscInt *size)
683: {
690: (*x->ops->getlocalsize)(x,size);
691: return(0);
692: }
696: /*@C
697: VecGetOwnershipRange - Returns the range of indices owned by
698: this processor, assuming that the vectors are laid out with the
699: first n1 elements on the first processor, next n2 elements on the
700: second, etc. For certain parallel layouts this range may not be
701: well defined.
703: Not Collective
705: Input Parameter:
706: . x - the vector
708: Output Parameters:
709: + low - the first local element, pass in NULL if not interested
710: - high - one more than the last local element, pass in NULL if not interested
712: Note:
713: The high argument is one more than the last element stored locally.
715: Fortran: NULL_INTEGER should be used instead of NULL
717: Level: beginner
719: Concepts: ownership^of vectors
720: Concepts: vector^ownership of elements
722: .seealso: MatGetOwnershipRange(), MatGetOwnershipRanges(), VecGetOwnershipRanges()
723: @*/
724: PetscErrorCode VecGetOwnershipRange(Vec x,PetscInt *low,PetscInt *high)
725: {
731: if (low) *low = x->map->rstart;
732: if (high) *high = x->map->rend;
733: return(0);
734: }
738: /*@C
739: VecGetOwnershipRanges - Returns the range of indices owned by EACH processor,
740: assuming that the vectors are laid out with the
741: first n1 elements on the first processor, next n2 elements on the
742: second, etc. For certain parallel layouts this range may not be
743: well defined.
745: Not Collective
747: Input Parameter:
748: . x - the vector
750: Output Parameters:
751: . range - array of length size+1 with the start and end+1 for each process
753: Note:
754: The high argument is one more than the last element stored locally.
756: Fortran: You must PASS in an array of length size+1
758: Level: beginner
760: Concepts: ownership^of vectors
761: Concepts: vector^ownership of elements
763: .seealso: MatGetOwnershipRange(), MatGetOwnershipRanges(), VecGetOwnershipRange()
764: @*/
765: PetscErrorCode VecGetOwnershipRanges(Vec x,const PetscInt *ranges[])
766: {
772: PetscLayoutGetRanges(x->map,ranges);
773: return(0);
774: }
778: /*@
779: VecSetOption - Sets an option for controling a vector's behavior.
781: Collective on Vec
783: Input Parameter:
784: + x - the vector
785: . op - the option
786: - flag - turn the option on or off
788: Supported Options:
789: + VEC_IGNORE_OFF_PROC_ENTRIES, which causes VecSetValues() to ignore
790: entries destined to be stored on a separate processor. This can be used
791: to eliminate the global reduction in the VecAssemblyXXXX() if you know
792: that you have only used VecSetValues() to set local elements
793: . VEC_IGNORE_NEGATIVE_INDICES, which means you can pass negative indices
794: in ix in calls to VecSetValues() or VecGetValues(). These rows are simply
795: ignored.
797: Level: intermediate
799: @*/
800: PetscErrorCode VecSetOption(Vec x,VecOption op,PetscBool flag)
801: {
807: if (x->ops->setoption) {
808: (*x->ops->setoption)(x,op,flag);
809: }
810: return(0);
811: }
815: /* Default routines for obtaining and releasing; */
816: /* may be used by any implementation */
817: PetscErrorCode VecDuplicateVecs_Default(Vec w,PetscInt m,Vec *V[])
818: {
820: PetscInt i;
825: if (m <= 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"m must be > 0: m = %D",m);
826: PetscMalloc1(m,V);
827: for (i=0; i<m; i++) {VecDuplicate(w,*V+i);}
828: return(0);
829: }
833: PetscErrorCode VecDestroyVecs_Default(PetscInt m,Vec v[])
834: {
836: PetscInt i;
840: if (m <= 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"m must be > 0: m = %D",m);
841: for (i=0; i<m; i++) {VecDestroy(&v[i]);}
842: PetscFree(v);
843: return(0);
844: }
848: /*@
849: VecResetArray - Resets a vector to use its default memory. Call this
850: after the use of VecPlaceArray().
852: Not Collective
854: Input Parameters:
855: . vec - the vector
857: Level: developer
859: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecPlaceArray()
861: @*/
862: PetscErrorCode VecResetArray(Vec vec)
863: {
869: if (vec->ops->resetarray) {
870: (*vec->ops->resetarray)(vec);
871: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Cannot reset array in this type of vector");
872: PetscObjectStateIncrease((PetscObject)vec);
873: return(0);
874: }
878: /*@C
879: VecLoad - Loads a vector that has been stored in binary or HDF5 format
880: with VecView().
882: Collective on PetscViewer
884: Input Parameters:
885: + newvec - the newly loaded vector, this needs to have been created with VecCreate() or
886: some related function before a call to VecLoad().
887: - viewer - binary file viewer, obtained from PetscViewerBinaryOpen() or
888: HDF5 file viewer, obtained from PetscViewerHDF5Open()
890: Level: intermediate
892: Notes:
893: Defaults to the standard Seq or MPI Vec, if you want some other type of Vec call VecSetFromOptions()
894: before calling this.
896: The input file must contain the full global vector, as
897: written by the routine VecView().
899: If the type or size of newvec is not set before a call to VecLoad, PETSc
900: sets the type and the local and global sizes.If type and/or
901: sizes are already set, then the same are used.
903: IF using HDF5, you must assign the Vec the same name as was used in the Vec
904: that was stored in the file using PetscObjectSetName(). Otherwise you will
905: get the error message: "Cannot H5DOpen2() with Vec name NAMEOFOBJECT"
907: Notes for advanced users:
908: Most users should not need to know the details of the binary storage
909: format, since VecLoad() and VecView() completely hide these details.
910: But for anyone who's interested, the standard binary matrix storage
911: format is
912: .vb
913: int VEC_FILE_CLASSID
914: int number of rows
915: PetscScalar *values of all entries
916: .ve
918: In addition, PETSc automatically does the byte swapping for
919: machines that store the bytes reversed, e.g. DEC alpha, freebsd,
920: linux, Windows and the paragon; thus if you write your own binary
921: read/write routines you have to swap the bytes; see PetscBinaryRead()
922: and PetscBinaryWrite() to see how this may be done.
924: Concepts: vector^loading from file
926: .seealso: PetscViewerBinaryOpen(), VecView(), MatLoad(), VecLoad()
927: @*/
928: PetscErrorCode VecLoad(Vec newvec, PetscViewer viewer)
929: {
931: PetscBool isbinary,ishdf5;
932: PetscViewerFormat format;
937: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
938: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5);
939: if (!isbinary && !ishdf5) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Invalid viewer; open viewer with PetscViewerBinaryOpen()");
941: PetscLogEventBegin(VEC_Load,viewer,0,0,0);
942: if (!((PetscObject)newvec)->type_name && !newvec->ops->create) {
943: VecSetType(newvec, VECSTANDARD);
944: }
945: PetscViewerGetFormat(viewer,&format);
946: if (format == PETSC_VIEWER_NATIVE && newvec->ops->loadnative) {
947: (*newvec->ops->loadnative)(newvec,viewer);
948: } else {
949: (*newvec->ops->load)(newvec,viewer);
950: }
951: PetscLogEventEnd(VEC_Load,viewer,0,0,0);
952: return(0);
953: }
958: /*@
959: VecReciprocal - Replaces each component of a vector by its reciprocal.
961: Logically Collective on Vec
963: Input Parameter:
964: . vec - the vector
966: Output Parameter:
967: . vec - the vector reciprocal
969: Level: intermediate
971: Concepts: vector^reciprocal
973: .seealso: VecLog(), VecExp(), VecSqrtAbs()
975: @*/
976: PetscErrorCode VecReciprocal(Vec vec)
977: {
983: if (vec->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
984: if (!vec->ops->reciprocal) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Vector does not support reciprocal operation");
985: (*vec->ops->reciprocal)(vec);
986: PetscObjectStateIncrease((PetscObject)vec);
987: return(0);
988: }
992: /*@C
993: VecSetOperation - Allows user to set a vector operation.
995: Logically Collective on Vec
997: Input Parameters:
998: + vec - the vector
999: . op - the name of the operation
1000: - f - the function that provides the operation.
1002: Level: advanced
1004: Usage:
1005: $ PetscErrorCode userview(Vec,PetscViewer);
1006: $ VecCreateMPI(comm,m,M,&x);
1007: $ VecSetOperation(x,VECOP_VIEW,(void(*)(void))userview);
1009: Notes:
1010: See the file include/petscvec.h for a complete list of matrix
1011: operations, which all have the form VECOP_<OPERATION>, where
1012: <OPERATION> is the name (in all capital letters) of the
1013: user interface routine (e.g., VecView() -> VECOP_VIEW).
1015: This function is not currently available from Fortran.
1017: .keywords: vector, set, operation
1019: .seealso: VecCreate(), MatShellSetOperation()
1020: @*/
1021: PetscErrorCode VecSetOperation(Vec vec,VecOperation op, void (*f)(void))
1022: {
1025: if (op == VECOP_VIEW && !vec->ops->viewnative) {
1026: vec->ops->viewnative = vec->ops->view;
1027: } else if (op == VECOP_LOAD && !vec->ops->loadnative) {
1028: vec->ops->loadnative = vec->ops->load;
1029: }
1030: (((void(**)(void))vec->ops)[(int)op]) = f;
1031: return(0);
1032: }
1037: /*@
1038: VecStashSetInitialSize - sets the sizes of the vec-stash, that is
1039: used during the assembly process to store values that belong to
1040: other processors.
1042: Not Collective, different processes can have different size stashes
1044: Input Parameters:
1045: + vec - the vector
1046: . size - the initial size of the stash.
1047: - bsize - the initial size of the block-stash(if used).
1049: Options Database Keys:
1050: + -vecstash_initial_size <size> or <size0,size1,...sizep-1>
1051: - -vecstash_block_initial_size <bsize> or <bsize0,bsize1,...bsizep-1>
1053: Level: intermediate
1055: Notes:
1056: The block-stash is used for values set with VecSetValuesBlocked() while
1057: the stash is used for values set with VecSetValues()
1059: Run with the option -info and look for output of the form
1060: VecAssemblyBegin_MPIXXX:Stash has MM entries, uses nn mallocs.
1061: to determine the appropriate value, MM, to use for size and
1062: VecAssemblyBegin_MPIXXX:Block-Stash has BMM entries, uses nn mallocs.
1063: to determine the value, BMM to use for bsize
1065: Concepts: vector^stash
1066: Concepts: stash^vector
1068: .seealso: VecSetBlockSize(), VecSetValues(), VecSetValuesBlocked(), VecStashView()
1070: @*/
1071: PetscErrorCode VecStashSetInitialSize(Vec vec,PetscInt size,PetscInt bsize)
1072: {
1077: VecStashSetInitialSize_Private(&vec->stash,size);
1078: VecStashSetInitialSize_Private(&vec->bstash,bsize);
1079: return(0);
1080: }
1084: /*@
1085: VecConjugate - Conjugates a vector.
1087: Logically Collective on Vec
1089: Input Parameters:
1090: . x - the vector
1092: Level: intermediate
1094: Concepts: vector^conjugate
1096: @*/
1097: PetscErrorCode VecConjugate(Vec x)
1098: {
1099: #if defined(PETSC_USE_COMPLEX)
1105: if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1106: (*x->ops->conjugate)(x);
1107: /* we need to copy norms here */
1108: PetscObjectStateIncrease((PetscObject)x);
1109: return(0);
1110: #else
1111: return(0);
1112: #endif
1113: }
1117: /*@
1118: VecPointwiseMult - Computes the componentwise multiplication w = x*y.
1120: Logically Collective on Vec
1122: Input Parameters:
1123: . x, y - the vectors
1125: Output Parameter:
1126: . w - the result
1128: Level: advanced
1130: Notes: any subset of the x, y, and w may be the same vector.
1132: Concepts: vector^pointwise multiply
1134: .seealso: VecPointwiseDivide(), VecPointwiseMax(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
1135: @*/
1136: PetscErrorCode VecPointwiseMult(Vec w, Vec x,Vec y)
1137: {
1149: 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");
1151: PetscLogEventBegin(VEC_PointwiseMult,x,y,w,0);
1152: (*w->ops->pointwisemult)(w,x,y);
1153: PetscLogEventEnd(VEC_PointwiseMult,x,y,w,0);
1154: PetscObjectStateIncrease((PetscObject)w);
1155: return(0);
1156: }
1160: /*@
1161: VecSetRandom - Sets all components of a vector to random numbers.
1163: Logically Collective on Vec
1165: Input Parameters:
1166: + x - the vector
1167: - rctx - the random number context, formed by PetscRandomCreate(), or NULL and
1168: it will create one internally.
1170: Output Parameter:
1171: . x - the vector
1173: Example of Usage:
1174: .vb
1175: PetscRandomCreate(PETSC_COMM_WORLD,&rctx);
1176: VecSetRandom(x,rctx);
1177: PetscRandomDestroy(rctx);
1178: .ve
1180: Level: intermediate
1182: Concepts: vector^setting to random
1183: Concepts: random^vector
1185: .seealso: VecSet(), VecSetValues(), PetscRandomCreate(), PetscRandomDestroy()
1186: @*/
1187: PetscErrorCode VecSetRandom(Vec x,PetscRandom rctx)
1188: {
1190: PetscRandom randObj = NULL;
1196: if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1198: if (!rctx) {
1199: MPI_Comm comm;
1200: PetscObjectGetComm((PetscObject)x,&comm);
1201: PetscRandomCreate(comm,&randObj);
1202: PetscRandomSetFromOptions(randObj);
1203: rctx = randObj;
1204: }
1206: PetscLogEventBegin(VEC_SetRandom,x,rctx,0,0);
1207: (*x->ops->setrandom)(x,rctx);
1208: PetscLogEventEnd(VEC_SetRandom,x,rctx,0,0);
1210: PetscRandomDestroy(&randObj);
1211: PetscObjectStateIncrease((PetscObject)x);
1212: return(0);
1213: }
1217: /*@
1218: VecZeroEntries - puts a 0.0 in each element of a vector
1220: Logically Collective on Vec
1222: Input Parameter:
1223: . vec - The vector
1225: Level: beginner
1227: Developer Note: This routine does not need to exist since the exact functionality is obtained with
1228: VecSet(vec,0); I guess someone added it to mirror the functionality of MatZeroEntries() but Mat is nothing
1229: like a Vec (one is an operator and one is an element of a vector space, yeah yeah dual blah blah blah) so
1230: this routine should not exist.
1232: .keywords: Vec, set, options, database
1233: .seealso: VecCreate(), VecSetOptionsPrefix(), VecSet(), VecSetValues()
1234: @*/
1235: PetscErrorCode VecZeroEntries(Vec vec)
1236: {
1240: VecSet(vec,0);
1241: return(0);
1242: }
1246: /*
1247: VecSetTypeFromOptions_Private - Sets the type of vector from user options. Defaults to a PETSc sequential vector on one
1248: processor and a PETSc MPI vector on more than one processor.
1250: Collective on Vec
1252: Input Parameter:
1253: . vec - The vector
1255: Level: intermediate
1257: .keywords: Vec, set, options, database, type
1258: .seealso: VecSetFromOptions(), VecSetType()
1259: */
1260: static PetscErrorCode VecSetTypeFromOptions_Private(Vec vec)
1261: {
1262: PetscBool opt;
1263: VecType defaultType;
1264: char typeName[256];
1265: PetscMPIInt size;
1269: if (((PetscObject)vec)->type_name) defaultType = ((PetscObject)vec)->type_name;
1270: else {
1271: MPI_Comm_size(PetscObjectComm((PetscObject)vec), &size);
1272: if (size > 1) defaultType = VECMPI;
1273: else defaultType = VECSEQ;
1274: }
1276: if (!VecRegisterAllCalled) {VecRegisterAll();}
1277: PetscOptionsFList("-vec_type","Vector type","VecSetType",VecList,defaultType,typeName,256,&opt);
1278: if (opt) {
1279: VecSetType(vec, typeName);
1280: } else {
1281: VecSetType(vec, defaultType);
1282: }
1283: return(0);
1284: }
1288: /*@
1289: VecSetFromOptions - Configures the vector from the options database.
1291: Collective on Vec
1293: Input Parameter:
1294: . vec - The vector
1296: Notes: To see all options, run your program with the -help option, or consult the users manual.
1297: Must be called after VecCreate() but before the vector is used.
1299: Level: beginner
1301: Concepts: vectors^setting options
1302: Concepts: vectors^setting type
1304: .keywords: Vec, set, options, database
1305: .seealso: VecCreate(), VecSetOptionsPrefix()
1306: @*/
1307: PetscErrorCode VecSetFromOptions(Vec vec)
1308: {
1314: PetscObjectOptionsBegin((PetscObject)vec);
1315: /* Handle vector type options */
1316: VecSetTypeFromOptions_Private(vec);
1318: /* Handle specific vector options */
1319: if (vec->ops->setfromoptions) {
1320: (*vec->ops->setfromoptions)(vec);
1321: }
1323: /* process any options handlers added with PetscObjectAddOptionsHandler() */
1324: PetscObjectProcessOptionsHandlers((PetscObject)vec);
1325: PetscOptionsEnd();
1326: return(0);
1327: }
1331: /*@
1332: VecSetSizes - Sets the local and global sizes, and checks to determine compatibility
1334: Collective on Vec
1336: Input Parameters:
1337: + v - the vector
1338: . n - the local size (or PETSC_DECIDE to have it set)
1339: - N - the global size (or PETSC_DECIDE)
1341: Notes:
1342: n and N cannot be both PETSC_DECIDE
1343: If one processor calls this with N of PETSC_DECIDE then all processors must, otherwise the program will hang.
1345: Level: intermediate
1347: .seealso: VecGetSize(), PetscSplitOwnership()
1348: @*/
1349: PetscErrorCode VecSetSizes(Vec v, PetscInt n, PetscInt N)
1350: {
1356: 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);
1357: 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);
1358: v->map->n = n;
1359: v->map->N = N;
1360: if (v->ops->create) {
1361: (*v->ops->create)(v);
1362: v->ops->create = 0;
1363: }
1364: return(0);
1365: }
1369: /*@
1370: VecSetBlockSize - Sets the blocksize for future calls to VecSetValuesBlocked()
1371: and VecSetValuesBlockedLocal().
1373: Logically Collective on Vec
1375: Input Parameter:
1376: + v - the vector
1377: - bs - the blocksize
1379: Notes:
1380: All vectors obtained by VecDuplicate() inherit the same blocksize.
1382: Level: advanced
1384: .seealso: VecSetValuesBlocked(), VecSetLocalToGlobalMapping(), VecGetBlockSize()
1386: Concepts: block size^vectors
1387: @*/
1388: PetscErrorCode VecSetBlockSize(Vec v,PetscInt bs)
1389: {
1394: if (bs < 0 || bs == v->map->bs) return(0);
1396: PetscLayoutSetBlockSize(v->map,bs);
1397: v->bstash.bs = bs; /* use the same blocksize for the vec's block-stash */
1398: return(0);
1399: }
1403: /*@
1404: VecGetBlockSize - Gets the blocksize for the vector, i.e. what is used for VecSetValuesBlocked()
1405: and VecSetValuesBlockedLocal().
1407: Not Collective
1409: Input Parameter:
1410: . v - the vector
1412: Output Parameter:
1413: . bs - the blocksize
1415: Notes:
1416: All vectors obtained by VecDuplicate() inherit the same blocksize.
1418: Level: advanced
1420: .seealso: VecSetValuesBlocked(), VecSetLocalToGlobalMapping(), VecSetBlockSize()
1422: Concepts: vector^block size
1423: Concepts: block^vector
1425: @*/
1426: PetscErrorCode VecGetBlockSize(Vec v,PetscInt *bs)
1427: {
1433: PetscLayoutGetBlockSize(v->map,bs);
1434: return(0);
1435: }
1439: /*@C
1440: VecSetOptionsPrefix - Sets the prefix used for searching for all
1441: Vec options in the database.
1443: Logically Collective on Vec
1445: Input Parameter:
1446: + v - the Vec context
1447: - prefix - the prefix to prepend to all option names
1449: Notes:
1450: A hyphen (-) must NOT be given at the beginning of the prefix name.
1451: The first character of all runtime options is AUTOMATICALLY the hyphen.
1453: Level: advanced
1455: .keywords: Vec, set, options, prefix, database
1457: .seealso: VecSetFromOptions()
1458: @*/
1459: PetscErrorCode VecSetOptionsPrefix(Vec v,const char prefix[])
1460: {
1465: PetscObjectSetOptionsPrefix((PetscObject)v,prefix);
1466: return(0);
1467: }
1471: /*@C
1472: VecAppendOptionsPrefix - Appends to the prefix used for searching for all
1473: Vec options in the database.
1475: Logically Collective on Vec
1477: Input Parameters:
1478: + v - the Vec context
1479: - prefix - the prefix to prepend to all option names
1481: Notes:
1482: A hyphen (-) must NOT be given at the beginning of the prefix name.
1483: The first character of all runtime options is AUTOMATICALLY the hyphen.
1485: Level: advanced
1487: .keywords: Vec, append, options, prefix, database
1489: .seealso: VecGetOptionsPrefix()
1490: @*/
1491: PetscErrorCode VecAppendOptionsPrefix(Vec v,const char prefix[])
1492: {
1497: PetscObjectAppendOptionsPrefix((PetscObject)v,prefix);
1498: return(0);
1499: }
1503: /*@C
1504: VecGetOptionsPrefix - Sets the prefix used for searching for all
1505: Vec options in the database.
1507: Not Collective
1509: Input Parameter:
1510: . v - the Vec context
1512: Output Parameter:
1513: . prefix - pointer to the prefix string used
1515: Notes: On the fortran side, the user should pass in a string 'prefix' of
1516: sufficient length to hold the prefix.
1518: Level: advanced
1520: .keywords: Vec, get, options, prefix, database
1522: .seealso: VecAppendOptionsPrefix()
1523: @*/
1524: PetscErrorCode VecGetOptionsPrefix(Vec v,const char *prefix[])
1525: {
1530: PetscObjectGetOptionsPrefix((PetscObject)v,prefix);
1531: return(0);
1532: }
1536: /*@
1537: VecSetUp - Sets up the internal vector data structures for the later use.
1539: Collective on Vec
1541: Input Parameters:
1542: . v - the Vec context
1544: Notes:
1545: For basic use of the Vec classes the user need not explicitly call
1546: VecSetUp(), since these actions will happen automatically.
1548: Level: advanced
1550: .keywords: Vec, setup
1552: .seealso: VecCreate(), VecDestroy()
1553: @*/
1554: PetscErrorCode VecSetUp(Vec v)
1555: {
1556: PetscMPIInt size;
1561: if (!((PetscObject)v)->type_name) {
1562: MPI_Comm_size(PetscObjectComm((PetscObject)v), &size);
1563: if (size == 1) {
1564: VecSetType(v, VECSEQ);
1565: } else {
1566: VecSetType(v, VECMPI);
1567: }
1568: }
1569: return(0);
1570: }
1572: /*
1573: These currently expose the PetscScalar/PetscReal in updating the
1574: cached norm. If we push those down into the implementation these
1575: will become independent of PetscScalar/PetscReal
1576: */
1580: /*@
1581: VecCopy - Copies a vector. y <- x
1583: Logically Collective on Vec
1585: Input Parameter:
1586: . x - the vector
1588: Output Parameter:
1589: . y - the copy
1591: Notes:
1592: For default parallel PETSc vectors, both x and y must be distributed in
1593: the same manner; local copies are done.
1595: Level: beginner
1597: .seealso: VecDuplicate()
1598: @*/
1599: PetscErrorCode VecCopy(Vec x,Vec y)
1600: {
1601: PetscBool flgs[4];
1602: PetscReal norms[4] = {0.0,0.0,0.0,0.0};
1604: PetscInt i;
1611: if (x == y) return(0);
1612: if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1613: 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);
1615: #if !defined(PETSC_USE_MIXED_PRECISION)
1616: for (i=0; i<4; i++) {
1617: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[i],norms[i],flgs[i]);
1618: }
1619: #endif
1621: PetscLogEventBegin(VEC_Copy,x,y,0,0);
1622: #if defined(PETSC_USE_MIXED_PRECISION)
1623: extern PetscErrorCode VecGetArray(Vec,double**);
1624: extern PetscErrorCode VecRestoreArray(Vec,double**);
1625: extern PetscErrorCode VecGetArray(Vec,float**);
1626: extern PetscErrorCode VecRestoreArray(Vec,float**);
1627: extern PetscErrorCode VecGetArrayRead(Vec,const double**);
1628: extern PetscErrorCode VecRestoreArrayRead(Vec,const double**);
1629: extern PetscErrorCode VecGetArrayRead(Vec,const float**);
1630: extern PetscErrorCode VecRestoreArrayRead(Vec,const float**);
1631: if ((((PetscObject)x)->precision == PETSC_PRECISION_SINGLE) && (((PetscObject)y)->precision == PETSC_PRECISION_DOUBLE)) {
1632: PetscInt i,n;
1633: const float *xx;
1634: double *yy;
1635: VecGetArrayRead(x,&xx);
1636: VecGetArray(y,&yy);
1637: VecGetLocalSize(x,&n);
1638: for (i=0; i<n; i++) yy[i] = xx[i];
1639: VecRestoreArrayRead(x,&xx);
1640: VecRestoreArray(y,&yy);
1641: } else if ((((PetscObject)x)->precision == PETSC_PRECISION_DOUBLE) && (((PetscObject)y)->precision == PETSC_PRECISION_SINGLE)) {
1642: PetscInt i,n;
1643: float *yy;
1644: const double *xx;
1645: VecGetArrayRead(x,&xx);
1646: VecGetArray(y,&yy);
1647: VecGetLocalSize(x,&n);
1648: for (i=0; i<n; i++) yy[i] = (float) xx[i];
1649: VecRestoreArrayRead(x,&xx);
1650: VecRestoreArray(y,&yy);
1651: } else {
1652: (*x->ops->copy)(x,y);
1653: }
1654: #else
1655: (*x->ops->copy)(x,y);
1656: #endif
1658: PetscObjectStateIncrease((PetscObject)y);
1659: #if !defined(PETSC_USE_MIXED_PRECISION)
1660: for (i=0; i<4; i++) {
1661: if (flgs[i]) {
1662: PetscObjectComposedDataSetReal((PetscObject)y,NormIds[i],norms[i]);
1663: }
1664: }
1665: #endif
1667: PetscLogEventEnd(VEC_Copy,x,y,0,0);
1668: return(0);
1669: }
1673: /*@
1674: VecSwap - Swaps the vectors x and y.
1676: Logically Collective on Vec
1678: Input Parameters:
1679: . x, y - the vectors
1681: Level: advanced
1683: Concepts: vector^swapping values
1685: @*/
1686: PetscErrorCode VecSwap(Vec x,Vec y)
1687: {
1688: PetscReal normxs[4]={0.0,0.0,0.0,0.0},normys[4]={0.0,0.0,0.0,0.0};
1689: PetscBool flgxs[4],flgys[4];
1691: PetscInt i;
1699: if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1700: if (y->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1701: if (x->map->N != y->map->N) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths");
1702: if (x->map->n != y->map->n) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");
1704: PetscLogEventBegin(VEC_Swap,x,y,0,0);
1705: for (i=0; i<4; i++) {
1706: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[i],normxs[i],flgxs[i]);
1707: PetscObjectComposedDataGetReal((PetscObject)y,NormIds[i],normys[i],flgys[i]);
1708: }
1709: (*x->ops->swap)(x,y);
1710: PetscObjectStateIncrease((PetscObject)x);
1711: PetscObjectStateIncrease((PetscObject)y);
1712: for (i=0; i<4; i++) {
1713: if (flgxs[i]) {
1714: PetscObjectComposedDataSetReal((PetscObject)y,NormIds[i],normxs[i]);
1715: }
1716: if (flgys[i]) {
1717: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[i],normys[i]);
1718: }
1719: }
1720: PetscLogEventEnd(VEC_Swap,x,y,0,0);
1721: return(0);
1722: }
1726: /*
1727: VecStashViewFromOptions - Processes command line options to determine if/how an VecStash object is to be viewed.
1729: Collective on VecStash
1731: Input Parameters:
1732: + obj - the VecStash object
1733: . prefix - prefix to use for viewing, or NULL to use prefix of 'mat'
1734: - optionname - option to activate viewing
1736: Level: intermediate
1738: Developer Note: This cannot use PetscObjectViewFromOptions() because it takes a Vec as an argument but does not use VecView
1740: */
1741: PetscErrorCode VecStashViewFromOptions(Vec obj,const char prefix[],const char optionname[])
1742: {
1743: PetscErrorCode ierr;
1744: PetscViewer viewer;
1745: PetscBool flg;
1746: static PetscBool incall = PETSC_FALSE;
1747: PetscViewerFormat format;
1750: if (incall) return(0);
1751: incall = PETSC_TRUE;
1752: PetscOptionsGetViewer(PetscObjectComm((PetscObject)obj),prefix,optionname,&viewer,&format,&flg);
1753: if (flg) {
1754: PetscViewerPushFormat(viewer,format);
1755: VecStashView(obj,viewer);
1756: PetscViewerPopFormat(viewer);
1757: PetscViewerDestroy(&viewer);
1758: }
1759: incall = PETSC_FALSE;
1760: return(0);
1761: }
1765: /*@
1766: VecStashView - Prints the entries in the vector stash and block stash.
1768: Collective on Vec
1770: Input Parameters:
1771: + v - the vector
1772: - viewer - the viewer
1774: Level: advanced
1776: Concepts: vector^stash
1777: Concepts: stash^vector
1779: .seealso: VecSetBlockSize(), VecSetValues(), VecSetValuesBlocked()
1781: @*/
1782: PetscErrorCode VecStashView(Vec v,PetscViewer viewer)
1783: {
1785: PetscMPIInt rank;
1786: PetscInt i,j;
1787: PetscBool match;
1788: VecStash *s;
1789: PetscScalar val;
1796: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&match);
1797: if (!match) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_SUP,"Stash viewer only works with ASCII viewer not %s\n",((PetscObject)v)->type_name);
1798: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
1799: MPI_Comm_rank(PetscObjectComm((PetscObject)v),&rank);
1800: s = &v->bstash;
1802: /* print block stash */
1803: PetscViewerASCIISynchronizedAllow(viewer,PETSC_TRUE);
1804: PetscViewerASCIISynchronizedPrintf(viewer,"[%d]Vector Block stash size %D block size %D\n",rank,s->n,s->bs);
1805: for (i=0; i<s->n; i++) {
1806: PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Element %D ",rank,s->idx[i]);
1807: for (j=0; j<s->bs; j++) {
1808: val = s->array[i*s->bs+j];
1809: #if defined(PETSC_USE_COMPLEX)
1810: PetscViewerASCIISynchronizedPrintf(viewer,"(%18.16e %18.16e) ",PetscRealPart(val),PetscImaginaryPart(val));
1811: #else
1812: PetscViewerASCIISynchronizedPrintf(viewer,"%18.16e ",val);
1813: #endif
1814: }
1815: PetscViewerASCIISynchronizedPrintf(viewer,"\n");
1816: }
1817: PetscViewerFlush(viewer);
1819: s = &v->stash;
1821: /* print basic stash */
1822: PetscViewerASCIISynchronizedPrintf(viewer,"[%d]Vector stash size %D\n",rank,s->n);
1823: for (i=0; i<s->n; i++) {
1824: val = s->array[i];
1825: #if defined(PETSC_USE_COMPLEX)
1826: PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Element %D (%18.16e %18.16e) ",rank,s->idx[i],PetscRealPart(val),PetscImaginaryPart(val));
1827: #else
1828: PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Element %D %18.16e\n",rank,s->idx[i],val);
1829: #endif
1830: }
1831: PetscViewerFlush(viewer);
1832: PetscViewerASCIISynchronizedAllow(viewer,PETSC_FALSE);
1834: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
1835: return(0);
1836: }
1840: PetscErrorCode PetscOptionsVec(const char key[],const char text[],const char man[],Vec v,PetscBool *set)
1841: {
1842: PetscInt i,N,rstart,rend;
1844: PetscScalar *xx;
1845: PetscReal *xreal;
1846: PetscBool iset;
1849: VecGetOwnershipRange(v,&rstart,&rend);
1850: VecGetSize(v,&N);
1851: PetscCalloc1(N,&xreal);
1852: PetscOptionsRealArray(key,text,man,xreal,&N,&iset);
1853: if (iset) {
1854: VecGetArray(v,&xx);
1855: for (i=rstart; i<rend; i++) xx[i-rstart] = xreal[i];
1856: VecRestoreArray(v,&xx);
1857: }
1858: PetscFree(xreal);
1859: if (set) *set = iset;
1860: return(0);
1861: }
1865: /*@
1866: VecGetLayout - get PetscLayout describing vector layout
1868: Not Collective
1870: Input Arguments:
1871: . x - the vector
1873: Output Arguments:
1874: . map - the layout
1876: Level: developer
1878: .seealso: VecGetSizes(), VecGetOwnershipRange(), VecGetOwnershipRanges()
1879: @*/
1880: PetscErrorCode VecGetLayout(Vec x,PetscLayout *map)
1881: {
1885: *map = x->map;
1886: return(0);
1887: }
1891: /*@
1892: VecSetLayout - set PetscLayout describing vector layout
1894: Not Collective
1896: Input Arguments:
1897: + x - the vector
1898: - map - the layout
1900: Notes:
1901: It is normally only valid to replace the layout with a layout known to be equivalent.
1903: Level: developer
1905: .seealso: VecGetLayout(), VecGetSizes(), VecGetOwnershipRange(), VecGetOwnershipRanges()
1906: @*/
1907: PetscErrorCode VecSetLayout(Vec x,PetscLayout map)
1908: {
1913: PetscLayoutReference(map,&x->map);
1914: return(0);
1915: }