Actual source code: vector.c

petsc-dev 2014-04-18
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>    /*I  "petscvec.h"   I*/

  8: /* Logging support */
  9: PetscClassId  VEC_CLASSID;
 10: PetscLogEvent VEC_View, VEC_Max, VEC_Min, VEC_DotBarrier, VEC_Dot, VEC_MDotBarrier, VEC_MDot, VEC_TDot;
 11: PetscLogEvent VEC_Norm, VEC_Normalize, VEC_Scale, VEC_Copy, VEC_Set, VEC_AXPY, VEC_AYPX, VEC_WAXPY;
 12: PetscLogEvent VEC_MTDot, VEC_NormBarrier, VEC_MAXPY, VEC_Swap, VEC_AssemblyBegin, VEC_ScatterBegin, VEC_ScatterEnd;
 13: PetscLogEvent VEC_AssemblyEnd, VEC_PointwiseMult, VEC_SetValues, VEC_Load, VEC_ScatterBarrier;
 14: PetscLogEvent VEC_SetRandom, VEC_ReduceArithmetic, VEC_ReduceBarrier, VEC_ReduceCommunication,VEC_ReduceBegin,VEC_ReduceEnd,VEC_Ops;
 15: PetscLogEvent VEC_DotNormBarrier, VEC_DotNorm, VEC_AXPBYPCZ, VEC_CUSPCopyFromGPU, VEC_CUSPCopyToGPU;
 16: PetscLogEvent VEC_CUSPCopyFromGPUSome, VEC_CUSPCopyToGPUSome;
 17: PetscLogEvent VEC_ViennaCLCopyFromGPU, VEC_ViennaCLCopyToGPU;

 19: extern PetscErrorCode VecStashGetInfo_Private(VecStash*,PetscInt*,PetscInt*);
 22: /*@
 23:    VecStashGetInfo - Gets how many values are currently in the vector stash, i.e. need
 24:        to be communicated to other processors during the VecAssemblyBegin/End() process

 26:     Not collective

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

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

 37:    Level: advanced

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

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

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

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

 59:    Logically Collective on Vec

 61:    Input Parameters:
 62: +  x - vector
 63: -  mapping - mapping created with ISLocalToGlobalMappingCreate() or ISLocalToGlobalMappingCreateIS()

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

 68:    Level: intermediate

 70:    Concepts: vector^setting values with local numbering

 72: seealso:  VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetValuesLocal(),
 73:            VecSetLocalToGlobalMappingBlock(), VecSetValuesBlockedLocal()
 74: @*/
 75: PetscErrorCode  VecSetLocalToGlobalMapping(Vec x,ISLocalToGlobalMapping mapping)
 76: {


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

 93: /*@
 94:    VecSetLocalToGlobalMappingBlock - Sets a local numbering to global numbering used
 95:    by the routine VecSetValuesBlockedLocal() to allow users to insert vector entries
 96:    using a local (per-processor) numbering.

 98:    Logically Collective on Vec

100:    Input Parameters:
101: +  x - vector
102: -  mapping - mapping created with ISLocalToGlobalMappingCreate() or ISLocalToGlobalMappingCreateIS()

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

107:    Level: intermediate

109:    Concepts: vector^setting values blocked with local numbering

111: .seealso:  VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetValuesLocal(),
112:            VecSetLocalToGlobalMapping(), VecSetValuesBlockedLocal()
113: @*/
114: PetscErrorCode  VecSetLocalToGlobalMappingBlock(Vec x,ISLocalToGlobalMapping mapping)
115: {


122:   PetscLayoutSetISLocalToGlobalMappingBlock(x->map,mapping);
123:   return(0);
124: }

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

131:    Not Collective

133:    Input Parameter:
134: .  X - the vector

136:    Output Parameter:
137: .  mapping - the mapping

139:    Level: advanced

141:    Concepts: vectors^local to global mapping
142:    Concepts: local to global mapping^for vectors

144: .seealso:  VecSetValuesLocal(), VecGetLocalToGlobalMappingBlock()
145: @*/
146: PetscErrorCode VecGetLocalToGlobalMapping(Vec X,ISLocalToGlobalMapping *mapping)
147: {
152:   *mapping = X->map->mapping;
153:   return(0);
154: }

158: /*@
159:    VecGetLocalToGlobalMappingBlock - Gets the local-to-global numbering set by VecSetLocalToGlobalMappingBlock()

161:    Not Collective

163:    Input Parameters:
164: .  X - the vector

166:    Output Parameters:
167: .  mapping - the mapping

169:    Level: advanced

171:    Concepts: vectors^local to global mapping blocked
172:    Concepts: local to global mapping^for vectors, blocked

174: .seealso:  VecSetValuesBlockedLocal(), VecGetLocalToGlobalMapping()
175: @*/
176: PetscErrorCode VecGetLocalToGlobalMappingBlock(Vec X,ISLocalToGlobalMapping *mapping)
177: {
182:   *mapping = X->map->bmapping;
183:   return(0);
184: }

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

192:    Collective on Vec

194:    Input Parameter:
195: .  vec - the vector

197:    Level: beginner

199:    Concepts: assembly^vectors

201: .seealso: VecAssemblyEnd(), VecSetValues()
202: @*/
203: PetscErrorCode  VecAssemblyBegin(Vec vec)
204: {

210:   VecStashViewFromOptions(vec,NULL,"-vec_view_stash");
211:   PetscLogEventBegin(VEC_AssemblyBegin,vec,0,0,0);
212:   if (vec->ops->assemblybegin) {
213:     (*vec->ops->assemblybegin)(vec);
214:   }
215:   PetscLogEventEnd(VEC_AssemblyBegin,vec,0,0,0);
216:   PetscObjectStateIncrease((PetscObject)vec);
217:   return(0);
218: }

222: /*@
223:    VecAssemblyEnd - Completes assembling the vector.  This routine should
224:    be called after VecAssemblyBegin().

226:    Collective on Vec

228:    Input Parameter:
229: .  vec - the vector

231:    Options Database Keys:
232: +  -vec_view - Prints vector in ASCII format
233: .  -vec_view ::ascii_matlab - Prints vector in ASCII MATLAB format to stdout
234: .  -vec_view matlab:filename - Prints vector in MATLAB format to matlaboutput.mat
235: .  -vec_view draw - Activates vector viewing using drawing tools
236: .  -display <name> - Sets display name (default is host)
237: .  -draw_pause <sec> - Sets number of seconds to pause after display
238: -  -vec_view socket - Activates vector viewing using a socket

240:    Level: beginner

242: .seealso: VecAssemblyBegin(), VecSetValues()
243: @*/
244: PetscErrorCode  VecAssemblyEnd(Vec vec)
245: {

250:   PetscLogEventBegin(VEC_AssemblyEnd,vec,0,0,0);
252:   if (vec->ops->assemblyend) {
253:     (*vec->ops->assemblyend)(vec);
254:   }
255:   PetscLogEventEnd(VEC_AssemblyEnd,vec,0,0,0);
256:   VecViewFromOptions(vec,NULL,"-vec_view");
257:   return(0);
258: }

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

265:    Logically Collective on Vec

267:    Input Parameters:
268: .  x, y  - the vectors

270:    Output Parameter:
271: .  w - the result

273:    Level: advanced

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

278:    Concepts: vector^pointwise multiply

280: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
281: @*/
282: PetscErrorCode  VecPointwiseMax(Vec w,Vec x,Vec y)
283: {

295:   if (x->map->N != y->map->N || x->map->N != w->map->N) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths");
296:   if (x->map->n != y->map->n || x->map->n != w->map->n) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");

298:   (*w->ops->pointwisemax)(w,x,y);
299:   PetscObjectStateIncrease((PetscObject)w);
300:   return(0);
301: }


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

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.
320:           For complex numbers compares only the real part

322:    Concepts: vector^pointwise multiply

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

339:   if (x->map->N != y->map->N || x->map->N != w->map->N) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths");
340:   if (x->map->n != y->map->n || x->map->n != w->map->n) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");

342:   (*w->ops->pointwisemin)(w,x,y);
343:   PetscObjectStateIncrease((PetscObject)w);
344:   return(0);
345: }

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

352:    Logically Collective on Vec

354:    Input Parameters:
355: .  x, y  - the vectors

357:    Output Parameter:
358: .  w - the result

360:    Level: advanced

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

364:    Concepts: vector^pointwise multiply

366: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMin(), VecPointwiseMax(), VecMaxPointwiseDivide()
367: @*/
368: PetscErrorCode  VecPointwiseMaxAbs(Vec w,Vec x,Vec y)
369: {

381:   if (x->map->N != y->map->N || x->map->N != w->map->N) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths");
382:   if (x->map->n != y->map->n || x->map->n != w->map->n) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");

384:   (*w->ops->pointwisemaxabs)(w,x,y);
385:   PetscObjectStateIncrease((PetscObject)w);
386:   return(0);
387: }

391: /*@
392:    VecPointwiseDivide - Computes the componentwise division w = x/y.

394:    Logically Collective on Vec

396:    Input Parameters:
397: .  x, y  - the vectors

399:    Output Parameter:
400: .  w - the result

402:    Level: advanced

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

406:    Concepts: vector^pointwise divide

408: .seealso: VecPointwiseMult(), VecPointwiseMax(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
409: @*/
410: PetscErrorCode  VecPointwiseDivide(Vec w,Vec x,Vec y)
411: {

423:   if (x->map->N != y->map->N || x->map->N != w->map->N) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths");
424:   if (x->map->n != y->map->n || x->map->n != w->map->n) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");

426:   (*w->ops->pointwisedivide)(w,x,y);
427:   PetscObjectStateIncrease((PetscObject)w);
428:   return(0);
429: }


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

437:    Collective on Vec

439:    Input Parameters:
440: .  v - a vector to mimic

442:    Output Parameter:
443: .  newv - location to put new vector

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

449:    Use VecDestroy() to free the space. Use VecDuplicateVecs() to get several
450:    vectors.

452:    Level: beginner

454: .seealso: VecDestroy(), VecDuplicateVecs(), VecCreate(), VecCopy()
455: @*/
456: PetscErrorCode  VecDuplicate(Vec v,Vec *newv)
457: {

464:   (*v->ops->duplicate)(v,newv);
465:   PetscObjectStateIncrease((PetscObject)*newv);
466:   return(0);
467: }

471: /*@
472:    VecDestroy - Destroys a vector.

474:    Collective on Vec

476:    Input Parameters:
477: .  v  - the vector

479:    Level: beginner

481: .seealso: VecDuplicate(), VecDestroyVecs()
482: @*/
483: PetscErrorCode  VecDestroy(Vec *v)
484: {

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

492:   PetscObjectSAWsViewOff((PetscObject)*v);
493:   /* destroy the internal part */
494:   if ((*v)->ops->destroy) {
495:     (*(*v)->ops->destroy)(*v);
496:   }
497:   /* destroy the external/common part */
498:   PetscLayoutDestroy(&(*v)->map);
499:   PetscHeaderDestroy(v);
500:   return(0);
501: }

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

508:    Collective on Vec

510:    Input Parameters:
511: +  m - the number of vectors to obtain
512: -  v - a vector to mimic

514:    Output Parameter:
515: .  V - location to put pointer to array of vectors

517:    Notes:
518:    Use VecDestroyVecs() to free the space. Use VecDuplicate() to form a single
519:    vector.

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

526:    Level: intermediate

528: .seealso:  VecDestroyVecs(), VecDuplicate(), VecCreate(), VecDuplicateVecsF90()
529: @*/
530: PetscErrorCode  VecDuplicateVecs(Vec v,PetscInt m,Vec *V[])
531: {

538:   (*v->ops->duplicatevecs)(v, m,V);
539:   return(0);
540: }

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

547:    Collective on Vec

549:    Input Parameters:
550: +  vv - pointer to pointer to array of vector pointers
551: -  m - the number of vectors previously obtained

553:    Fortran Note:
554:    The Fortran interface is slightly different from that given below.
555:    See the Fortran chapter of the users manual

557:    Level: intermediate

559: .seealso: VecDuplicateVecs(), VecDestroyVecsf90()
560: @*/
561: PetscErrorCode  VecDestroyVecs(PetscInt m,Vec *vv[])
562: {

567:   if (!*vv) return(0);
570:   if (m < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Trying to destroy negative number of vectors %D",m);
571:   (*(**vv)->ops->destroyvecs)(m,*vv);
572:   *vv  = 0;
573:   return(0);
574: }

578: /*@C
579:    VecView - Views a vector object.

581:    Collective on Vec

583:    Input Parameters:
584: +  vec - the vector
585: -  viewer - an optional visualization context

587:    Notes:
588:    The available visualization contexts include
589: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
590: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
591:          output where only the first processor opens
592:          the file.  All other processors send their
593:          data to the first processor to print.

595:    You can change the format the vector is printed using the
596:    option PetscViewerSetFormat().

598:    The user can open alternative visualization contexts with
599: +    PetscViewerASCIIOpen() - Outputs vector to a specified file
600: .    PetscViewerBinaryOpen() - Outputs vector in binary to a
601:          specified file; corresponding input uses VecLoad()
602: .    PetscViewerDrawOpen() - Outputs vector to an X window display
603: -    PetscViewerSocketOpen() - Outputs vector to Socket viewer

605:    The user can call PetscViewerSetFormat() to specify the output
606:    format of ASCII printed objects (when using PETSC_VIEWER_STDOUT_SELF,
607:    PETSC_VIEWER_STDOUT_WORLD and PetscViewerASCIIOpen).  Available formats include
608: +    PETSC_VIEWER_DEFAULT - default, prints vector contents
609: .    PETSC_VIEWER_ASCII_MATLAB - prints vector contents in MATLAB format
610: .    PETSC_VIEWER_ASCII_INDEX - prints vector contents, including indices of vector elements
611: -    PETSC_VIEWER_ASCII_COMMON - prints vector contents, using a
612:          format common among all vector types

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

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

622:    Level: beginner

624:    Concepts: vector^printing
625:    Concepts: vector^saving to disk

627: .seealso: PetscViewerASCIIOpen(), PetscViewerDrawOpen(), PetscDrawLGCreate(),
628:           PetscViewerSocketOpen(), PetscViewerBinaryOpen(), VecLoad(), PetscViewerCreate(),
629:           PetscRealView(), PetscScalarView(), PetscIntView()
630: @*/
631: PetscErrorCode  VecView(Vec vec,PetscViewer viewer)
632: {
633:   PetscErrorCode    ierr;
634:   PetscBool         iascii;
635:   PetscViewerFormat format;

640:   if (!viewer) {
641:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)vec),&viewer);
642:   }
645:   if (vec->stash.n || vec->bstash.n) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call VecAssemblyBegin/End() before viewing this vector");

647:   PetscLogEventBegin(VEC_View,vec,viewer,0,0);
648:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
649:   if (iascii) {
650:     PetscInt rows,bs;

652:     PetscObjectPrintClassNamePrefixType((PetscObject)vec,viewer);
653:     PetscViewerGetFormat(viewer,&format);
654:     if (format == PETSC_VIEWER_ASCII_INFO || format == PETSC_VIEWER_ASCII_INFO_DETAIL) {
655:       PetscViewerASCIIPushTab(viewer);
656:       VecGetSize(vec,&rows);
657:       VecGetBlockSize(vec,&bs);
658:       if (bs != 1) {
659:         PetscViewerASCIIPrintf(viewer,"length=%D, bs=%D\n",rows,bs);
660:       } else {
661:         PetscViewerASCIIPrintf(viewer,"length=%D\n",rows);
662:       }
663:       PetscViewerASCIIPopTab(viewer);
664:     }
665:   }
666:   (*vec->ops->view)(vec,viewer);
667:   PetscLogEventEnd(VEC_View,vec,viewer,0,0);
668:   return(0);
669: }

671: #if defined(PETSC_USE_DEBUG)
672: #include <../src/sys/totalview/tv_data_display.h>
673: PETSC_UNUSED static int TV_display_type(const struct _p_Vec *v)
674: {
675:   const PetscScalar *values;
676:   char              type[32];
677:   PetscErrorCode    ierr;


680:   TV_add_row("Local rows", "int", &v->map->n);
681:   TV_add_row("Global rows", "int", &v->map->N);
682:   TV_add_row("Typename", TV_ascii_string_type, ((PetscObject)v)->type_name);
683:   VecGetArrayRead((Vec)v,&values);
684:   PetscSNPrintf(type,32,"double[%d]",v->map->n);
685:   TV_add_row("values",type, values);
686:   VecRestoreArrayRead((Vec)v,&values);
687:   return TV_format_OK;
688: }
689: #endif

693: /*@
694:    VecGetSize - Returns the global number of elements of the vector.

696:    Not Collective

698:    Input Parameter:
699: .  x - the vector

701:    Output Parameters:
702: .  size - the global length of the vector

704:    Level: beginner

706:    Concepts: vector^local size

708: .seealso: VecGetLocalSize()
709: @*/
710: PetscErrorCode  VecGetSize(Vec x,PetscInt *size)
711: {

718:   (*x->ops->getsize)(x,size);
719:   return(0);
720: }

724: /*@
725:    VecGetLocalSize - Returns the number of elements of the vector stored
726:    in local memory. This routine may be implementation dependent, so use
727:    with care.

729:    Not Collective

731:    Input Parameter:
732: .  x - the vector

734:    Output Parameter:
735: .  size - the length of the local piece of the vector

737:    Level: beginner

739:    Concepts: vector^size

741: .seealso: VecGetSize()
742: @*/
743: PetscErrorCode  VecGetLocalSize(Vec x,PetscInt *size)
744: {

751:   (*x->ops->getlocalsize)(x,size);
752:   return(0);
753: }

757: /*@C
758:    VecGetOwnershipRange - Returns the range of indices owned by
759:    this processor, assuming that the vectors are laid out with the
760:    first n1 elements on the first processor, next n2 elements on the
761:    second, etc.  For certain parallel layouts this range may not be
762:    well defined.

764:    Not Collective

766:    Input Parameter:
767: .  x - the vector

769:    Output Parameters:
770: +  low - the first local element, pass in NULL if not interested
771: -  high - one more than the last local element, pass in NULL if not interested

773:    Note:
774:    The high argument is one more than the last element stored locally.

776:    Fortran: NULL_INTEGER should be used instead of NULL

778:    Level: beginner

780:    Concepts: ownership^of vectors
781:    Concepts: vector^ownership of elements

783: .seealso:   MatGetOwnershipRange(), MatGetOwnershipRanges(), VecGetOwnershipRanges()
784: @*/
785: PetscErrorCode  VecGetOwnershipRange(Vec x,PetscInt *low,PetscInt *high)
786: {
792:   if (low)  *low  = x->map->rstart;
793:   if (high) *high = x->map->rend;
794:   return(0);
795: }

799: /*@C
800:    VecGetOwnershipRanges - Returns the range of indices owned by EACH processor,
801:    assuming that the vectors are laid out with the
802:    first n1 elements on the first processor, next n2 elements on the
803:    second, etc.  For certain parallel layouts this range may not be
804:    well defined.

806:    Not Collective

808:    Input Parameter:
809: .  x - the vector

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

814:    Note:
815:    The high argument is one more than the last element stored locally.

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

819:    Level: beginner

821:    Concepts: ownership^of vectors
822:    Concepts: vector^ownership of elements

824: .seealso:   MatGetOwnershipRange(), MatGetOwnershipRanges(), VecGetOwnershipRange()
825: @*/
826: PetscErrorCode  VecGetOwnershipRanges(Vec x,const PetscInt *ranges[])
827: {

833:   PetscLayoutGetRanges(x->map,ranges);
834:   return(0);
835: }

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

842:    Collective on Vec

844:    Input Parameter:
845: +  x - the vector
846: .  op - the option
847: -  flag - turn the option on or off

849:    Supported Options:
850: +     VEC_IGNORE_OFF_PROC_ENTRIES, which causes VecSetValues() to ignore
851:           entries destined to be stored on a separate processor. This can be used
852:           to eliminate the global reduction in the VecAssemblyXXXX() if you know
853:           that you have only used VecSetValues() to set local elements
854: .     VEC_IGNORE_NEGATIVE_INDICES, which means you can pass negative indices
855:           in ix in calls to VecSetValues() or VecGetValues(). These rows are simply
856:           ignored.

858:    Level: intermediate

860: @*/
861: PetscErrorCode  VecSetOption(Vec x,VecOption op,PetscBool flag)
862: {

868:   if (x->ops->setoption) {
869:     (*x->ops->setoption)(x,op,flag);
870:   }
871:   return(0);
872: }

876: /* Default routines for obtaining and releasing; */
877: /* may be used by any implementation */
878: PetscErrorCode VecDuplicateVecs_Default(Vec w,PetscInt m,Vec *V[])
879: {
881:   PetscInt       i;

886:   if (m <= 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"m must be > 0: m = %D",m);
887:   PetscMalloc1(m,V);
888:   for (i=0; i<m; i++) {VecDuplicate(w,*V+i);}
889:   return(0);
890: }

894: PetscErrorCode VecDestroyVecs_Default(PetscInt m,Vec v[])
895: {
897:   PetscInt       i;

901:   if (m <= 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"m must be > 0: m = %D",m);
902:   for (i=0; i<m; i++) {VecDestroy(&v[i]);}
903:   PetscFree(v);
904:   return(0);
905: }

909: /*@
910:    VecResetArray - Resets a vector to use its default memory. Call this
911:    after the use of VecPlaceArray().

913:    Not Collective

915:    Input Parameters:
916: .  vec - the vector

918:    Level: developer

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

922: @*/
923: PetscErrorCode  VecResetArray(Vec vec)
924: {

930:   if (vec->ops->resetarray) {
931:     (*vec->ops->resetarray)(vec);
932:   } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Cannot reset array in this type of vector");
933:   PetscObjectStateIncrease((PetscObject)vec);
934:   return(0);
935: }

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

943:   Collective on PetscViewer

945:   Input Parameters:
946: + newvec - the newly loaded vector, this needs to have been created with VecCreate() or
947:            some related function before a call to VecLoad().
948: - viewer - binary file viewer, obtained from PetscViewerBinaryOpen() or
949:            HDF5 file viewer, obtained from PetscViewerHDF5Open()

951:    Level: intermediate

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

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

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

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

968:   Notes for advanced users:
969:   Most users should not need to know the details of the binary storage
970:   format, since VecLoad() and VecView() completely hide these details.
971:   But for anyone who's interested, the standard binary matrix storage
972:   format is
973: .vb
974:      int    VEC_FILE_CLASSID
975:      int    number of rows
976:      PetscScalar *values of all entries
977: .ve

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

985:   Concepts: vector^loading from file

987: .seealso: PetscViewerBinaryOpen(), VecView(), MatLoad(), VecLoad()
988: @*/
989: PetscErrorCode  VecLoad(Vec newvec, PetscViewer viewer)
990: {
992:   PetscBool      isbinary,ishdf5;

997:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
998:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5);
999:   if (!isbinary && !ishdf5) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Invalid viewer; open viewer with PetscViewerBinaryOpen()");

1001:   PetscLogEventBegin(VEC_Load,viewer,0,0,0);
1002:   if (!((PetscObject)newvec)->type_name && !newvec->ops->create) {
1003:     VecSetType(newvec, VECSTANDARD);
1004:   }
1005:   (*newvec->ops->load)(newvec,viewer);
1006:   PetscLogEventEnd(VEC_Load,viewer,0,0,0);
1007:   return(0);
1008: }


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

1016:    Logically Collective on Vec

1018:    Input Parameter:
1019: .  vec - the vector

1021:    Output Parameter:
1022: .  vec - the vector reciprocal

1024:    Level: intermediate

1026:    Concepts: vector^reciprocal

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

1030: @*/
1031: PetscErrorCode  VecReciprocal(Vec vec)
1032: {

1038:   if (vec->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1039:   if (!vec->ops->reciprocal) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Vector does not support reciprocal operation");
1040:   (*vec->ops->reciprocal)(vec);
1041:   PetscObjectStateIncrease((PetscObject)vec);
1042:   return(0);
1043: }

1047: PetscErrorCode  VecSetOperation(Vec vec,VecOperation op, void (*f)(void))
1048: {
1051:   (((void(**)(void))vec->ops)[(int)op]) = f;
1052:   return(0);
1053: }


1058: /*@
1059:    VecStashSetInitialSize - sets the sizes of the vec-stash, that is
1060:    used during the assembly process to store values that belong to
1061:    other processors.

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

1065:    Input Parameters:
1066: +  vec   - the vector
1067: .  size  - the initial size of the stash.
1068: -  bsize - the initial size of the block-stash(if used).

1070:    Options Database Keys:
1071: +   -vecstash_initial_size <size> or <size0,size1,...sizep-1>
1072: -   -vecstash_block_initial_size <bsize> or <bsize0,bsize1,...bsizep-1>

1074:    Level: intermediate

1076:    Notes:
1077:      The block-stash is used for values set with VecSetValuesBlocked() while
1078:      the stash is used for values set with VecSetValues()

1080:      Run with the option -info and look for output of the form
1081:      VecAssemblyBegin_MPIXXX:Stash has MM entries, uses nn mallocs.
1082:      to determine the appropriate value, MM, to use for size and
1083:      VecAssemblyBegin_MPIXXX:Block-Stash has BMM entries, uses nn mallocs.
1084:      to determine the value, BMM to use for bsize

1086:    Concepts: vector^stash
1087:    Concepts: stash^vector

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

1091: @*/
1092: PetscErrorCode  VecStashSetInitialSize(Vec vec,PetscInt size,PetscInt bsize)
1093: {

1098:   VecStashSetInitialSize_Private(&vec->stash,size);
1099:   VecStashSetInitialSize_Private(&vec->bstash,bsize);
1100:   return(0);
1101: }

1105: /*@
1106:    VecConjugate - Conjugates a vector.

1108:    Logically Collective on Vec

1110:    Input Parameters:
1111: .  x - the vector

1113:    Level: intermediate

1115:    Concepts: vector^conjugate

1117: @*/
1118: PetscErrorCode  VecConjugate(Vec x)
1119: {
1120: #if defined(PETSC_USE_COMPLEX)

1126:   if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1127:   (*x->ops->conjugate)(x);
1128:   /* we need to copy norms here */
1129:   PetscObjectStateIncrease((PetscObject)x);
1130:   return(0);
1131: #else
1132:   return(0);
1133: #endif
1134: }

1138: /*@
1139:    VecPointwiseMult - Computes the componentwise multiplication w = x*y.

1141:    Logically Collective on Vec

1143:    Input Parameters:
1144: .  x, y  - the vectors

1146:    Output Parameter:
1147: .  w - the result

1149:    Level: advanced

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

1153:    Concepts: vector^pointwise multiply

1155: .seealso: VecPointwiseDivide(), VecPointwiseMax(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
1156: @*/
1157: PetscErrorCode  VecPointwiseMult(Vec w, Vec x,Vec y)
1158: {

1170:   if (x->map->n != y->map->n || x->map->n != w->map->n) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");

1172:   PetscLogEventBegin(VEC_PointwiseMult,x,y,w,0);
1173:   (*w->ops->pointwisemult)(w,x,y);
1174:   PetscLogEventEnd(VEC_PointwiseMult,x,y,w,0);
1175:   PetscObjectStateIncrease((PetscObject)w);
1176:   return(0);
1177: }

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

1184:    Logically Collective on Vec

1186:    Input Parameters:
1187: +  x  - the vector
1188: -  rctx - the random number context, formed by PetscRandomCreate(), or NULL and
1189:           it will create one internally.

1191:    Output Parameter:
1192: .  x  - the vector

1194:    Example of Usage:
1195: .vb
1196:      PetscRandomCreate(PETSC_COMM_WORLD,&rctx);
1197:      VecSetRandom(x,rctx);
1198:      PetscRandomDestroy(rctx);
1199: .ve

1201:    Level: intermediate

1203:    Concepts: vector^setting to random
1204:    Concepts: random^vector

1206: .seealso: VecSet(), VecSetValues(), PetscRandomCreate(), PetscRandomDestroy()
1207: @*/
1208: PetscErrorCode  VecSetRandom(Vec x,PetscRandom rctx)
1209: {
1211:   PetscRandom    randObj = NULL;

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

1219:   if (!rctx) {
1220:     MPI_Comm comm;
1221:     PetscObjectGetComm((PetscObject)x,&comm);
1222:     PetscRandomCreate(comm,&randObj);
1223:     PetscRandomSetFromOptions(randObj);
1224:     rctx = randObj;
1225:   }

1227:   PetscLogEventBegin(VEC_SetRandom,x,rctx,0,0);
1228:   (*x->ops->setrandom)(x,rctx);
1229:   PetscLogEventEnd(VEC_SetRandom,x,rctx,0,0);

1231:   PetscRandomDestroy(&randObj);
1232:   PetscObjectStateIncrease((PetscObject)x);
1233:   return(0);
1234: }

1238: /*@
1239:   VecZeroEntries - puts a 0.0 in each element of a vector

1241:   Logically Collective on Vec

1243:   Input Parameter:
1244: . vec - The vector

1246:   Level: beginner

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

1253: .keywords: Vec, set, options, database
1254: .seealso: VecCreate(),  VecSetOptionsPrefix(), VecSet(), VecSetValues()
1255: @*/
1256: PetscErrorCode  VecZeroEntries(Vec vec)
1257: {

1261:   VecSet(vec,0);
1262:   return(0);
1263: }

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

1271:   Collective on Vec

1273:   Input Parameter:
1274: . vec - The vector

1276:   Level: intermediate

1278: .keywords: Vec, set, options, database, type
1279: .seealso: VecSetFromOptions(), VecSetType()
1280: */
1281: static PetscErrorCode VecSetTypeFromOptions_Private(Vec vec)
1282: {
1283:   PetscBool      opt;
1284:   VecType        defaultType;
1285:   char           typeName[256];
1286:   PetscMPIInt    size;

1290:   if (((PetscObject)vec)->type_name) defaultType = ((PetscObject)vec)->type_name;
1291:   else {
1292:     MPI_Comm_size(PetscObjectComm((PetscObject)vec), &size);
1293:     if (size > 1) defaultType = VECMPI;
1294:     else defaultType = VECSEQ;
1295:   }

1297:   if (!VecRegisterAllCalled) {VecRegisterAll();}
1298:   PetscOptionsFList("-vec_type","Vector type","VecSetType",VecList,defaultType,typeName,256,&opt);
1299:   if (opt) {
1300:     VecSetType(vec, typeName);
1301:   } else {
1302:     VecSetType(vec, defaultType);
1303:   }
1304:   return(0);
1305: }

1309: /*@
1310:   VecSetFromOptions - Configures the vector from the options database.

1312:   Collective on Vec

1314:   Input Parameter:
1315: . vec - The vector

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

1320:   Level: beginner

1322:   Concepts: vectors^setting options
1323:   Concepts: vectors^setting type

1325: .keywords: Vec, set, options, database
1326: .seealso: VecCreate(), VecSetOptionsPrefix()
1327: @*/
1328: PetscErrorCode  VecSetFromOptions(Vec vec)
1329: {


1335:   PetscObjectOptionsBegin((PetscObject)vec);
1336:   /* Handle vector type options */
1337:   VecSetTypeFromOptions_Private(vec);

1339:   /* Handle specific vector options */
1340:   if (vec->ops->setfromoptions) {
1341:     (*vec->ops->setfromoptions)(vec);
1342:   }

1344:   /* process any options handlers added with PetscObjectAddOptionsHandler() */
1345:   PetscObjectProcessOptionsHandlers((PetscObject)vec);
1346:   PetscOptionsEnd();
1347:   return(0);
1348: }

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

1355:   Collective on Vec

1357:   Input Parameters:
1358: + v - the vector
1359: . n - the local size (or PETSC_DECIDE to have it set)
1360: - N - the global size (or PETSC_DECIDE)

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

1366:   Level: intermediate

1368: .seealso: VecGetSize(), PetscSplitOwnership()
1369: @*/
1370: PetscErrorCode  VecSetSizes(Vec v, PetscInt n, PetscInt N)
1371: {

1377:   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);
1378:   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);
1379:   v->map->n = n;
1380:   v->map->N = N;
1381:   if (v->ops->create) {
1382:     (*v->ops->create)(v);
1383:     v->ops->create = 0;
1384:   }
1385:   return(0);
1386: }

1390: /*@
1391:    VecSetBlockSize - Sets the blocksize for future calls to VecSetValuesBlocked()
1392:    and VecSetValuesBlockedLocal().

1394:    Logically Collective on Vec

1396:    Input Parameter:
1397: +  v - the vector
1398: -  bs - the blocksize

1400:    Notes:
1401:    All vectors obtained by VecDuplicate() inherit the same blocksize.

1403:    Level: advanced

1405: .seealso: VecSetValuesBlocked(), VecSetLocalToGlobalMappingBlock(), VecGetBlockSize()

1407:   Concepts: block size^vectors
1408: @*/
1409: PetscErrorCode  VecSetBlockSize(Vec v,PetscInt bs)
1410: {

1415:   if (bs < 0 || bs == v->map->bs) return(0);
1417:   PetscLayoutSetBlockSize(v->map,bs);
1418:   v->bstash.bs = bs; /* use the same blocksize for the vec's block-stash */
1419:   return(0);
1420: }

1424: /*@
1425:    VecGetBlockSize - Gets the blocksize for the vector, i.e. what is used for VecSetValuesBlocked()
1426:    and VecSetValuesBlockedLocal().

1428:    Not Collective

1430:    Input Parameter:
1431: .  v - the vector

1433:    Output Parameter:
1434: .  bs - the blocksize

1436:    Notes:
1437:    All vectors obtained by VecDuplicate() inherit the same blocksize.

1439:    Level: advanced

1441: .seealso: VecSetValuesBlocked(), VecSetLocalToGlobalMappingBlock(), VecSetBlockSize()

1443:    Concepts: vector^block size
1444:    Concepts: block^vector

1446: @*/
1447: PetscErrorCode  VecGetBlockSize(Vec v,PetscInt *bs)
1448: {

1454:   PetscLayoutGetBlockSize(v->map,bs);
1455:   return(0);
1456: }

1460: /*@C
1461:    VecSetOptionsPrefix - Sets the prefix used for searching for all
1462:    Vec options in the database.

1464:    Logically Collective on Vec

1466:    Input Parameter:
1467: +  v - the Vec context
1468: -  prefix - the prefix to prepend to all option names

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

1474:    Level: advanced

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

1478: .seealso: VecSetFromOptions()
1479: @*/
1480: PetscErrorCode  VecSetOptionsPrefix(Vec v,const char prefix[])
1481: {

1486:   PetscObjectSetOptionsPrefix((PetscObject)v,prefix);
1487:   return(0);
1488: }

1492: /*@C
1493:    VecAppendOptionsPrefix - Appends to the prefix used for searching for all
1494:    Vec options in the database.

1496:    Logically Collective on Vec

1498:    Input Parameters:
1499: +  v - the Vec context
1500: -  prefix - the prefix to prepend to all option names

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

1506:    Level: advanced

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

1510: .seealso: VecGetOptionsPrefix()
1511: @*/
1512: PetscErrorCode  VecAppendOptionsPrefix(Vec v,const char prefix[])
1513: {

1518:   PetscObjectAppendOptionsPrefix((PetscObject)v,prefix);
1519:   return(0);
1520: }

1524: /*@C
1525:    VecGetOptionsPrefix - Sets the prefix used for searching for all
1526:    Vec options in the database.

1528:    Not Collective

1530:    Input Parameter:
1531: .  v - the Vec context

1533:    Output Parameter:
1534: .  prefix - pointer to the prefix string used

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

1539:    Level: advanced

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

1543: .seealso: VecAppendOptionsPrefix()
1544: @*/
1545: PetscErrorCode  VecGetOptionsPrefix(Vec v,const char *prefix[])
1546: {

1551:   PetscObjectGetOptionsPrefix((PetscObject)v,prefix);
1552:   return(0);
1553: }

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

1560:    Collective on Vec

1562:    Input Parameters:
1563: .  v - the Vec context

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

1569:    Level: advanced

1571: .keywords: Vec, setup

1573: .seealso: VecCreate(), VecDestroy()
1574: @*/
1575: PetscErrorCode  VecSetUp(Vec v)
1576: {
1577:   PetscMPIInt    size;

1582:   if (!((PetscObject)v)->type_name) {
1583:     MPI_Comm_size(PetscObjectComm((PetscObject)v), &size);
1584:     if (size == 1) {
1585:       VecSetType(v, VECSEQ);
1586:     } else {
1587:       VecSetType(v, VECMPI);
1588:     }
1589:   }
1590:   return(0);
1591: }

1593: /*
1594:     These currently expose the PetscScalar/PetscReal in updating the
1595:     cached norm. If we push those down into the implementation these
1596:     will become independent of PetscScalar/PetscReal
1597: */

1601: /*@
1602:    VecCopy - Copies a vector. y <- x

1604:    Logically Collective on Vec

1606:    Input Parameter:
1607: .  x - the vector

1609:    Output Parameter:
1610: .  y - the copy

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

1616:    Level: beginner

1618: .seealso: VecDuplicate()
1619: @*/
1620: PetscErrorCode  VecCopy(Vec x,Vec y)
1621: {
1622:   PetscBool      flgs[4];
1623:   PetscReal      norms[4] = {0.0,0.0,0.0,0.0};
1625:   PetscInt       i;

1632:   if (x == y) return(0);
1633:   if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1634:   if (x->map->n != y->map->n) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths %d != %d", x->map->n, y->map->n);

1636: #if !defined(PETSC_USE_MIXED_PRECISION)
1637:   for (i=0; i<4; i++) {
1638:     PetscObjectComposedDataGetReal((PetscObject)x,NormIds[i],norms[i],flgs[i]);
1639:   }
1640: #endif

1642:   PetscLogEventBegin(VEC_Copy,x,y,0,0);
1643: #if defined(PETSC_USE_MIXED_PRECISION)
1644:   extern PetscErrorCode VecGetArray(Vec,double**);
1645:   extern PetscErrorCode VecRestoreArray(Vec,double**);
1646:   extern PetscErrorCode VecGetArray(Vec,float**);
1647:   extern PetscErrorCode VecRestoreArray(Vec,float**);
1648:   extern PetscErrorCode VecGetArrayRead(Vec,const double**);
1649:   extern PetscErrorCode VecRestoreArrayRead(Vec,const double**);
1650:   extern PetscErrorCode VecGetArrayRead(Vec,const float**);
1651:   extern PetscErrorCode VecRestoreArrayRead(Vec,const float**);
1652:   if ((((PetscObject)x)->precision == PETSC_PRECISION_SINGLE) && (((PetscObject)y)->precision == PETSC_PRECISION_DOUBLE)) {
1653:     PetscInt    i,n;
1654:     const float *xx;
1655:     double      *yy;
1656:     VecGetArrayRead(x,&xx);
1657:     VecGetArray(y,&yy);
1658:     VecGetLocalSize(x,&n);
1659:     for (i=0; i<n; i++) yy[i] = xx[i];
1660:     VecRestoreArrayRead(x,&xx);
1661:     VecRestoreArray(y,&yy);
1662:   } else if ((((PetscObject)x)->precision == PETSC_PRECISION_DOUBLE) && (((PetscObject)y)->precision == PETSC_PRECISION_SINGLE)) {
1663:     PetscInt     i,n;
1664:     float        *yy;
1665:     const double *xx;
1666:     VecGetArrayRead(x,&xx);
1667:     VecGetArray(y,&yy);
1668:     VecGetLocalSize(x,&n);
1669:     for (i=0; i<n; i++) yy[i] = (float) xx[i];
1670:     VecRestoreArrayRead(x,&xx);
1671:     VecRestoreArray(y,&yy);
1672:   } else {
1673:     (*x->ops->copy)(x,y);
1674:   }
1675: #else
1676:   (*x->ops->copy)(x,y);
1677: #endif

1679:   PetscObjectStateIncrease((PetscObject)y);
1680: #if !defined(PETSC_USE_MIXED_PRECISION)
1681:   for (i=0; i<4; i++) {
1682:     if (flgs[i]) {
1683:       PetscObjectComposedDataSetReal((PetscObject)y,NormIds[i],norms[i]);
1684:     }
1685:   }
1686: #endif

1688:   PetscLogEventEnd(VEC_Copy,x,y,0,0);
1689:   return(0);
1690: }

1694: /*@
1695:    VecSwap - Swaps the vectors x and y.

1697:    Logically Collective on Vec

1699:    Input Parameters:
1700: .  x, y  - the vectors

1702:    Level: advanced

1704:    Concepts: vector^swapping values

1706: @*/
1707: PetscErrorCode  VecSwap(Vec x,Vec y)
1708: {
1709:   PetscReal      normxs[4]={0.0,0.0,0.0,0.0},normys[4]={0.0,0.0,0.0,0.0};
1710:   PetscBool      flgxs[4],flgys[4];
1712:   PetscInt       i;

1720:   if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1721:   if (y->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1722:   if (x->map->N != y->map->N) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths");
1723:   if (x->map->n != y->map->n) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");

1725:   PetscLogEventBegin(VEC_Swap,x,y,0,0);
1726:   for (i=0; i<4; i++) {
1727:     PetscObjectComposedDataGetReal((PetscObject)x,NormIds[i],normxs[i],flgxs[i]);
1728:     PetscObjectComposedDataGetReal((PetscObject)y,NormIds[i],normys[i],flgys[i]);
1729:   }
1730:   (*x->ops->swap)(x,y);
1731:   PetscObjectStateIncrease((PetscObject)x);
1732:   PetscObjectStateIncrease((PetscObject)y);
1733:   for (i=0; i<4; i++) {
1734:     if (flgxs[i]) {
1735:       PetscObjectComposedDataSetReal((PetscObject)y,NormIds[i],normxs[i]);
1736:     }
1737:     if (flgys[i]) {
1738:       PetscObjectComposedDataSetReal((PetscObject)x,NormIds[i],normys[i]);
1739:     }
1740:   }
1741:   PetscLogEventEnd(VEC_Swap,x,y,0,0);
1742:   return(0);
1743: }

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

1750:   Collective on VecStash

1752:   Input Parameters:
1753: + obj   - the VecStash object
1754: . prefix - prefix to use for viewing, or NULL to use prefix of 'mat'
1755: - optionname - option to activate viewing

1757:   Level: intermediate

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

1761: */
1762: PetscErrorCode VecStashViewFromOptions(Vec obj,const char prefix[],const char optionname[])
1763: {
1764:   PetscErrorCode    ierr;
1765:   PetscViewer       viewer;
1766:   PetscBool         flg;
1767:   static PetscBool  incall = PETSC_FALSE;
1768:   PetscViewerFormat format;

1771:   if (incall) return(0);
1772:   incall = PETSC_TRUE;
1773:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)obj),prefix,optionname,&viewer,&format,&flg);
1774:   if (flg) {
1775:     PetscViewerPushFormat(viewer,format);
1776:     VecStashView(obj,viewer);
1777:     PetscViewerPopFormat(viewer);
1778:     PetscViewerDestroy(&viewer);
1779:   }
1780:   incall = PETSC_FALSE;
1781:   return(0);
1782: }

1786: /*@
1787:    VecStashView - Prints the entries in the vector stash and block stash.

1789:    Collective on Vec

1791:    Input Parameters:
1792: +  v - the vector
1793: -  viewer - the viewer

1795:    Level: advanced

1797:    Concepts: vector^stash
1798:    Concepts: stash^vector

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

1802: @*/
1803: PetscErrorCode  VecStashView(Vec v,PetscViewer viewer)
1804: {
1806:   PetscMPIInt    rank;
1807:   PetscInt       i,j;
1808:   PetscBool      match;
1809:   VecStash       *s;
1810:   PetscScalar    val;


1817:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&match);
1818:   if (!match) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_SUP,"Stash viewer only works with ASCII viewer not %s\n",((PetscObject)v)->type_name);
1819:   PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
1820:   MPI_Comm_rank(PetscObjectComm((PetscObject)v),&rank);
1821:   s    = &v->bstash;

1823:   /* print block stash */
1824:   PetscViewerASCIISynchronizedAllow(viewer,PETSC_TRUE);
1825:   PetscViewerASCIISynchronizedPrintf(viewer,"[%d]Vector Block stash size %D block size %D\n",rank,s->n,s->bs);
1826:   for (i=0; i<s->n; i++) {
1827:     PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Element %D ",rank,s->idx[i]);
1828:     for (j=0; j<s->bs; j++) {
1829:       val = s->array[i*s->bs+j];
1830: #if defined(PETSC_USE_COMPLEX)
1831:       PetscViewerASCIISynchronizedPrintf(viewer,"(%18.16e %18.16e) ",PetscRealPart(val),PetscImaginaryPart(val));
1832: #else
1833:       PetscViewerASCIISynchronizedPrintf(viewer,"%18.16e ",val);
1834: #endif
1835:     }
1836:     PetscViewerASCIISynchronizedPrintf(viewer,"\n");
1837:   }
1838:   PetscViewerFlush(viewer);

1840:   s = &v->stash;

1842:   /* print basic stash */
1843:   PetscViewerASCIISynchronizedPrintf(viewer,"[%d]Vector stash size %D\n",rank,s->n);
1844:   for (i=0; i<s->n; i++) {
1845:     val = s->array[i];
1846: #if defined(PETSC_USE_COMPLEX)
1847:     PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Element %D (%18.16e %18.16e) ",rank,s->idx[i],PetscRealPart(val),PetscImaginaryPart(val));
1848: #else
1849:     PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Element %D %18.16e\n",rank,s->idx[i],val);
1850: #endif
1851:   }
1852:   PetscViewerFlush(viewer);
1853:   PetscViewerASCIISynchronizedAllow(viewer,PETSC_FALSE);

1855:   PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
1856:   return(0);
1857: }

1861: PetscErrorCode PetscOptionsVec(const char key[],const char text[],const char man[],Vec v,PetscBool *set)
1862: {
1863:   PetscInt       i,N,rstart,rend;
1865:   PetscScalar    *xx;
1866:   PetscReal      *xreal;
1867:   PetscBool      iset;

1870:   VecGetOwnershipRange(v,&rstart,&rend);
1871:   VecGetSize(v,&N);
1872:   PetscCalloc1(N,&xreal);
1873:   PetscOptionsRealArray(key,text,man,xreal,&N,&iset);
1874:   if (iset) {
1875:     VecGetArray(v,&xx);
1876:     for (i=rstart; i<rend; i++) xx[i-rstart] = xreal[i];
1877:     VecRestoreArray(v,&xx);
1878:   }
1879:   PetscFree(xreal);
1880:   if (set) *set = iset;
1881:   return(0);
1882: }

1886: /*@
1887:    VecGetLayout - get PetscLayout describing vector layout

1889:    Not Collective

1891:    Input Arguments:
1892: .  x - the vector

1894:    Output Arguments:
1895: .  map - the layout

1897:    Level: developer

1899: .seealso: VecGetSizes(), VecGetOwnershipRange(), VecGetOwnershipRanges()
1900: @*/
1901: PetscErrorCode VecGetLayout(Vec x,PetscLayout *map)
1902: {

1906:   *map = x->map;
1907:   return(0);
1908: }

1912: /*@
1913:    VecSetLayout - set PetscLayout describing vector layout

1915:    Not Collective

1917:    Input Arguments:
1918: +  x - the vector
1919: -  map - the layout

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

1924:    Level: developer

1926: .seealso: VecGetLayout(), VecGetSizes(), VecGetOwnershipRange(), VecGetOwnershipRanges()
1927: @*/
1928: PetscErrorCode VecSetLayout(Vec x,PetscLayout map)
1929: {

1934:   PetscLayoutReference(map,&x->map);
1935:   return(0);
1936: }