Actual source code: vector.c

petsc-3.9.2 2018-05-20
Report Typos and Errors

  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>

  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;
 16: PetscLogEvent VEC_ViennaCLCopyFromGPU, VEC_ViennaCLCopyToGPU;
 17: PetscLogEvent VEC_CUDACopyFromGPU, VEC_CUDACopyToGPU;
 18: PetscLogEvent VEC_CUDACopyFromGPUSome, VEC_CUDACopyToGPUSome;

 20: /*@
 21:    VecStashGetInfo - Gets how many values are currently in the vector stash, i.e. need
 22:        to be communicated to other processors during the VecAssemblyBegin/End() process

 24:     Not collective

 26:    Input Parameter:
 27: .   vec - the vector

 29:    Output Parameters:
 30: +   nstash   - the size of the stash
 31: .   reallocs - the number of additional mallocs incurred.
 32: .   bnstash   - the size of the block stash
 33: -   breallocs - the number of additional mallocs incurred.in the block stash

 35:    Level: advanced

 37: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), Vec, VecStashSetInitialSize(), VecStashView()

 39: @*/
 40: PetscErrorCode  VecStashGetInfo(Vec vec,PetscInt *nstash,PetscInt *reallocs,PetscInt *bnstash,PetscInt *breallocs)
 41: {

 45:   VecStashGetInfo_Private(&vec->stash,nstash,reallocs);
 46:   VecStashGetInfo_Private(&vec->bstash,bnstash,breallocs);
 47:   return(0);
 48: }

 50: /*@
 51:    VecSetLocalToGlobalMapping - Sets a local numbering to global numbering used
 52:    by the routine VecSetValuesLocal() to allow users to insert vector entries
 53:    using a local (per-processor) numbering.

 55:    Logically Collective on Vec

 57:    Input Parameters:
 58: +  x - vector
 59: -  mapping - mapping created with ISLocalToGlobalMappingCreate() or ISLocalToGlobalMappingCreateIS()

 61:    Notes:
 62:    All vectors obtained with VecDuplicate() from this vector inherit the same mapping.

 64:    Level: intermediate

 66:    Concepts: vector^setting values with local numbering

 68: seealso:  VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetValuesLocal(),
 69:            VecSetLocalToGlobalMapping(), VecSetValuesBlockedLocal()
 70: @*/
 71: PetscErrorCode  VecSetLocalToGlobalMapping(Vec x,ISLocalToGlobalMapping mapping)
 72: {


 79:   if (x->ops->setlocaltoglobalmapping) {
 80:     (*x->ops->setlocaltoglobalmapping)(x,mapping);
 81:   } else {
 82:     PetscLayoutSetISLocalToGlobalMapping(x->map,mapping);
 83:   }
 84:   return(0);
 85: }

 87: /*@
 88:    VecGetLocalToGlobalMapping - Gets the local-to-global numbering set by VecSetLocalToGlobalMapping()

 90:    Not Collective

 92:    Input Parameter:
 93: .  X - the vector

 95:    Output Parameter:
 96: .  mapping - the mapping

 98:    Level: advanced

100:    Concepts: vectors^local to global mapping
101:    Concepts: local to global mapping^for vectors

103: .seealso:  VecSetValuesLocal()
104: @*/
105: PetscErrorCode VecGetLocalToGlobalMapping(Vec X,ISLocalToGlobalMapping *mapping)
106: {
111:   *mapping = X->map->mapping;
112:   return(0);
113: }

115: /*@
116:    VecAssemblyBegin - Begins assembling the vector.  This routine should
117:    be called after completing all calls to VecSetValues().

119:    Collective on Vec

121:    Input Parameter:
122: .  vec - the vector

124:    Level: beginner

126:    Concepts: assembly^vectors

128: .seealso: VecAssemblyEnd(), VecSetValues()
129: @*/
130: PetscErrorCode  VecAssemblyBegin(Vec vec)
131: {

137:   VecStashViewFromOptions(vec,NULL,"-vec_view_stash");
138:   PetscLogEventBegin(VEC_AssemblyBegin,vec,0,0,0);
139:   if (vec->ops->assemblybegin) {
140:     (*vec->ops->assemblybegin)(vec);
141:   }
142:   PetscLogEventEnd(VEC_AssemblyBegin,vec,0,0,0);
143:   PetscObjectStateIncrease((PetscObject)vec);
144:   return(0);
145: }

147: /*@
148:    VecAssemblyEnd - Completes assembling the vector.  This routine should
149:    be called after VecAssemblyBegin().

151:    Collective on Vec

153:    Input Parameter:
154: .  vec - the vector

156:    Options Database Keys:
157: +  -vec_view - Prints vector in ASCII format
158: .  -vec_view ::ascii_matlab - Prints vector in ASCII MATLAB format to stdout
159: .  -vec_view matlab:filename - Prints vector in MATLAB format to matlaboutput.mat
160: .  -vec_view draw - Activates vector viewing using drawing tools
161: .  -display <name> - Sets display name (default is host)
162: .  -draw_pause <sec> - Sets number of seconds to pause after display
163: -  -vec_view socket - Activates vector viewing using a socket

165:    Level: beginner

167: .seealso: VecAssemblyBegin(), VecSetValues()
168: @*/
169: PetscErrorCode  VecAssemblyEnd(Vec vec)
170: {

175:   PetscLogEventBegin(VEC_AssemblyEnd,vec,0,0,0);
177:   if (vec->ops->assemblyend) {
178:     (*vec->ops->assemblyend)(vec);
179:   }
180:   PetscLogEventEnd(VEC_AssemblyEnd,vec,0,0,0);
181:   VecViewFromOptions(vec,NULL,"-vec_view");
182:   return(0);
183: }

185: /*@
186:    VecPointwiseMax - Computes the componentwise maximum w_i = max(x_i, y_i).

188:    Logically Collective on Vec

190:    Input Parameters:
191: .  x, y  - the vectors

193:    Output Parameter:
194: .  w - the result

196:    Level: advanced

198:    Notes: any subset of the x, y, and w may be the same vector.
199:           For complex numbers compares only the real part

201:    Concepts: vector^pointwise multiply

203: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
204: @*/
205: PetscErrorCode  VecPointwiseMax(Vec w,Vec x,Vec y)
206: {

218:   VecCheckSameSize(w,1,x,2);
219:   VecCheckSameSize(w,1,y,3);
220:   (*w->ops->pointwisemax)(w,x,y);
221:   PetscObjectStateIncrease((PetscObject)w);
222:   return(0);
223: }


226: /*@
227:    VecPointwiseMin - Computes the componentwise minimum w_i = min(x_i, y_i).

229:    Logically Collective on Vec

231:    Input Parameters:
232: .  x, y  - the vectors

234:    Output Parameter:
235: .  w - the result

237:    Level: advanced

239:    Notes: any subset of the x, y, and w may be the same vector.
240:           For complex numbers compares only the real part

242:    Concepts: vector^pointwise multiply

244: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
245: @*/
246: PetscErrorCode  VecPointwiseMin(Vec w,Vec x,Vec y)
247: {

259:   VecCheckSameSize(w,1,x,2);
260:   VecCheckSameSize(w,1,y,3);
261:   (*w->ops->pointwisemin)(w,x,y);
262:   PetscObjectStateIncrease((PetscObject)w);
263:   return(0);
264: }

266: /*@
267:    VecPointwiseMaxAbs - Computes the componentwise maximum of the absolute values w_i = max(abs(x_i), abs(y_i)).

269:    Logically Collective on Vec

271:    Input Parameters:
272: .  x, y  - the vectors

274:    Output Parameter:
275: .  w - the result

277:    Level: advanced

279:    Notes: any subset of the x, y, and w may be the same vector.

281:    Concepts: vector^pointwise multiply

283: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMin(), VecPointwiseMax(), VecMaxPointwiseDivide()
284: @*/
285: PetscErrorCode  VecPointwiseMaxAbs(Vec w,Vec x,Vec y)
286: {

298:   VecCheckSameSize(w,1,x,2);
299:   VecCheckSameSize(w,1,y,3);
300:   (*w->ops->pointwisemaxabs)(w,x,y);
301:   PetscObjectStateIncrease((PetscObject)w);
302:   return(0);
303: }

305: /*@
306:    VecPointwiseDivide - Computes the componentwise division w = x/y.

308:    Logically Collective on Vec

310:    Input Parameters:
311: .  x, y  - the vectors

313:    Output Parameter:
314: .  w - the result

316:    Level: advanced

318:    Notes: any subset of the x, y, and w may be the same vector.

320:    Concepts: vector^pointwise divide

322: .seealso: VecPointwiseMult(), VecPointwiseMax(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
323: @*/
324: PetscErrorCode  VecPointwiseDivide(Vec w,Vec x,Vec y)
325: {

337:   VecCheckSameSize(w,1,x,2);
338:   VecCheckSameSize(w,1,y,3);
339:   (*w->ops->pointwisedivide)(w,x,y);
340:   PetscObjectStateIncrease((PetscObject)w);
341:   return(0);
342: }


345: /*@
346:    VecDuplicate - Creates a new vector of the same type as an existing vector.

348:    Collective on Vec

350:    Input Parameters:
351: .  v - a vector to mimic

353:    Output Parameter:
354: .  newv - location to put new vector

356:    Notes:
357:    VecDuplicate() DOES NOT COPY the vector entries, but rather allocates storage
358:    for the new vector.  Use VecCopy() to copy a vector.

360:    Use VecDestroy() to free the space. Use VecDuplicateVecs() to get several
361:    vectors.

363:    Level: beginner

365: .seealso: VecDestroy(), VecDuplicateVecs(), VecCreate(), VecCopy()
366: @*/
367: PetscErrorCode  VecDuplicate(Vec v,Vec *newv)
368: {

375:   (*v->ops->duplicate)(v,newv);
376:   PetscObjectStateIncrease((PetscObject)*newv);
377:   return(0);
378: }

380: /*@
381:    VecDestroy - Destroys a vector.

383:    Collective on Vec

385:    Input Parameters:
386: .  v  - the vector

388:    Level: beginner

390: .seealso: VecDuplicate(), VecDestroyVecs()
391: @*/
392: PetscErrorCode  VecDestroy(Vec *v)
393: {

397:   if (!*v) return(0);
399:   if (--((PetscObject)(*v))->refct > 0) {*v = 0; return(0);}

401:   PetscObjectSAWsViewOff((PetscObject)*v);
402:   /* destroy the internal part */
403:   if ((*v)->ops->destroy) {
404:     (*(*v)->ops->destroy)(*v);
405:   }
406:   /* destroy the external/common part */
407:   PetscLayoutDestroy(&(*v)->map);
408:   PetscHeaderDestroy(v);
409:   return(0);
410: }

412: /*@C
413:    VecDuplicateVecs - Creates several vectors of the same type as an existing vector.

415:    Collective on Vec

417:    Input Parameters:
418: +  m - the number of vectors to obtain
419: -  v - a vector to mimic

421:    Output Parameter:
422: .  V - location to put pointer to array of vectors

424:    Notes:
425:    Use VecDestroyVecs() to free the space. Use VecDuplicate() to form a single
426:    vector.

428:    Fortran Note:
429:    The Fortran interface is slightly different from that given below, it
430:    requires one to pass in V a Vec (integer) array of size at least m.
431:    See the Fortran chapter of the users manual and petsc/src/vec/vec/examples for details.

433:    Level: intermediate

435: .seealso:  VecDestroyVecs(), VecDuplicate(), VecCreate(), VecDuplicateVecsF90()
436: @*/
437: PetscErrorCode  VecDuplicateVecs(Vec v,PetscInt m,Vec *V[])
438: {

445:   (*v->ops->duplicatevecs)(v,m,V);
446:   return(0);
447: }

449: /*@C
450:    VecDestroyVecs - Frees a block of vectors obtained with VecDuplicateVecs().

452:    Collective on Vec

454:    Input Parameters:
455: +  vv - pointer to pointer to array of vector pointers, if NULL no vectors are destroyed
456: -  m - the number of vectors previously obtained, if zero no vectors are destroyed

458:    Fortran Note:
459:    The Fortran interface is slightly different from that given below.
460:    See the Fortran chapter of the users manual

462:    Level: intermediate

464: .seealso: VecDuplicateVecs(), VecDestroyVecsf90()
465: @*/
466: PetscErrorCode  VecDestroyVecs(PetscInt m,Vec *vv[])
467: {

472:   if (m < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Trying to destroy negative number of vectors %D",m);
473:   if (!m || !*vv) {*vv  = NULL; return(0);}
476:   (*(**vv)->ops->destroyvecs)(m,*vv);
477:   *vv  = NULL;
478:   return(0);
479: }

481: /*@C
482:    VecView - Views a vector object.

484:    Collective on Vec

486:    Input Parameters:
487: +  vec - the vector
488: -  viewer - an optional visualization context

490:    Notes:
491:    The available visualization contexts include
492: +     PETSC_VIEWER_STDOUT_SELF - for sequential vectors
493: .     PETSC_VIEWER_STDOUT_WORLD - for parallel vectors created on PETSC_COMM_WORLD
494: -     PETSC_VIEWER_STDOUT_(comm) - for parallel vectors created on MPI communicator comm

496:    You can change the format the vector is printed using the
497:    option PetscViewerPushFormat().

499:    The user can open alternative visualization contexts with
500: +    PetscViewerASCIIOpen() - Outputs vector to a specified file
501: .    PetscViewerBinaryOpen() - Outputs vector in binary to a
502:          specified file; corresponding input uses VecLoad()
503: .    PetscViewerDrawOpen() - Outputs vector to an X window display
504: -    PetscViewerSocketOpen() - Outputs vector to Socket viewer

506:    The user can call PetscViewerPushFormat() to specify the output
507:    format of ASCII printed objects (when using PETSC_VIEWER_STDOUT_SELF,
508:    PETSC_VIEWER_STDOUT_WORLD and PetscViewerASCIIOpen).  Available formats include
509: +    PETSC_VIEWER_DEFAULT - default, prints vector contents
510: .    PETSC_VIEWER_ASCII_MATLAB - prints vector contents in MATLAB format
511: .    PETSC_VIEWER_ASCII_INDEX - prints vector contents, including indices of vector elements
512: -    PETSC_VIEWER_ASCII_COMMON - prints vector contents, using a
513:          format common among all vector types

515:    Notes: You can pass any number of vector objects, or other PETSc objects to the same viewer.

517:    Notes for binary viewer: If you pass multiply vectors to a binary viewer you can read them back in in the same order
518: $     with VecLoad().
519: $
520: $    If the blocksize of the vector is greater than one then you must provide a unique prefix to
521: $    the vector with PetscObjectSetOptionsPrefix((PetscObject)vec,"uniqueprefix"); BEFORE calling VecView() on the
522: $    vector to be stored and then set that same unique prefix on the vector that you pass to VecLoad(). The blocksize
523: $    information is stored in an ASCII file with the same name as the binary file plus a ".info" appended to the
524: $    filename. If you copy the binary file, make sure you copy the associated .info file with it.

526:    Notes for HDF5 Viewer: the name of the Vec (given with PetscObjectSetName() is the name that is used
527: $    for the object in the HDF5 file. If you wish to store the same vector to the HDF5 viewer (with different values,
528: $    obviously) several times, you must change its name each time before calling the VecView(). The name you use
529: $    here should equal the name that you use in the Vec object that you use with VecLoad().

531:    See the manual page for VecLoad() on the exact format the binary viewer stores
532:    the values in the file.

534:    Level: beginner

536:    Concepts: vector^printing
537:    Concepts: vector^saving to disk

539: .seealso: PetscViewerASCIIOpen(), PetscViewerDrawOpen(), PetscDrawLGCreate(),
540:           PetscViewerSocketOpen(), PetscViewerBinaryOpen(), VecLoad(), PetscViewerCreate(),
541:           PetscRealView(), PetscScalarView(), PetscIntView()
542: @*/
543: PetscErrorCode  VecView(Vec vec,PetscViewer viewer)
544: {
545:   PetscErrorCode    ierr;
546:   PetscBool         iascii;
547:   PetscViewerFormat format;
548:   PetscMPIInt       size;

553:   if (!viewer) {
554:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)vec),&viewer);
555:   }
557:   PetscViewerGetFormat(viewer,&format);
558:   MPI_Comm_size(PetscObjectComm((PetscObject)vec),&size);
559:   if (size == 1 && format == PETSC_VIEWER_LOAD_BALANCE) return(0);

561:   if (vec->stash.n || vec->bstash.n) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call VecAssemblyBegin/End() before viewing this vector");

563:   PetscLogEventBegin(VEC_View,vec,viewer,0,0);
564:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
565:   if (iascii) {
566:     PetscInt rows,bs;

568:     PetscObjectPrintClassNamePrefixType((PetscObject)vec,viewer);
569:     if (format == PETSC_VIEWER_ASCII_INFO || format == PETSC_VIEWER_ASCII_INFO_DETAIL) {
570:       PetscViewerASCIIPushTab(viewer);
571:       VecGetSize(vec,&rows);
572:       VecGetBlockSize(vec,&bs);
573:       if (bs != 1) {
574:         PetscViewerASCIIPrintf(viewer,"length=%D, bs=%D\n",rows,bs);
575:       } else {
576:         PetscViewerASCIIPrintf(viewer,"length=%D\n",rows);
577:       }
578:       PetscViewerASCIIPopTab(viewer);
579:     }
580:   }
581:   VecLockPush(vec);
582:   if ((format == PETSC_VIEWER_NATIVE || format == PETSC_VIEWER_LOAD_BALANCE) && vec->ops->viewnative) {
583:     (*vec->ops->viewnative)(vec,viewer);
584:   } else {
585:     (*vec->ops->view)(vec,viewer);
586:   }
587:   VecLockPop(vec);
588:   PetscLogEventEnd(VEC_View,vec,viewer,0,0);
589:   return(0);
590: }

592: #if defined(PETSC_USE_DEBUG)
593: #include <../src/sys/totalview/tv_data_display.h>
594: PETSC_UNUSED static int TV_display_type(const struct _p_Vec *v)
595: {
596:   const PetscScalar *values;
597:   char              type[32];
598:   PetscErrorCode    ierr;


601:   TV_add_row("Local rows", "int", &v->map->n);
602:   TV_add_row("Global rows", "int", &v->map->N);
603:   TV_add_row("Typename", TV_ascii_string_type, ((PetscObject)v)->type_name);
604:   VecGetArrayRead((Vec)v,&values);
605:   PetscSNPrintf(type,32,"double[%d]",v->map->n);
606:   TV_add_row("values",type, values);
607:   VecRestoreArrayRead((Vec)v,&values);
608:   return TV_format_OK;
609: }
610: #endif

612: /*@
613:    VecGetSize - Returns the global number of elements of the vector.

615:    Not Collective

617:    Input Parameter:
618: .  x - the vector

620:    Output Parameters:
621: .  size - the global length of the vector

623:    Level: beginner

625:    Concepts: vector^local size

627: .seealso: VecGetLocalSize()
628: @*/
629: PetscErrorCode  VecGetSize(Vec x,PetscInt *size)
630: {

637:   (*x->ops->getsize)(x,size);
638:   return(0);
639: }

641: /*@
642:    VecGetLocalSize - Returns the number of elements of the vector stored
643:    in local memory. This routine may be implementation dependent, so use
644:    with care.

646:    Not Collective

648:    Input Parameter:
649: .  x - the vector

651:    Output Parameter:
652: .  size - the length of the local piece of the vector

654:    Level: beginner

656:    Concepts: vector^size

658: .seealso: VecGetSize()
659: @*/
660: PetscErrorCode  VecGetLocalSize(Vec x,PetscInt *size)
661: {

668:   (*x->ops->getlocalsize)(x,size);
669:   return(0);
670: }

672: /*@C
673:    VecGetOwnershipRange - Returns the range of indices owned by
674:    this processor, assuming that the vectors are laid out with the
675:    first n1 elements on the first processor, next n2 elements on the
676:    second, etc.  For certain parallel layouts this range may not be
677:    well defined.

679:    Not Collective

681:    Input Parameter:
682: .  x - the vector

684:    Output Parameters:
685: +  low - the first local element, pass in NULL if not interested
686: -  high - one more than the last local element, pass in NULL if not interested

688:    Note:
689:    The high argument is one more than the last element stored locally.

691:    Fortran: NULL_INTEGER should be used instead of NULL

693:    Level: beginner

695:    Concepts: ownership^of vectors
696:    Concepts: vector^ownership of elements

698: .seealso:   MatGetOwnershipRange(), MatGetOwnershipRanges(), VecGetOwnershipRanges()
699: @*/
700: PetscErrorCode  VecGetOwnershipRange(Vec x,PetscInt *low,PetscInt *high)
701: {
707:   if (low)  *low  = x->map->rstart;
708:   if (high) *high = x->map->rend;
709:   return(0);
710: }

712: /*@C
713:    VecGetOwnershipRanges - Returns the range of indices owned by EACH processor,
714:    assuming that the vectors are laid out with the
715:    first n1 elements on the first processor, next n2 elements on the
716:    second, etc.  For certain parallel layouts this range may not be
717:    well defined.

719:    Not Collective

721:    Input Parameter:
722: .  x - the vector

724:    Output Parameters:
725: .  range - array of length size+1 with the start and end+1 for each process

727:    Note:
728:    The high argument is one more than the last element stored locally.

730:    Fortran: You must PASS in an array of length size+1

732:    Level: beginner

734:    Concepts: ownership^of vectors
735:    Concepts: vector^ownership of elements

737: .seealso:   MatGetOwnershipRange(), MatGetOwnershipRanges(), VecGetOwnershipRange()
738: @*/
739: PetscErrorCode  VecGetOwnershipRanges(Vec x,const PetscInt *ranges[])
740: {

746:   PetscLayoutGetRanges(x->map,ranges);
747:   return(0);
748: }

750: /*@
751:    VecSetOption - Sets an option for controling a vector's behavior.

753:    Collective on Vec

755:    Input Parameter:
756: +  x - the vector
757: .  op - the option
758: -  flag - turn the option on or off

760:    Supported Options:
761: +     VEC_IGNORE_OFF_PROC_ENTRIES, which causes VecSetValues() to ignore
762:           entries destined to be stored on a separate processor. This can be used
763:           to eliminate the global reduction in the VecAssemblyXXXX() if you know
764:           that you have only used VecSetValues() to set local elements
765: .     VEC_IGNORE_NEGATIVE_INDICES, which means you can pass negative indices
766:           in ix in calls to VecSetValues() or VecGetValues(). These rows are simply
767:           ignored.
768: -     VEC_SUBSET_OFF_PROC_ENTRIES, which causes VecAssemblyBegin() to assume that the off-process
769:           entries will always be a subset (possibly equal) of the off-process entries set on the
770:           first assembly.  This reuses the communication pattern, thus avoiding a global reduction.
771:           Subsequent assemblies setting off-process values should use the same InsertMode as the
772:           first assembly.

774:    Developer Note:
775:    The InsertMode restriction could be removed by packing the stash messages out of place.

777:    Level: intermediate

779: @*/
780: PetscErrorCode  VecSetOption(Vec x,VecOption op,PetscBool flag)
781: {

787:   if (x->ops->setoption) {
788:     (*x->ops->setoption)(x,op,flag);
789:   }
790:   return(0);
791: }

793: /* Default routines for obtaining and releasing; */
794: /* may be used by any implementation */
795: PetscErrorCode VecDuplicateVecs_Default(Vec w,PetscInt m,Vec *V[])
796: {
798:   PetscInt       i;

803:   if (m <= 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"m must be > 0: m = %D",m);
804:   PetscMalloc1(m,V);
805:   for (i=0; i<m; i++) {VecDuplicate(w,*V+i);}
806:   return(0);
807: }

809: PetscErrorCode VecDestroyVecs_Default(PetscInt m,Vec v[])
810: {
812:   PetscInt       i;

816:   for (i=0; i<m; i++) {VecDestroy(&v[i]);}
817:   PetscFree(v);
818:   return(0);
819: }

821: /*@
822:    VecResetArray - Resets a vector to use its default memory. Call this
823:    after the use of VecPlaceArray().

825:    Not Collective

827:    Input Parameters:
828: .  vec - the vector

830:    Level: developer

832: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecPlaceArray()

834: @*/
835: PetscErrorCode  VecResetArray(Vec vec)
836: {

842:   if (vec->ops->resetarray) {
843:     (*vec->ops->resetarray)(vec);
844:   } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Cannot reset array in this type of vector");
845:   PetscObjectStateIncrease((PetscObject)vec);
846:   return(0);
847: }

849: /*@C
850:   VecLoad - Loads a vector that has been stored in binary or HDF5 format
851:   with VecView().

853:   Collective on PetscViewer

855:   Input Parameters:
856: + newvec - the newly loaded vector, this needs to have been created with VecCreate() or
857:            some related function before a call to VecLoad().
858: - viewer - binary file viewer, obtained from PetscViewerBinaryOpen() or
859:            HDF5 file viewer, obtained from PetscViewerHDF5Open()

861:    Level: intermediate

863:   Notes:
864:   Defaults to the standard Seq or MPI Vec, if you want some other type of Vec call VecSetFromOptions()
865:   before calling this.

867:   The input file must contain the full global vector, as
868:   written by the routine VecView().

870:   If the type or size of newvec is not set before a call to VecLoad, PETSc
871:   sets the type and the local and global sizes. If type and/or
872:   sizes are already set, then the same are used.

874:   If using binary and the blocksize of the vector is greater than one then you must provide a unique prefix to
875:   the vector with PetscObjectSetOptionsPrefix((PetscObject)vec,"uniqueprefix"); BEFORE calling VecView() on the
876:   vector to be stored and then set that same unique prefix on the vector that you pass to VecLoad(). The blocksize
877:   information is stored in an ASCII file with the same name as the binary file plus a ".info" appended to the
878:   filename. If you copy the binary file, make sure you copy the associated .info file with it.

880:   If using HDF5, you must assign the Vec the same name as was used in the Vec
881:   that was stored in the file using PetscObjectSetName(). Otherwise you will
882:   get the error message: "Cannot H5DOpen2() with Vec name NAMEOFOBJECT"

884:   Notes for advanced users:
885:   Most users should not need to know the details of the binary storage
886:   format, since VecLoad() and VecView() completely hide these details.
887:   But for anyone who's interested, the standard binary vector storage
888:   format is
889: .vb
890:      int    VEC_FILE_CLASSID
891:      int    number of rows
892:      PetscScalar *values of all entries
893: .ve

895:    In addition, PETSc automatically does the byte swapping for
896: machines that store the bytes reversed, e.g.  DEC alpha, freebsd,
897: linux, Windows and the paragon; thus if you write your own binary
898: read/write routines you have to swap the bytes; see PetscBinaryRead()
899: and PetscBinaryWrite() to see how this may be done.

901:   Concepts: vector^loading from file

903: .seealso: PetscViewerBinaryOpen(), VecView(), MatLoad(), VecLoad()
904: @*/
905: PetscErrorCode  VecLoad(Vec newvec, PetscViewer viewer)
906: {
908:   PetscBool      isbinary,ishdf5;
909:   PetscViewerFormat format;

914:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
915:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5);
916:   if (!isbinary && !ishdf5) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Invalid viewer; open viewer with PetscViewerBinaryOpen()");

918:   PetscLogEventBegin(VEC_Load,viewer,0,0,0);
919:   if (!((PetscObject)newvec)->type_name && !newvec->ops->create) {
920:     VecSetType(newvec, VECSTANDARD);
921:   }
922:   PetscViewerGetFormat(viewer,&format);
923:   if (format == PETSC_VIEWER_NATIVE && newvec->ops->loadnative) {
924:     (*newvec->ops->loadnative)(newvec,viewer);
925:   } else {
926:     (*newvec->ops->load)(newvec,viewer);
927:   }
928:   PetscLogEventEnd(VEC_Load,viewer,0,0,0);
929:   return(0);
930: }


933: /*@
934:    VecReciprocal - Replaces each component of a vector by its reciprocal.

936:    Logically Collective on Vec

938:    Input Parameter:
939: .  vec - the vector

941:    Output Parameter:
942: .  vec - the vector reciprocal

944:    Level: intermediate

946:    Concepts: vector^reciprocal

948: .seealso: VecLog(), VecExp(), VecSqrtAbs()

950: @*/
951: PetscErrorCode  VecReciprocal(Vec vec)
952: {

958:   if (vec->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
959:   if (!vec->ops->reciprocal) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Vector does not support reciprocal operation");
960:   (*vec->ops->reciprocal)(vec);
961:   PetscObjectStateIncrease((PetscObject)vec);
962:   return(0);
963: }

965: /*@C
966:     VecSetOperation - Allows user to set a vector operation.

968:    Logically Collective on Vec

970:     Input Parameters:
971: +   vec - the vector
972: .   op - the name of the operation
973: -   f - the function that provides the operation.

975:    Level: advanced

977:     Usage:
978: $      PetscErrorCode userview(Vec,PetscViewer);
979: $      VecCreateMPI(comm,m,M,&x);
980: $      VecSetOperation(x,VECOP_VIEW,(void(*)(void))userview);

982:     Notes:
983:     See the file include/petscvec.h for a complete list of matrix
984:     operations, which all have the form VECOP_<OPERATION>, where
985:     <OPERATION> is the name (in all capital letters) of the
986:     user interface routine (e.g., VecView() -> VECOP_VIEW).

988:     This function is not currently available from Fortran.

990: .keywords: vector, set, operation

992: .seealso: VecCreate(), MatShellSetOperation()
993: @*/
994: PetscErrorCode VecSetOperation(Vec vec,VecOperation op, void (*f)(void))
995: {
998:   if (op == VECOP_VIEW && !vec->ops->viewnative) {
999:     vec->ops->viewnative = vec->ops->view;
1000:   } else if (op == VECOP_LOAD && !vec->ops->loadnative) {
1001:     vec->ops->loadnative = vec->ops->load;
1002:   }
1003:   (((void(**)(void))vec->ops)[(int)op]) = f;
1004:   return(0);
1005: }


1008: /*@
1009:    VecStashSetInitialSize - sets the sizes of the vec-stash, that is
1010:    used during the assembly process to store values that belong to
1011:    other processors.

1013:    Not Collective, different processes can have different size stashes

1015:    Input Parameters:
1016: +  vec   - the vector
1017: .  size  - the initial size of the stash.
1018: -  bsize - the initial size of the block-stash(if used).

1020:    Options Database Keys:
1021: +   -vecstash_initial_size <size> or <size0,size1,...sizep-1>
1022: -   -vecstash_block_initial_size <bsize> or <bsize0,bsize1,...bsizep-1>

1024:    Level: intermediate

1026:    Notes:
1027:      The block-stash is used for values set with VecSetValuesBlocked() while
1028:      the stash is used for values set with VecSetValues()

1030:      Run with the option -info and look for output of the form
1031:      VecAssemblyBegin_MPIXXX:Stash has MM entries, uses nn mallocs.
1032:      to determine the appropriate value, MM, to use for size and
1033:      VecAssemblyBegin_MPIXXX:Block-Stash has BMM entries, uses nn mallocs.
1034:      to determine the value, BMM to use for bsize

1036:    Concepts: vector^stash
1037:    Concepts: stash^vector

1039: .seealso: VecSetBlockSize(), VecSetValues(), VecSetValuesBlocked(), VecStashView()

1041: @*/
1042: PetscErrorCode  VecStashSetInitialSize(Vec vec,PetscInt size,PetscInt bsize)
1043: {

1048:   VecStashSetInitialSize_Private(&vec->stash,size);
1049:   VecStashSetInitialSize_Private(&vec->bstash,bsize);
1050:   return(0);
1051: }

1053: /*@
1054:    VecConjugate - Conjugates a vector.

1056:    Logically Collective on Vec

1058:    Input Parameters:
1059: .  x - the vector

1061:    Level: intermediate

1063:    Concepts: vector^conjugate

1065: @*/
1066: PetscErrorCode  VecConjugate(Vec x)
1067: {
1068: #if defined(PETSC_USE_COMPLEX)

1074:   if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1075:   (*x->ops->conjugate)(x);
1076:   /* we need to copy norms here */
1077:   PetscObjectStateIncrease((PetscObject)x);
1078:   return(0);
1079: #else
1080:   return(0);
1081: #endif
1082: }

1084: /*@
1085:    VecPointwiseMult - Computes the componentwise multiplication w = x*y.

1087:    Logically Collective on Vec

1089:    Input Parameters:
1090: .  x, y  - the vectors

1092:    Output Parameter:
1093: .  w - the result

1095:    Level: advanced

1097:    Notes: any subset of the x, y, and w may be the same vector.

1099:    Concepts: vector^pointwise multiply

1101: .seealso: VecPointwiseDivide(), VecPointwiseMax(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
1102: @*/
1103: PetscErrorCode  VecPointwiseMult(Vec w, Vec x,Vec y)
1104: {

1116:   VecCheckSameSize(w,1,x,2);
1117:   VecCheckSameSize(w,2,y,3);
1118:   PetscLogEventBegin(VEC_PointwiseMult,x,y,w,0);
1119:   (*w->ops->pointwisemult)(w,x,y);
1120:   PetscLogEventEnd(VEC_PointwiseMult,x,y,w,0);
1121:   PetscObjectStateIncrease((PetscObject)w);
1122:   return(0);
1123: }

1125: /*@C
1126:    VecSetRandom - Sets all components of a vector to random numbers.

1128:    Logically Collective on Vec

1130:    Input Parameters:
1131: +  x  - the vector
1132: -  rctx - the random number context, formed by PetscRandomCreate(), or NULL and
1133:           it will create one internally.

1135:    Output Parameter:
1136: .  x  - the vector

1138:    Example of Usage:
1139: .vb
1140:      PetscRandomCreate(PETSC_COMM_WORLD,&rctx);
1141:      VecSetRandom(x,rctx);
1142:      PetscRandomDestroy(rctx);
1143: .ve

1145:    Level: intermediate

1147:    Concepts: vector^setting to random
1148:    Concepts: random^vector

1150: .seealso: VecSet(), VecSetValues(), PetscRandomCreate(), PetscRandomDestroy()
1151: @*/
1152: PetscErrorCode  VecSetRandom(Vec x,PetscRandom rctx)
1153: {
1155:   PetscRandom    randObj = NULL;

1161:   if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");

1163:   if (!rctx) {
1164:     MPI_Comm comm;
1165:     PetscObjectGetComm((PetscObject)x,&comm);
1166:     PetscRandomCreate(comm,&randObj);
1167:     PetscRandomSetFromOptions(randObj);
1168:     rctx = randObj;
1169:   }

1171:   PetscLogEventBegin(VEC_SetRandom,x,rctx,0,0);
1172:   (*x->ops->setrandom)(x,rctx);
1173:   PetscLogEventEnd(VEC_SetRandom,x,rctx,0,0);

1175:   PetscRandomDestroy(&randObj);
1176:   PetscObjectStateIncrease((PetscObject)x);
1177:   return(0);
1178: }

1180: /*@
1181:   VecZeroEntries - puts a 0.0 in each element of a vector

1183:   Logically Collective on Vec

1185:   Input Parameter:
1186: . vec - The vector

1188:   Level: beginner

1190:   Developer Note: This routine does not need to exist since the exact functionality is obtained with
1191:      VecSet(vec,0);  I guess someone added it to mirror the functionality of MatZeroEntries() but Mat is nothing
1192:      like a Vec (one is an operator and one is an element of a vector space, yeah yeah dual blah blah blah) so
1193:      this routine should not exist.

1195: .keywords: Vec, set, options, database
1196: .seealso: VecCreate(),  VecSetOptionsPrefix(), VecSet(), VecSetValues()
1197: @*/
1198: PetscErrorCode  VecZeroEntries(Vec vec)
1199: {

1203:   VecSet(vec,0);
1204:   return(0);
1205: }

1207: /*
1208:   VecSetTypeFromOptions_Private - Sets the type of vector from user options. Defaults to a PETSc sequential vector on one
1209:   processor and a PETSc MPI vector on more than one processor.

1211:   Collective on Vec

1213:   Input Parameter:
1214: . vec - The vector

1216:   Level: intermediate

1218: .keywords: Vec, set, options, database, type
1219: .seealso: VecSetFromOptions(), VecSetType()
1220: */
1221: static PetscErrorCode VecSetTypeFromOptions_Private(PetscOptionItems *PetscOptionsObject,Vec vec)
1222: {
1223:   PetscBool      opt;
1224:   VecType        defaultType;
1225:   char           typeName[256];
1226:   PetscMPIInt    size;

1230:   if (((PetscObject)vec)->type_name) defaultType = ((PetscObject)vec)->type_name;
1231:   else {
1232:     MPI_Comm_size(PetscObjectComm((PetscObject)vec), &size);
1233:     if (size > 1) defaultType = VECMPI;
1234:     else defaultType = VECSEQ;
1235:   }

1237:   VecRegisterAll();
1238:   PetscOptionsFList("-vec_type","Vector type","VecSetType",VecList,defaultType,typeName,256,&opt);
1239:   if (opt) {
1240:     VecSetType(vec, typeName);
1241:   } else {
1242:     VecSetType(vec, defaultType);
1243:   }
1244:   return(0);
1245: }

1247: /*@
1248:   VecSetFromOptions - Configures the vector from the options database.

1250:   Collective on Vec

1252:   Input Parameter:
1253: . vec - The vector

1255:   Notes:  To see all options, run your program with the -help option, or consult the users manual.
1256:           Must be called after VecCreate() but before the vector is used.

1258:   Level: beginner

1260:   Concepts: vectors^setting options
1261:   Concepts: vectors^setting type

1263: .keywords: Vec, set, options, database
1264: .seealso: VecCreate(), VecSetOptionsPrefix()
1265: @*/
1266: PetscErrorCode  VecSetFromOptions(Vec vec)
1267: {


1273:   PetscObjectOptionsBegin((PetscObject)vec);
1274:   /* Handle vector type options */
1275:   VecSetTypeFromOptions_Private(PetscOptionsObject,vec);

1277:   /* Handle specific vector options */
1278:   if (vec->ops->setfromoptions) {
1279:     (*vec->ops->setfromoptions)(PetscOptionsObject,vec);
1280:   }

1282:   /* process any options handlers added with PetscObjectAddOptionsHandler() */
1283:   PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject)vec);
1284:   PetscOptionsEnd();
1285:   return(0);
1286: }

1288: /*@
1289:   VecSetSizes - Sets the local and global sizes, and checks to determine compatibility

1291:   Collective on Vec

1293:   Input Parameters:
1294: + v - the vector
1295: . n - the local size (or PETSC_DECIDE to have it set)
1296: - N - the global size (or PETSC_DECIDE)

1298:   Notes:
1299:   n and N cannot be both PETSC_DECIDE
1300:   If one processor calls this with N of PETSC_DECIDE then all processors must, otherwise the program will hang.

1302:   Level: intermediate

1304: .seealso: VecGetSize(), PetscSplitOwnership()
1305: @*/
1306: PetscErrorCode  VecSetSizes(Vec v, PetscInt n, PetscInt N)
1307: {

1313:   if (N >= 0 && n > N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Local size %D cannot be larger than global size %D",n,N);
1314:   if ((v->map->n >= 0 || v->map->N >= 0) && (v->map->n != n || v->map->N != N)) SETERRQ4(PETSC_COMM_SELF,PETSC_ERR_SUP,"Cannot change/reset vector sizes to %D local %D global after previously setting them to %D local %D global",n,N,v->map->n,v->map->N);
1315:   v->map->n = n;
1316:   v->map->N = N;
1317:   if (v->ops->create) {
1318:     (*v->ops->create)(v);
1319:     v->ops->create = 0;
1320:   }
1321:   return(0);
1322: }

1324: /*@
1325:    VecSetBlockSize - Sets the blocksize for future calls to VecSetValuesBlocked()
1326:    and VecSetValuesBlockedLocal().

1328:    Logically Collective on Vec

1330:    Input Parameter:
1331: +  v - the vector
1332: -  bs - the blocksize

1334:    Notes:
1335:    All vectors obtained by VecDuplicate() inherit the same blocksize.

1337:    Level: advanced

1339: .seealso: VecSetValuesBlocked(), VecSetLocalToGlobalMapping(), VecGetBlockSize()

1341:   Concepts: block size^vectors
1342: @*/
1343: PetscErrorCode  VecSetBlockSize(Vec v,PetscInt bs)
1344: {

1349:   if (bs < 0 || bs == v->map->bs) return(0);
1351:   PetscLayoutSetBlockSize(v->map,bs);
1352:   v->bstash.bs = bs; /* use the same blocksize for the vec's block-stash */
1353:   return(0);
1354: }

1356: /*@
1357:    VecGetBlockSize - Gets the blocksize for the vector, i.e. what is used for VecSetValuesBlocked()
1358:    and VecSetValuesBlockedLocal().

1360:    Not Collective

1362:    Input Parameter:
1363: .  v - the vector

1365:    Output Parameter:
1366: .  bs - the blocksize

1368:    Notes:
1369:    All vectors obtained by VecDuplicate() inherit the same blocksize.

1371:    Level: advanced

1373: .seealso: VecSetValuesBlocked(), VecSetLocalToGlobalMapping(), VecSetBlockSize()

1375:    Concepts: vector^block size
1376:    Concepts: block^vector

1378: @*/
1379: PetscErrorCode  VecGetBlockSize(Vec v,PetscInt *bs)
1380: {

1386:   PetscLayoutGetBlockSize(v->map,bs);
1387:   return(0);
1388: }

1390: /*@C
1391:    VecSetOptionsPrefix - Sets the prefix used for searching for all
1392:    Vec options in the database.

1394:    Logically Collective on Vec

1396:    Input Parameter:
1397: +  v - the Vec context
1398: -  prefix - the prefix to prepend to all option names

1400:    Notes:
1401:    A hyphen (-) must NOT be given at the beginning of the prefix name.
1402:    The first character of all runtime options is AUTOMATICALLY the hyphen.

1404:    Level: advanced

1406: .keywords: Vec, set, options, prefix, database

1408: .seealso: VecSetFromOptions()
1409: @*/
1410: PetscErrorCode  VecSetOptionsPrefix(Vec v,const char prefix[])
1411: {

1416:   PetscObjectSetOptionsPrefix((PetscObject)v,prefix);
1417:   return(0);
1418: }

1420: /*@C
1421:    VecAppendOptionsPrefix - Appends to the prefix used for searching for all
1422:    Vec options in the database.

1424:    Logically Collective on Vec

1426:    Input Parameters:
1427: +  v - the Vec context
1428: -  prefix - the prefix to prepend to all option names

1430:    Notes:
1431:    A hyphen (-) must NOT be given at the beginning of the prefix name.
1432:    The first character of all runtime options is AUTOMATICALLY the hyphen.

1434:    Level: advanced

1436: .keywords: Vec, append, options, prefix, database

1438: .seealso: VecGetOptionsPrefix()
1439: @*/
1440: PetscErrorCode  VecAppendOptionsPrefix(Vec v,const char prefix[])
1441: {

1446:   PetscObjectAppendOptionsPrefix((PetscObject)v,prefix);
1447:   return(0);
1448: }

1450: /*@C
1451:    VecGetOptionsPrefix - Sets the prefix used for searching for all
1452:    Vec options in the database.

1454:    Not Collective

1456:    Input Parameter:
1457: .  v - the Vec context

1459:    Output Parameter:
1460: .  prefix - pointer to the prefix string used

1462:    Notes: On the fortran side, the user should pass in a string 'prefix' of
1463:    sufficient length to hold the prefix.

1465:    Level: advanced

1467: .keywords: Vec, get, options, prefix, database

1469: .seealso: VecAppendOptionsPrefix()
1470: @*/
1471: PetscErrorCode  VecGetOptionsPrefix(Vec v,const char *prefix[])
1472: {

1477:   PetscObjectGetOptionsPrefix((PetscObject)v,prefix);
1478:   return(0);
1479: }

1481: /*@
1482:    VecSetUp - Sets up the internal vector data structures for the later use.

1484:    Collective on Vec

1486:    Input Parameters:
1487: .  v - the Vec context

1489:    Notes:
1490:    For basic use of the Vec classes the user need not explicitly call
1491:    VecSetUp(), since these actions will happen automatically.

1493:    Level: advanced

1495: .keywords: Vec, setup

1497: .seealso: VecCreate(), VecDestroy()
1498: @*/
1499: PetscErrorCode  VecSetUp(Vec v)
1500: {
1501:   PetscMPIInt    size;

1506:   if (!((PetscObject)v)->type_name) {
1507:     MPI_Comm_size(PetscObjectComm((PetscObject)v), &size);
1508:     if (size == 1) {
1509:       VecSetType(v, VECSEQ);
1510:     } else {
1511:       VecSetType(v, VECMPI);
1512:     }
1513:   }
1514:   return(0);
1515: }

1517: /*
1518:     These currently expose the PetscScalar/PetscReal in updating the
1519:     cached norm. If we push those down into the implementation these
1520:     will become independent of PetscScalar/PetscReal
1521: */

1523: /*@
1524:    VecCopy - Copies a vector. y <- x

1526:    Logically Collective on Vec

1528:    Input Parameter:
1529: .  x - the vector

1531:    Output Parameter:
1532: .  y - the copy

1534:    Notes:
1535:    For default parallel PETSc vectors, both x and y must be distributed in
1536:    the same manner; local copies are done.

1538:    Developer Notes:
1540:    of the vectors to be sequential and one to be parallel so long as both have the same
1541:    local sizes. This is used in some internal functions in PETSc.

1543:    Level: beginner

1545: .seealso: VecDuplicate()
1546: @*/
1547: PetscErrorCode  VecCopy(Vec x,Vec y)
1548: {
1549:   PetscBool      flgs[4];
1550:   PetscReal      norms[4] = {0.0,0.0,0.0,0.0};
1552:   PetscInt       i;

1559:   if (x == y) return(0);
1560:   VecCheckSameLocalSize(x,1,y,2);
1561:   if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1562:   VecLocked(y,2);

1564: #if !defined(PETSC_USE_MIXED_PRECISION)
1565:   for (i=0; i<4; i++) {
1566:     PetscObjectComposedDataGetReal((PetscObject)x,NormIds[i],norms[i],flgs[i]);
1567:   }
1568: #endif

1570:   PetscLogEventBegin(VEC_Copy,x,y,0,0);
1571: #if defined(PETSC_USE_MIXED_PRECISION)
1572:   extern PetscErrorCode VecGetArray(Vec,double**);
1573:   extern PetscErrorCode VecRestoreArray(Vec,double**);
1574:   extern PetscErrorCode VecGetArray(Vec,float**);
1575:   extern PetscErrorCode VecRestoreArray(Vec,float**);
1576:   extern PetscErrorCode VecGetArrayRead(Vec,const double**);
1577:   extern PetscErrorCode VecRestoreArrayRead(Vec,const double**);
1578:   extern PetscErrorCode VecGetArrayRead(Vec,const float**);
1579:   extern PetscErrorCode VecRestoreArrayRead(Vec,const float**);
1580:   if ((((PetscObject)x)->precision == PETSC_PRECISION_SINGLE) && (((PetscObject)y)->precision == PETSC_PRECISION_DOUBLE)) {
1581:     PetscInt    i,n;
1582:     const float *xx;
1583:     double      *yy;
1584:     VecGetArrayRead(x,&xx);
1585:     VecGetArray(y,&yy);
1586:     VecGetLocalSize(x,&n);
1587:     for (i=0; i<n; i++) yy[i] = xx[i];
1588:     VecRestoreArrayRead(x,&xx);
1589:     VecRestoreArray(y,&yy);
1590:   } else if ((((PetscObject)x)->precision == PETSC_PRECISION_DOUBLE) && (((PetscObject)y)->precision == PETSC_PRECISION_SINGLE)) {
1591:     PetscInt     i,n;
1592:     float        *yy;
1593:     const double *xx;
1594:     VecGetArrayRead(x,&xx);
1595:     VecGetArray(y,&yy);
1596:     VecGetLocalSize(x,&n);
1597:     for (i=0; i<n; i++) yy[i] = (float) xx[i];
1598:     VecRestoreArrayRead(x,&xx);
1599:     VecRestoreArray(y,&yy);
1600:   } else {
1601:     (*x->ops->copy)(x,y);
1602:   }
1603: #else
1604:   (*x->ops->copy)(x,y);
1605: #endif

1607:   PetscObjectStateIncrease((PetscObject)y);
1608: #if !defined(PETSC_USE_MIXED_PRECISION)
1609:   for (i=0; i<4; i++) {
1610:     if (flgs[i]) {
1611:       PetscObjectComposedDataSetReal((PetscObject)y,NormIds[i],norms[i]);
1612:     }
1613:   }
1614: #endif

1616:   PetscLogEventEnd(VEC_Copy,x,y,0,0);
1617:   return(0);
1618: }

1620: /*@
1621:    VecSwap - Swaps the vectors x and y.

1623:    Logically Collective on Vec

1625:    Input Parameters:
1626: .  x, y  - the vectors

1628:    Level: advanced

1630:    Concepts: vector^swapping values

1632: @*/
1633: PetscErrorCode  VecSwap(Vec x,Vec y)
1634: {
1635:   PetscReal      normxs[4]={0.0,0.0,0.0,0.0},normys[4]={0.0,0.0,0.0,0.0};
1636:   PetscBool      flgxs[4],flgys[4];
1638:   PetscInt       i;

1646:   VecCheckSameSize(x,1,y,2);
1647:   if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1648:   if (y->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");

1650:   PetscLogEventBegin(VEC_Swap,x,y,0,0);
1651:   for (i=0; i<4; i++) {
1652:     PetscObjectComposedDataGetReal((PetscObject)x,NormIds[i],normxs[i],flgxs[i]);
1653:     PetscObjectComposedDataGetReal((PetscObject)y,NormIds[i],normys[i],flgys[i]);
1654:   }
1655:   (*x->ops->swap)(x,y);
1656:   PetscObjectStateIncrease((PetscObject)x);
1657:   PetscObjectStateIncrease((PetscObject)y);
1658:   for (i=0; i<4; i++) {
1659:     if (flgxs[i]) {
1660:       PetscObjectComposedDataSetReal((PetscObject)y,NormIds[i],normxs[i]);
1661:     }
1662:     if (flgys[i]) {
1663:       PetscObjectComposedDataSetReal((PetscObject)x,NormIds[i],normys[i]);
1664:     }
1665:   }
1666:   PetscLogEventEnd(VEC_Swap,x,y,0,0);
1667:   return(0);
1668: }

1670: /*
1671:   VecStashViewFromOptions - Processes command line options to determine if/how an VecStash object is to be viewed.

1673:   Collective on VecStash

1675:   Input Parameters:
1676: + obj   - the VecStash object
1677: . bobj - optional other object that provides the prefix
1678: - optionname - option to activate viewing

1680:   Level: intermediate

1682:   Developer Note: This cannot use PetscObjectViewFromOptions() because it takes a Vec as an argument but does not use VecView

1684: */
1685: PetscErrorCode VecStashViewFromOptions(Vec obj,PetscObject bobj,const char optionname[])
1686: {
1687:   PetscErrorCode    ierr;
1688:   PetscViewer       viewer;
1689:   PetscBool         flg;
1690:   PetscViewerFormat format;
1691:   char              *prefix;

1694:   prefix = bobj ? bobj->prefix : ((PetscObject)obj)->prefix;
1695:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)obj),prefix,optionname,&viewer,&format,&flg);
1696:   if (flg) {
1697:     PetscViewerPushFormat(viewer,format);
1698:     VecStashView(obj,viewer);
1699:     PetscViewerPopFormat(viewer);
1700:     PetscViewerDestroy(&viewer);
1701:   }
1702:   return(0);
1703: }

1705: /*@
1706:    VecStashView - Prints the entries in the vector stash and block stash.

1708:    Collective on Vec

1710:    Input Parameters:
1711: +  v - the vector
1712: -  viewer - the viewer

1714:    Level: advanced

1716:    Concepts: vector^stash
1717:    Concepts: stash^vector

1719: .seealso: VecSetBlockSize(), VecSetValues(), VecSetValuesBlocked()

1721: @*/
1722: PetscErrorCode  VecStashView(Vec v,PetscViewer viewer)
1723: {
1725:   PetscMPIInt    rank;
1726:   PetscInt       i,j;
1727:   PetscBool      match;
1728:   VecStash       *s;
1729:   PetscScalar    val;


1736:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&match);
1737:   if (!match) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_SUP,"Stash viewer only works with ASCII viewer not %s\n",((PetscObject)v)->type_name);
1738:   PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
1739:   MPI_Comm_rank(PetscObjectComm((PetscObject)v),&rank);
1740:   s    = &v->bstash;

1742:   /* print block stash */
1743:   PetscViewerASCIIPushSynchronized(viewer);
1744:   PetscViewerASCIISynchronizedPrintf(viewer,"[%d]Vector Block stash size %D block size %D\n",rank,s->n,s->bs);
1745:   for (i=0; i<s->n; i++) {
1746:     PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Element %D ",rank,s->idx[i]);
1747:     for (j=0; j<s->bs; j++) {
1748:       val = s->array[i*s->bs+j];
1749: #if defined(PETSC_USE_COMPLEX)
1750:       PetscViewerASCIISynchronizedPrintf(viewer,"(%18.16e %18.16e) ",PetscRealPart(val),PetscImaginaryPart(val));
1751: #else
1752:       PetscViewerASCIISynchronizedPrintf(viewer,"%18.16e ",val);
1753: #endif
1754:     }
1755:     PetscViewerASCIISynchronizedPrintf(viewer,"\n");
1756:   }
1757:   PetscViewerFlush(viewer);

1759:   s = &v->stash;

1761:   /* print basic stash */
1762:   PetscViewerASCIISynchronizedPrintf(viewer,"[%d]Vector stash size %D\n",rank,s->n);
1763:   for (i=0; i<s->n; i++) {
1764:     val = s->array[i];
1765: #if defined(PETSC_USE_COMPLEX)
1766:     PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Element %D (%18.16e %18.16e) ",rank,s->idx[i],PetscRealPart(val),PetscImaginaryPart(val));
1767: #else
1768:     PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Element %D %18.16e\n",rank,s->idx[i],val);
1769: #endif
1770:   }
1771:   PetscViewerFlush(viewer);
1772:   PetscViewerASCIIPopSynchronized(viewer);
1773:   PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
1774:   return(0);
1775: }

1777: PetscErrorCode PetscOptionsGetVec(PetscOptions options,const char prefix[],const char key[],Vec v,PetscBool *set)
1778: {
1779:   PetscInt       i,N,rstart,rend;
1781:   PetscScalar    *xx;
1782:   PetscReal      *xreal;
1783:   PetscBool      iset;

1786:   VecGetOwnershipRange(v,&rstart,&rend);
1787:   VecGetSize(v,&N);
1788:   PetscCalloc1(N,&xreal);
1789:   PetscOptionsGetRealArray(options,prefix,key,xreal,&N,&iset);
1790:   if (iset) {
1791:     VecGetArray(v,&xx);
1792:     for (i=rstart; i<rend; i++) xx[i-rstart] = xreal[i];
1793:     VecRestoreArray(v,&xx);
1794:   }
1795:   PetscFree(xreal);
1796:   if (set) *set = iset;
1797:   return(0);
1798: }

1800: /*@
1801:    VecGetLayout - get PetscLayout describing vector layout

1803:    Not Collective

1805:    Input Arguments:
1806: .  x - the vector

1808:    Output Arguments:
1809: .  map - the layout

1811:    Level: developer

1813: .seealso: VecGetSizes(), VecGetOwnershipRange(), VecGetOwnershipRanges()
1814: @*/
1815: PetscErrorCode VecGetLayout(Vec x,PetscLayout *map)
1816: {

1820:   *map = x->map;
1821:   return(0);
1822: }

1824: /*@
1825:    VecSetLayout - set PetscLayout describing vector layout

1827:    Not Collective

1829:    Input Arguments:
1830: +  x - the vector
1831: -  map - the layout

1833:    Notes:
1834:    It is normally only valid to replace the layout with a layout known to be equivalent.

1836:    Level: developer

1838: .seealso: VecGetLayout(), VecGetSizes(), VecGetOwnershipRange(), VecGetOwnershipRanges()
1839: @*/
1840: PetscErrorCode VecSetLayout(Vec x,PetscLayout map)
1841: {

1846:   PetscLayoutReference(map,&x->map);
1847:   return(0);
1848: }

1850: PetscErrorCode VecSetInf(Vec xin)
1851: {
1852:   PetscInt       i,n = xin->map->n;
1853:   PetscScalar    *xx;
1854:   PetscScalar    zero=0.0,one=1.0,inf=one/zero;

1858:   VecGetArray(xin,&xx);
1859:   for (i=0; i<n; i++) xx[i] = inf;
1860:   VecRestoreArray(xin,&xx);
1861:   return(0);
1862: }