Actual source code: vector.c

petsc-3.8.4 2018-03-24
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, VEC_CUSPCopyFromGPU, VEC_CUSPCopyToGPU;
 16: PetscLogEvent VEC_CUSPCopyFromGPUSome, VEC_CUSPCopyToGPUSome;
 17: PetscLogEvent VEC_ViennaCLCopyFromGPU, VEC_ViennaCLCopyToGPU;
 18: PetscLogEvent VEC_CUDACopyFromGPU, VEC_CUDACopyToGPU;
 19: PetscLogEvent VEC_CUDACopyFromGPUSome, VEC_CUDACopyToGPUSome;

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

 25:     Not collective

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

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

 36:    Level: advanced

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

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

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

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

 56:    Logically Collective on Vec

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

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

 65:    Level: intermediate

 67:    Concepts: vector^setting values with local numbering

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


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

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

 91:    Not Collective

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

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

 99:    Level: advanced

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

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

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

120:    Collective on Vec

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

125:    Level: beginner

127:    Concepts: assembly^vectors

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

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

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

152:    Collective on Vec

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

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

166:    Level: beginner

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

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

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

189:    Logically Collective on Vec

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

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

197:    Level: advanced

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

202:    Concepts: vector^pointwise multiply

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

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


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

230:    Logically Collective on Vec

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

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

238:    Level: advanced

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

243:    Concepts: vector^pointwise multiply

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

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

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

270:    Logically Collective on Vec

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

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

278:    Level: advanced

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

282:    Concepts: vector^pointwise multiply

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

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

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

309:    Logically Collective on Vec

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

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

317:    Level: advanced

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

321:    Concepts: vector^pointwise divide

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

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


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

349:    Collective on Vec

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

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

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

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

364:    Level: beginner

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

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

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

384:    Collective on Vec

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

389:    Level: beginner

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

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

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

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

416:    Collective on Vec

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

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

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

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

434:    Level: intermediate

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

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

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

453:    Collective on Vec

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

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

463:    Level: intermediate

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

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

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

485:    Collective on Vec

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

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

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

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

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

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

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

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

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

535:    Level: beginner

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

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

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

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

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

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

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


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

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

616:    Not Collective

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

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

624:    Level: beginner

626:    Concepts: vector^local size

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

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

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

647:    Not Collective

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

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

655:    Level: beginner

657:    Concepts: vector^size

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

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

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

680:    Not Collective

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

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

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

692:    Fortran: NULL_INTEGER should be used instead of NULL

694:    Level: beginner

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

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

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

720:    Not Collective

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

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

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

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

733:    Level: beginner

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

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

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

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

754:    Collective on Vec

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

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

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

778:    Level: intermediate

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

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

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

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

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

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

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

826:    Not Collective

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

831:    Level: developer

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

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

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

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

854:   Collective on PetscViewer

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

862:    Level: intermediate

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

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

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

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

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

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

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

902:   Concepts: vector^loading from file

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

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

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


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

937:    Logically Collective on Vec

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

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

945:    Level: intermediate

947:    Concepts: vector^reciprocal

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

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

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

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

969:    Logically Collective on Vec

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

976:    Level: advanced

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

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

989:     This function is not currently available from Fortran.

991: .keywords: vector, set, operation

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


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

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

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

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

1025:    Level: intermediate

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

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

1037:    Concepts: vector^stash
1038:    Concepts: stash^vector

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

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

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

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

1057:    Logically Collective on Vec

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

1062:    Level: intermediate

1064:    Concepts: vector^conjugate

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

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

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

1088:    Logically Collective on Vec

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

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

1096:    Level: advanced

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

1100:    Concepts: vector^pointwise multiply

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

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

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

1129:    Logically Collective on Vec

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

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

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

1146:    Level: intermediate

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

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

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

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

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

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

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

1184:   Logically Collective on Vec

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

1189:   Level: beginner

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

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

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

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

1212:   Collective on Vec

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

1217:   Level: intermediate

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

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

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

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

1251:   Collective on Vec

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

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

1259:   Level: beginner

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

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


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

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

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

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

1292:   Collective on Vec

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

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

1303:   Level: intermediate

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

1314:   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);
1315:   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);
1316:   v->map->n = n;
1317:   v->map->N = N;
1318:   if (v->ops->create) {
1319:     (*v->ops->create)(v);
1320:     v->ops->create = 0;
1321:   }
1322:   return(0);
1323: }

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

1329:    Logically Collective on Vec

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

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

1338:    Level: advanced

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

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

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

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

1361:    Not Collective

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

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

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

1372:    Level: advanced

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

1376:    Concepts: vector^block size
1377:    Concepts: block^vector

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

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

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

1395:    Logically Collective on Vec

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

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

1405:    Level: advanced

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

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

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

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

1425:    Logically Collective on Vec

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

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

1435:    Level: advanced

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

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

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

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

1455:    Not Collective

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

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

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

1466:    Level: advanced

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

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

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

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

1485:    Collective on Vec

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

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

1494:    Level: advanced

1496: .keywords: Vec, setup

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

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

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

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

1527:    Logically Collective on Vec

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

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

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

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

1544:    Level: beginner

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

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

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

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

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

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

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

1624:    Logically Collective on Vec

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

1629:    Level: advanced

1631:    Concepts: vector^swapping values

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

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

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

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

1674:   Collective on VecStash

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

1681:   Level: intermediate

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

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

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

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

1709:    Collective on Vec

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

1715:    Level: advanced

1717:    Concepts: vector^stash
1718:    Concepts: stash^vector

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

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


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

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

1760:   s = &v->stash;

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

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

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

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

1804:    Not Collective

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

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

1812:    Level: developer

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

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

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

1828:    Not Collective

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

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

1837:    Level: developer

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

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

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

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