Actual source code: vector.c

  1: #define PETSCVEC_DLL

  3: /*
  4:      Provides the interface functions for all vector operations.
  5:    These are the vector functions the user calls.
  6: */
 7:  #include vecimpl.h

  9: /* Logging support */
 10: PetscCookie PETSCVEC_DLLEXPORT VEC_COOKIE = 0;
 11: PetscEvent  VEC_View = 0, VEC_Max = 0, VEC_Min = 0, VEC_DotBarrier = 0, VEC_Dot = 0, VEC_MDotBarrier = 0, VEC_MDot = 0, VEC_TDot = 0;
 12: PetscEvent  VEC_Norm = 0, VEC_Normalize = 0, VEC_Scale = 0, VEC_Copy = 0, VEC_Set = 0, VEC_AXPY = 0, VEC_AYPX = 0, VEC_WAXPY = 0;
 13: PetscEvent  VEC_MTDot = 0, VEC_NormBarrier = 0, VEC_MAXPY = 0, VEC_Swap = 0, VEC_AssemblyBegin = 0, VEC_ScatterBegin = 0, VEC_ScatterEnd = 0;
 14: PetscEvent  VEC_AssemblyEnd = 0, VEC_PointwiseMult = 0, VEC_SetValues = 0, VEC_Load = 0, VEC_ScatterBarrier = 0;
 15: PetscEvent  VEC_SetRandom = 0, VEC_ReduceArithmetic = 0, VEC_ReduceBarrier = 0, VEC_ReduceCommunication = 0;

 17: /* ugly globals for VecSetValue() and VecSetValueLocal() */
 18: PetscInt    PETSCVEC_DLLEXPORT VecSetValue_Row = 0;
 19: PetscScalar PETSCVEC_DLLEXPORT VecSetValue_Value = 0.0;

 21: /*@
 22:   VecZeroEntries - puts a 0.0 in each element of a vector

 24:   Collective on Vec

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

 29:   Level: beginner

 31: .keywords: Vec, set, options, database
 32: .seealso: VecCreate(), VecPrintHelp(), VecSetOptionsPrefix(), VecSet(), VecSetValues()
 33: @*/
 34: PetscErrorCode PETSCVEC_DLLEXPORT VecZeroEntries (Vec vec)
 35: {
 38:   VecSet(vec,0.0);
 39:   return(0);
 40: }

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

 48:   Collective on Vec

 50:   Input Parameter:
 51: . vec - The vector

 53:   Level: intermediate

 55: .keywords: Vec, set, options, database, type
 56: .seealso: VecSetFromOptions(), VecSetType()
 57: */
 58: static PetscErrorCode VecSetTypeFromOptions_Private(Vec vec)
 59: {
 60:   PetscTruth     opt;
 61:   const char     *defaultType;
 62:   char           typeName[256];
 63:   PetscMPIInt    size;

 67:   if (vec->type_name) {
 68:     defaultType = vec->type_name;
 69:   } else {
 70:     MPI_Comm_size(vec->comm, &size);
 71:     if (size > 1) {
 72:       defaultType = VECMPI;
 73:     } else {
 74:       defaultType = VECSEQ;
 75:     }
 76:   }

 78:   if (!VecRegisterAllCalled) {
 79:     VecRegisterAll(PETSC_NULL);
 80:   }
 81:   PetscOptionsList("-vec_type","Vector type","VecSetType",VecList,defaultType,typeName,256,&opt);
 82:   if (opt) {
 83:     VecSetType(vec, typeName);
 84:   } else {
 85:     VecSetType(vec, defaultType);
 86:   }
 87:   return(0);
 88: }

 92: /*@C
 93:   VecSetFromOptions - Configures the vector from the options database.

 95:   Collective on Vec

 97:   Input Parameter:
 98: . vec - The vector

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

103:   Level: beginner

105:   Concepts: vectors^setting options
106:   Concepts: vectors^setting type

108: .keywords: Vec, set, options, database
109: .seealso: VecCreate(), VecPrintHelp(), VecSetOptionsPrefix()
110: @*/
111: PetscErrorCode PETSCVEC_DLLEXPORT VecSetFromOptions(Vec vec)
112: {
113:   PetscTruth     opt;


119:   PetscOptionsBegin(vec->comm, vec->prefix, "Vector options", "Vec");

121:   /* Handle generic vector options */
122:   PetscOptionsHasName(PETSC_NULL, "-help", &opt);
123:   if (opt) {
124:     VecPrintHelp(vec);
125:   }

127:   /* Handle vector type options */
128:   VecSetTypeFromOptions_Private(vec);

130:   /* Handle specific vector options */
131:   if (vec->ops->setfromoptions) {
132:     (*vec->ops->setfromoptions)(vec);
133:   }
134:   PetscOptionsEnd();

136:   VecViewFromOptions(vec, vec->name);
137:   return(0);
138: }

142: /*@
143:   VecPrintHelp - Prints some options for the Vec.

145:   Input Parameter:
146: . vec - The vector

148:   Options Database Keys:
149: $  -help, -h

151:   Level: intermediate

153: .keywords: Vec, help
154: .seealso: VecSetFromOptions()
155: @*/
156: PetscErrorCode PETSCVEC_DLLEXPORT VecPrintHelp(Vec vec)
157: {
160:   return(0);
161: }

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

168:   Collective on Vec

170:   Input Parameters:
171: + v - the vector
172: . n - the local size (or PETSC_DECIDE to have it set)
173: - N - the global size (or PETSC_DECIDE)

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

179:   Level: intermediate

181: .seealso: VecGetSize(), PetscSplitOwnership()
182: @*/
183: PetscErrorCode PETSCVEC_DLLEXPORT VecSetSizes(Vec v, PetscInt n, PetscInt N)
184: {
187:   if (N > 0 && n > N) SETERRQ2(PETSC_ERR_ARG_INCOMP,"Local size %D cannot be larger than global size %D",n,N);
188:   if ((v->n >= 0 || v->N >= 0) && (v->n != n || v->N != N)) SETERRQ4(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->n,v->N);
189:   v->n = n;
190:   v->N = N;
191:   return(0);
192: }

196: /*@
197:    VecSetBlockSize - Sets the blocksize for future calls to VecSetValuesBlocked()
198:    and VecSetValuesBlockedLocal().

200:    Collective on Vec

202:    Input Parameter:
203: +  v - the vector
204: -  bs - the blocksize

206:    Notes:
207:    All vectors obtained by VecDuplicate() inherit the same blocksize.

209:    Level: advanced

211: .seealso: VecSetValuesBlocked(), VecSetLocalToGlobalMappingBlocked(), VecGetBlockSize()

213:   Concepts: block size^vectors
214: @*/
215: PetscErrorCode PETSCVEC_DLLEXPORT VecSetBlockSize(Vec v,PetscInt bs)
216: {
219:   if (bs <= 0) bs = 1;
220:   if (bs == v->bs) return(0);
221:   if (v->bs != -1) SETERRQ2(PETSC_ERR_ARG_WRONGSTATE,"Cannot reset blocksize. Current size %D new %D",v->bs,bs);
222:   if (v->N != -1 && v->N % bs) SETERRQ2(PETSC_ERR_ARG_OUTOFRANGE,"Vector length not divisible by blocksize %D %D",v->N,bs);
223:   if (v->n != -1 && v->n % bs) SETERRQ2(PETSC_ERR_ARG_OUTOFRANGE,"Local vector length not divisible by blocksize %D %D\n\
224:    Try setting blocksize before setting the vector type",v->n,bs);
225: 
226:   v->bs        = bs;
227:   v->bstash.bs = bs; /* use the same blocksize for the vec's block-stash */
228:   return(0);
229: }

233: /*@
234:    VecGetBlockSize - Gets the blocksize for the vector, i.e. what is used for VecSetValuesBlocked()
235:    and VecSetValuesBlockedLocal().

237:    Collective on Vec

239:    Input Parameter:
240: .  v - the vector

242:    Output Parameter:
243: .  bs - the blocksize

245:    Notes:
246:    All vectors obtained by VecDuplicate() inherit the same blocksize.

248:    Level: advanced

250: .seealso: VecSetValuesBlocked(), VecSetLocalToGlobalMappingBlocked(), VecSetBlockSize()

252:    Concepts: vector^block size
253:    Concepts: block^vector

255: @*/
256: PetscErrorCode PETSCVEC_DLLEXPORT VecGetBlockSize(Vec v,PetscInt *bs)
257: {
261:   *bs = v->bs;
262:   return(0);
263: }

267: /*@
268:    VecValid - Checks whether a vector object is valid.

270:    Not Collective

272:    Input Parameter:
273: .  v - the object to check

275:    Output Parameter:
276:    flg - flag indicating vector status, either
277:    PETSC_TRUE if vector is valid, or PETSC_FALSE otherwise.

279:    Level: developer

281: @*/
282: PetscErrorCode PETSCVEC_DLLEXPORT VecValid(Vec v,PetscTruth *flg)
283: {
286:   if (!v)                           *flg = PETSC_FALSE;
287:   else if (v->cookie != VEC_COOKIE) *flg = PETSC_FALSE;
288:   else                              *flg = PETSC_TRUE;
289:   return(0);
290: }

294: /*@C
295:    VecSetOptionsPrefix - Sets the prefix used for searching for all 
296:    Vec options in the database.

298:    Collective on Vec

300:    Input Parameter:
301: +  v - the Vec context
302: -  prefix - the prefix to prepend to all option names

304:    Notes:
305:    A hyphen (-) must NOT be given at the beginning of the prefix name.
306:    The first character of all runtime options is AUTOVECICALLY the hyphen.

308:    Level: advanced

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

312: .seealso: VecSetFromOptions()
313: @*/
314: PetscErrorCode PETSCVEC_DLLEXPORT VecSetOptionsPrefix(Vec v,const char prefix[])
315: {

320:   PetscObjectSetOptionsPrefix((PetscObject)v,prefix);
321:   return(0);
322: }

326: /*@C
327:    VecAppendOptionsPrefix - Appends to the prefix used for searching for all 
328:    Vec options in the database.

330:    Collective on Vec

332:    Input Parameters:
333: +  v - the Vec context
334: -  prefix - the prefix to prepend to all option names

336:    Notes:
337:    A hyphen (-) must NOT be given at the beginning of the prefix name.
338:    The first character of all runtime options is AUTOVECICALLY the hyphen.

340:    Level: advanced

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

344: .seealso: VecGetOptionsPrefix()
345: @*/
346: PetscErrorCode PETSCVEC_DLLEXPORT VecAppendOptionsPrefix(Vec v,const char prefix[])
347: {
349: 
352:   PetscObjectAppendOptionsPrefix((PetscObject)v,prefix);
353:   return(0);
354: }

358: /*@C
359:    VecGetOptionsPrefix - Sets the prefix used for searching for all 
360:    Vec options in the database.

362:    Not Collective

364:    Input Parameter:
365: .  A - the Vec context

367:    Output Parameter:
368: .  prefix - pointer to the prefix string used

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

373:    Level: advanced

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

377: .seealso: VecAppendOptionsPrefix()
378: @*/
379: PetscErrorCode PETSCVEC_DLLEXPORT VecGetOptionsPrefix(Vec v,const char *prefix[])
380: {

385:   PetscObjectGetOptionsPrefix((PetscObject)v,prefix);
386:   return(0);
387: }

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

394:    Collective on Vec

396:    Input Parameters:
397: .  v - the Vec context

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

403:    Level: advanced

405: .keywords: Vec, setup

407: .seealso: VecCreate(), VecDestroy()
408: @*/
409: PetscErrorCode PETSCVEC_DLLEXPORT VecSetUp(Vec v)
410: {

415:   VecSetFromOptions(v);
416:   return(0);
417: }

421: /*@
422:    VecDot - Computes the vector dot product.

424:    Collective on Vec

426:    Input Parameters:
427: .  x, y - the vectors

429:    Output Parameter:
430: .  alpha - the dot product

432:    Performance Issues:
433: +    per-processor memory bandwidth
434: .    interprocessor latency
435: -    work load inbalance that causes certain processes to arrive much earlier than
436:      others

438:    Notes for Users of Complex Numbers:
439:    For complex vectors, VecDot() computes 
440: $     val = (x,y) = y^H x,
441:    where y^H denotes the conjugate transpose of y.

443:    Use VecTDot() for the indefinite form
444: $     val = (x,y) = y^T x,
445:    where y^T denotes the transpose of y.

447:    Level: intermediate

449:    Concepts: inner product
450:    Concepts: vector^inner product

452: .seealso: VecMDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd()
453: @*/
454: PetscErrorCode PETSCVEC_DLLEXPORT VecDot(Vec x,Vec y,PetscScalar *val)
455: {

465:   if (x->N != y->N) SETERRQ(PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths");
466:   if (x->n != y->n) SETERRQ(PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");

468:   PetscLogEventBarrierBegin(VEC_DotBarrier,x,y,0,0,x->comm);
469:   (*x->ops->dot)(x,y,val);
470:   PetscLogEventBarrierEnd(VEC_DotBarrier,x,y,0,0,x->comm);
471:   return(0);
472: }

476: /*@
477:    VecNorm  - Computes the vector norm.

479:    Collective on Vec

481:    Input Parameters:
482: +  x - the vector
483: -  type - one of NORM_1, NORM_2, NORM_INFINITY.  Also available
484:           NORM_1_AND_2, which computes both norms and stores them
485:           in a two element array.

487:    Output Parameter:
488: .  val - the norm 

490:    Notes:
491: $     NORM_1 denotes sum_i |x_i|
492: $     NORM_2 denotes sqrt(sum_i (x_i)^2)
493: $     NORM_INFINITY denotes max_i |x_i|

495:    Level: intermediate

497:    Performance Issues:
498: +    per-processor memory bandwidth
499: .    interprocessor latency
500: -    work load inbalance that causes certain processes to arrive much earlier than
501:      others

503:    Compile Option:
504:    PETSC_HAVE_SLOW_BLAS_NORM2 will cause a C (loop unrolled) version of the norm to be used, rather
505:  than the BLAS. This should probably only be used when one is using the FORTRAN BLAS routines 
506:  (as opposed to vendor provided) because the FORTRAN BLAS NRM2() routine is very slow. 

508:    Concepts: norm
509:    Concepts: vector^norm

511: .seealso: VecDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), 
512:           VecNormBegin(), VecNormEnd()

514: @*/
515: PetscErrorCode PETSCVEC_DLLEXPORT VecNorm(Vec x,NormType type,PetscReal *val)
516: {
517:   PetscTruth     flg;
519:   PetscInt       type_id;


526:   /*
527:    * Cached data?
528:    */
529:   if (type!=NORM_1_AND_2) {
530:     VecNormComposedDataID(type,&type_id);
531:     PetscObjectComposedDataGetReal((PetscObject)x,type_id,*val,flg);
532:     if (flg) return(0);
533:   }
534: 

536:   PetscLogEventBarrierBegin(VEC_NormBarrier,x,0,0,0,x->comm);
537:   (*x->ops->norm)(x,type,val);
538:   PetscLogEventBarrierEnd(VEC_NormBarrier,x,0,0,0,x->comm);

540:   if (type!=NORM_1_AND_2) {
541:     PetscObjectComposedDataSetReal((PetscObject)x,type_id,*val);
542:   }

544:   return(0);
545: }

549: PetscErrorCode PETSCVEC_DLLEXPORT VecNormComposedDataID(NormType type,PetscInt *type_id)
550: {
551:   static PetscInt id_norm1=-1,id_norm2=-1,id_normInf=-1,id_normF=-1,id_norm12=-1;
552:   PetscErrorCode  ierr;

555:   switch (type) {
556:   case NORM_1 :
557:     if (id_norm1==-1) {
558:       PetscObjectComposedDataRegister(&id_norm1);}
559:     *type_id = id_norm1; break;
560:   case NORM_2 :
561:     if (id_norm2==-1) {
562:       PetscObjectComposedDataRegister(&id_norm2);}
563:     *type_id = id_norm2; break;
564:   case NORM_1_AND_2 :
565:     /* we don't handle this one yet */
566:     if (id_norm1==-1) {
567:       PetscObjectComposedDataRegister(&id_norm1);}
568:     if (id_norm2==-1) {
569:       PetscObjectComposedDataRegister(&id_norm2);}
570:     *type_id = id_norm12; break;
571:   case NORM_INFINITY :
572:     if (id_normInf==-1) {
573:       PetscObjectComposedDataRegister(&id_normInf);}
574:     *type_id = id_normInf; break;
575:   case NORM_FROBENIUS :
576:     if (id_normF==-1) {
577:       PetscObjectComposedDataRegister(&id_normF);}
578:     *type_id = id_normF; break;
579:   }
580:   return(0);
581: }

585: /*@
586:    VecNormalize - Normalizes a vector by 2-norm. 

588:    Collective on Vec

590:    Input Parameters:
591: +  x - the vector

593:    Output Parameter:
594: .  x - the normalized vector
595: -  val - the vector norm before normalization

597:    Level: intermediate

599:    Concepts: vector^normalizing
600:    Concepts: normalizing^vector

602: @*/
603: PetscErrorCode PETSCVEC_DLLEXPORT VecNormalize(Vec x,PetscReal *val)
604: {

611:   PetscLogEventBegin(VEC_Normalize,x,0,0,0);
612:   VecNorm(x,NORM_2,val);
613:   if (!*val) {
614:     PetscLogInfo((x,"VecNormalize:Vector of zero norm can not be normalized; Returning only the zero norm\n"));
615:   } else {
616:     PetscScalar tmp = 1.0/(*val);
617:     VecScale(x,tmp);
618:   }

620:   PetscLogEventEnd(VEC_Normalize,x,0,0,0);
621:   return(0);
622: }

626: /*@C
627:    VecMax - Determines the maximum vector component and its location.

629:    Collective on Vec

631:    Input Parameter:
632: .  x - the vector

634:    Output Parameters:
635: +  val - the maximum component
636: -  p - the location of val

638:    Notes:
639:    Returns the value PETSC_MIN and p = -1 if the vector is of length 0.

641:    Level: intermediate

643:    Concepts: maximum^of vector
644:    Concepts: vector^maximum value

646: .seealso: VecNorm(), VecMin()
647: @*/
648: PetscErrorCode PETSCVEC_DLLEXPORT VecMax(Vec x,PetscInt *p,PetscReal *val)
649: {

656:   PetscLogEventBegin(VEC_Max,x,0,0,0);
657:   (*x->ops->max)(x,p,val);
658:   PetscLogEventEnd(VEC_Max,x,0,0,0);
659:   return(0);
660: }

664: /*@
665:    VecMin - Determines the minimum vector component and its location.

667:    Collective on Vec

669:    Input Parameters:
670: .  x - the vector

672:    Output Parameter:
673: +  val - the minimum component
674: -  p - the location of val

676:    Level: intermediate

678:    Notes:
679:    Returns the value PETSC_MAX and p = -1 if the vector is of length 0.

681:    Concepts: minimum^of vector
682:    Concepts: vector^minimum entry

684: .seealso: VecMax()
685: @*/
686: PetscErrorCode PETSCVEC_DLLEXPORT VecMin(Vec x,PetscInt *p,PetscReal *val)
687: {

694:   PetscLogEventBegin(VEC_Min,x,0,0,0);
695:   (*x->ops->min)(x,p,val);
696:   PetscLogEventEnd(VEC_Min,x,0,0,0);
697:   return(0);
698: }

702: /*@
703:    VecTDot - Computes an indefinite vector dot product. That is, this
704:    routine does NOT use the complex conjugate.

706:    Collective on Vec

708:    Input Parameters:
709: .  x, y - the vectors

711:    Output Parameter:
712: .  val - the dot product

714:    Notes for Users of Complex Numbers:
715:    For complex vectors, VecTDot() computes the indefinite form
716: $     val = (x,y) = y^T x,
717:    where y^T denotes the transpose of y.

719:    Use VecDot() for the inner product
720: $     val = (x,y) = y^H x,
721:    where y^H denotes the conjugate transpose of y.

723:    Level: intermediate

725:    Concepts: inner product^non-Hermitian
726:    Concepts: vector^inner product
727:    Concepts: non-Hermitian inner product

729: .seealso: VecDot(), VecMTDot()
730: @*/
731: PetscErrorCode PETSCVEC_DLLEXPORT VecTDot(Vec x,Vec y,PetscScalar *val)
732: {

742:   if (x->N != y->N) SETERRQ(PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths");
743:   if (x->n != y->n) SETERRQ(PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");

745:   PetscLogEventBegin(VEC_TDot,x,y,0,0);
746:   (*x->ops->tdot)(x,y,val);
747:   PetscLogEventEnd(VEC_TDot,x,y,0,0);
748:   return(0);
749: }

753: /*@
754:    VecScale - Scales a vector. 

756:    Collective on Vec

758:    Input Parameters:
759: +  x - the vector
760: -  alpha - the scalar

762:    Output Parameter:
763: .  x - the scaled vector

765:    Note:
766:    For a vector with n components, VecScale() computes 
767: $      x[i] = alpha * x[i], for i=1,...,n.

769:    Level: intermediate

771:    Concepts: vector^scaling
772:    Concepts: scaling^vector

774: @*/
775: PetscErrorCode PETSCVEC_DLLEXPORT VecScale (Vec x, PetscScalar alpha)
776: {
777:   PetscReal      scale,norm1=0.0,norm2=0.0,normInf=0.0,normF=0.0;
778:   PetscTruth     flg1,flg2,flgInf,flgF;
780:   PetscInt       type_id1,type_id2,type_idInf,type_idF;

785:   if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
786:   PetscLogEventBegin(VEC_Scale,x,0,0,0);
787:   (*x->ops->scale)(x,alpha);

789:   /*
790:    * Update cached data
791:    */
792:   /* see if we have cached norms */
793:   /* 1 */
794:   VecNormComposedDataID(NORM_1,&type_id1);
795:   PetscObjectComposedDataGetReal((PetscObject)x,type_id1,norm1,flg1);
796:   /* 2 */
797:   VecNormComposedDataID(NORM_2,&type_id2);
798:   PetscObjectComposedDataGetReal((PetscObject)x,type_id2,norm2,flg2);
799:   /* inf */
800:   VecNormComposedDataID(NORM_INFINITY,&type_idInf);
801:   PetscObjectComposedDataGetReal((PetscObject)x,type_idInf,normInf,flgInf);
802:   /* frobenius */
803:   VecNormComposedDataID(NORM_FROBENIUS,&type_idF);
804:   PetscObjectComposedDataGetReal((PetscObject)x,type_idF,normF,flgF);

806:   /* in general we consider this object touched */
807:   PetscObjectStateIncrease((PetscObject)x);

809:   /* however, norms can be simply updated */
810:   scale = PetscAbsScalar(alpha);
811:   /* 1 */
812:   if (flg1) {
813:     PetscObjectComposedDataSetReal((PetscObject)x,type_id1,scale*norm1);
814:   }
815:   /* 2 */
816:   if (flg2) {
817:     PetscObjectComposedDataSetReal((PetscObject)x,type_id2,scale*norm2);
818:   }
819:   /* inf */
820:   if (flgInf) {
821:     PetscObjectComposedDataSetReal((PetscObject)x,type_idInf,scale*normInf);
822:   }
823:   /* frobenius */
824:   if (flgF) {
825:     PetscObjectComposedDataSetReal((PetscObject)x,type_idF,scale*normF);
826:   }

828:   PetscLogEventEnd(VEC_Scale,x,0,0,0);
829:   return(0);
830: }

834: /*@
835:    VecCopy - Copies a vector. 

837:    Collective on Vec

839:    Input Parameter:
840: .  x - the vector

842:    Output Parameter:
843: .  y - the copy

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

849:    Level: beginner

851: .seealso: VecDuplicate()
852: @*/
853: PetscErrorCode PETSCVEC_DLLEXPORT VecCopy(Vec x,Vec y)
854: {
855:   PetscTruth     flg;
856:   PetscReal      norm=0.0;
858:   PetscInt       type_id;

866:   if (x->N != y->N) SETERRQ(PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths");
867:   if (x->n != y->n) SETERRQ(PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");

869:   PetscLogEventBegin(VEC_Copy,x,y,0,0);
870:   (*x->ops->copy)(x,y);

872:   /*
873:    * Update cached data
874:    */
875:   /* in general we consider this object touched */
876:   PetscObjectStateIncrease((PetscObject)y);
877:   /* however, norms can be simply copied over */
878:   /* 2 */
879:   VecNormComposedDataID(NORM_2,&type_id);
880:   PetscObjectComposedDataGetReal((PetscObject)x,type_id,norm,flg);
881:   if (flg) {
882:     PetscObjectComposedDataSetReal((PetscObject)y,type_id,norm);
883:   }
884:   /* 1 */
885:   VecNormComposedDataID(NORM_1,&type_id);
886:   PetscObjectComposedDataGetReal((PetscObject)x,type_id,norm,flg);
887:   if (flg) {
888:     PetscObjectComposedDataSetReal((PetscObject)y,type_id,norm);
889:   }
890:   /* inf */
891:   VecNormComposedDataID(NORM_INFINITY,&type_id);
892:   PetscObjectComposedDataGetReal((PetscObject)x,type_id,norm,flg);
893:   if (flg) {
894:     PetscObjectComposedDataSetReal((PetscObject)y,type_id,norm);
895:   }
896:   /* frobenius */
897:   VecNormComposedDataID(NORM_FROBENIUS,&type_id);
898:   PetscObjectComposedDataGetReal((PetscObject)x,type_id,norm,flg);
899:   if (flg) {
900:     PetscObjectComposedDataSetReal((PetscObject)y,type_id,norm);
901:   }

903:   PetscLogEventEnd(VEC_Copy,x,y,0,0);
904:   return(0);
905: }

909: /*@
910:    VecSet - Sets all components of a vector to a single scalar value. 

912:    Collective on Vec

914:    Input Parameters:
915: +  x  - the vector
916: -  alpha - the scalar

918:    Output Parameter:
919: .  x  - the vector

921:    Note:
922:    For a vector of dimension n, VecSet() computes
923: $     x[i] = alpha, for i=1,...,n,
924:    so that all vector entries then equal the identical
925:    scalar value, alpha.  Use the more general routine
926:    VecSetValues() to set different vector entries.

928:    You CANNOT call this after you have called VecSetValues() but before you call 
929:    VecAssemblyBegin/End().

931:    Level: beginner

933: .seealso VecSetValues(), VecSetValuesBlocked(), VecSetRandom()

935:    Concepts: vector^setting to constant

937: @*/
938: PetscErrorCode PETSCVEC_DLLEXPORT VecSet(Vec x,PetscScalar alpha)
939: {
940:   PetscReal      val;
942:   PetscInt       type_id;

947:   if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"You cannot call this after you have called VecSetValues() but\n before you have called VecAssemblyBegin/End()");
948: #if defined (PETSC_USE_DEBUG)
949:  {
950:    PetscReal alpha_local,alpha_max;
951:    alpha_local = PetscAbsScalar(alpha);
952:    MPI_Allreduce(&alpha_local,&alpha_max,1,MPIU_REAL,MPI_MAX,x->comm);
953:    if (alpha_local != alpha_max) SETERRQ(PETSC_ERR_ARG_WRONG,"Same value should be used across all processors");
954:  }
955: #endif

957:   PetscLogEventBegin(VEC_Set,x,0,0,0);
958:   (*x->ops->set)(x,alpha);
959:   PetscLogEventEnd(VEC_Set,x,0,0,0);

961:   /*
962:    * Update cached data
963:    */
964:   /* in general we consider this object touched */
965:   PetscObjectStateIncrease((PetscObject)x);
966:   /* however, norms can be simply set */
967:   /* 1 */
968:   val = PetscAbsScalar(alpha);
969:   VecNormComposedDataID(NORM_1,&type_id);
970:   PetscObjectComposedDataSetReal((PetscObject)x,type_id,x->N * val);
971:   /* inf */
972:   VecNormComposedDataID(NORM_INFINITY,&type_id);
973:   PetscObjectComposedDataSetReal((PetscObject)x,type_id,val);
974:   /* 2 */
975:   val = sqrt((double)x->N) * val;
976:   VecNormComposedDataID(NORM_2,&type_id);
977:   PetscObjectComposedDataSetReal((PetscObject)x,type_id,val);
978:   /* frobenius */
979:   VecNormComposedDataID(NORM_FROBENIUS,&type_id);
980:   PetscObjectComposedDataSetReal((PetscObject)x,type_id,val);

982:   return(0);
983: }

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

990:    Collective on Vec

992:    Input Parameters:
993: +  x  - the vector
994: -  rctx - the random number context, formed by PetscRandomCreate(), or PETSC_NULL and
995:           it will create one internally.

997:    Output Parameter:
998: .  x  - the vector

1000:    Example of Usage:
1001: .vb
1002:      PetscRandomCreate(PETSC_COMM_WORLD,RANDOM_DEFAULT,&rctx);
1003:      VecSetRandom(x,rctx);
1004:      PetscRandomDestroy(rctx);
1005: .ve

1007:    Level: intermediate

1009:    Concepts: vector^setting to random
1010:    Concepts: random^vector

1012: .seealso: VecSet(), VecSetValues(), PetscRandomCreate(), PetscRandomDestroy()
1013: @*/
1014: PetscErrorCode PETSCVEC_DLLEXPORT VecSetRandom(Vec x,PetscRandom rctx)
1015: {
1017:   PetscRandom    randObj = PETSC_NULL;

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

1025:   if (!rctx) {
1026:     MPI_Comm    comm;
1027:     PetscObjectGetComm((PetscObject)x,&comm);
1028:     PetscRandomCreate(comm,RANDOM_DEFAULT,&randObj);
1029:     rctx = randObj;
1030:   }

1032:   PetscLogEventBegin(VEC_SetRandom,x,rctx,0,0);
1033:   (*x->ops->setrandom)(x,rctx);
1034:   PetscLogEventEnd(VEC_SetRandom,x,rctx,0,0);
1035: 
1036:   if (randObj) {
1037:     PetscRandomDestroy(randObj);
1038:   }
1039:   PetscObjectStateIncrease((PetscObject)x);
1040:   return(0);
1041: }

1045: /*@
1046:    VecAXPY - Computes y = alpha x + y. 

1048:    Collective on Vec

1050:    Input Parameters:
1051: +  alpha - the scalar
1052: -  x, y  - the vectors

1054:    Output Parameter:
1055: .  y - output vector

1057:    Level: intermediate

1059:    Concepts: vector^BLAS
1060:    Concepts: BLAS

1062: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY()
1063: @*/
1064: PetscErrorCode PETSCVEC_DLLEXPORT VecAXPY(Vec y,PetscScalar alpha,Vec x)
1065: {

1074:   if (x->N != y->N) SETERRQ(PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths");
1075:   if (x->n != y->n) SETERRQ(PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");

1077:   PetscLogEventBegin(VEC_AXPY,x,y,0,0);
1078:   (*y->ops->axpy)(y,alpha,x);
1079:   PetscLogEventEnd(VEC_AXPY,x,y,0,0);
1080:   PetscObjectStateIncrease((PetscObject)y);
1081:   return(0);
1082: }

1086: /*@
1087:    VecAXPBY - Computes y = alpha x + beta y. 

1089:    Collective on Vec

1091:    Input Parameters:
1092: +  alpha,beta - the scalars
1093: -  x, y  - the vectors

1095:    Output Parameter:
1096: .  y - output vector

1098:    Level: intermediate

1100:    Concepts: BLAS
1101:    Concepts: vector^BLAS

1103: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPY()
1104: @*/
1105: PetscErrorCode PETSCVEC_DLLEXPORT VecAXPBY(Vec y,PetscScalar alpha,PetscScalar beta,Vec x)
1106: {

1115:   if (x->N != y->N) SETERRQ(PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths");
1116:   if (x->n != y->n) SETERRQ(PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");

1118:   PetscLogEventBegin(VEC_AXPY,x,y,0,0);
1119:   (*y->ops->axpby)(y,alpha,beta,x);
1120:   PetscLogEventEnd(VEC_AXPY,x,y,0,0);
1121:   PetscObjectStateIncrease((PetscObject)y);
1122:   return(0);
1123: }

1127: /*@
1128:    VecAYPX - Computes y = x + alpha y.

1130:    Collective on Vec

1132:    Input Parameters:
1133: +  alpha - the scalar
1134: -  x, y  - the vectors

1136:    Output Parameter:
1137: .  y - output vector

1139:    Level: intermediate

1141:    Concepts: vector^BLAS
1142:    Concepts: BLAS

1144: .seealso: VecAXPY(), VecWAXPY()
1145: @*/
1146: PetscErrorCode PETSCVEC_DLLEXPORT VecAYPX(Vec y,PetscScalar alpha,Vec x)
1147: {

1156:   if (x->N != y->N) SETERRQ(PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths");
1157:   if (x->n != y->n) SETERRQ(PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");

1159:   PetscLogEventBegin(VEC_AYPX,x,y,0,0);
1160:    (*y->ops->aypx)(y,alpha,x);
1161:   PetscLogEventEnd(VEC_AYPX,x,y,0,0);
1162:   PetscObjectStateIncrease((PetscObject)y);
1163:   return(0);
1164: }

1168: /*@
1169:    VecSwap - Swaps the vectors x and y.

1171:    Collective on Vec

1173:    Input Parameters:
1174: .  x, y  - the vectors

1176:    Level: advanced

1178:    Concepts: vector^swapping values

1180: @*/
1181: PetscErrorCode PETSCVEC_DLLEXPORT VecSwap(Vec x,Vec y)
1182: {
1183:   PetscReal      norm1x=0.0,norm2x=0.0,normInfx=0.0,normFx=0.0;
1184:   PetscReal      norm1y=0.0,norm2y=0.0,normInfy=0.0,normFy=0.0;
1185:   PetscTruth     flg1x,flg2x,flgInfx,flgFx;
1186:   PetscTruth     flg1y,flg2y,flgInfy,flgFy;
1187:   PetscInt       type_id1,type_id2,type_idInf,type_idF;

1196:   if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1197:   if (y->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
1198:   if (x->N != y->N) SETERRQ(PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths");
1199:   if (x->n != y->n) SETERRQ(PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");

1201:   PetscLogEventBegin(VEC_Swap,x,y,0,0);

1203:   /* See if we have cached norms */
1204:   /* 1 */
1205:   VecNormComposedDataID(NORM_1,&type_id1);
1206:   PetscObjectComposedDataGetReal((PetscObject)x,type_id1,norm1x,flg1x);
1207:   PetscObjectComposedDataGetReal((PetscObject)y,type_id1,norm1y,flg1y);
1208:   /* 2 */
1209:   VecNormComposedDataID(NORM_2,&type_id2);
1210:   PetscObjectComposedDataGetReal((PetscObject)x,type_id2,norm2x,flg2x);
1211:   PetscObjectComposedDataGetReal((PetscObject)y,type_id2,norm2y,flg2y);
1212:   /* inf */
1213:   VecNormComposedDataID(NORM_INFINITY,&type_idInf);
1214:   PetscObjectComposedDataGetReal((PetscObject)x,type_idInf,normInfx,flgInfx);
1215:   PetscObjectComposedDataGetReal((PetscObject)y,type_idInf,normInfy,flgInfy);
1216:   /* frobenius */
1217:   VecNormComposedDataID(NORM_FROBENIUS,&type_idF);
1218:   PetscObjectComposedDataGetReal((PetscObject)x,type_idF,normFx,flgFx);
1219:   PetscObjectComposedDataGetReal((PetscObject)y,type_idF,normFy,flgFy);

1221:   /* Do the actual swap */
1222:   (*x->ops->swap)(x,y);
1223:   PetscObjectStateIncrease((PetscObject)x);
1224:   PetscObjectStateIncrease((PetscObject)y);

1226:   /* Swap any cached norms */
1227:   /* 1 */
1228:   if (flg1x) {
1229:     PetscObjectComposedDataSetReal((PetscObject)y,type_id1,norm1x);
1230:   }
1231:   if (flg1y) {
1232:     PetscObjectComposedDataSetReal((PetscObject)x,type_id1,norm1y);
1233:   }
1234:   /* 2 */
1235:   if (flg2x) {
1236:     PetscObjectComposedDataSetReal((PetscObject)y,type_id2,norm2x);
1237:   }
1238:   if (flg2y) {
1239:     PetscObjectComposedDataSetReal((PetscObject)x,type_id2,norm2y);
1240:   }
1241:   /* inf */
1242:   if (flgInfx) {
1243:     PetscObjectComposedDataSetReal((PetscObject)y,type_idInf,normInfx);
1244:   }
1245:   if (flgInfy) {
1246:     PetscObjectComposedDataSetReal((PetscObject)x,type_idInf,normInfy);
1247:   }
1248:   /* frobenius */
1249:   if (flgFx) {
1250:     PetscObjectComposedDataSetReal((PetscObject)y,type_idF,normFx);
1251:   }
1252:   if (flgFy) {
1253:     PetscObjectComposedDataSetReal((PetscObject)x,type_idF,normFy);
1254:   }

1256:   PetscLogEventEnd(VEC_Swap,x,y,0,0);

1258:   return(0);
1259: }

1263: /*@
1264:    VecWAXPY - Computes w = alpha x + y.

1266:    Collective on Vec

1268:    Input Parameters:
1269: +  alpha - the scalar
1270: -  x, y  - the vectors

1272:    Output Parameter:
1273: .  w - the result

1275:    Level: intermediate

1277:    Concepts: vector^BLAS
1278:    Concepts: BLAS

1280: .seealso: VecAXPY(), VecAYPX()
1281: @*/
1282: PetscErrorCode PETSCVEC_DLLEXPORT VecWAXPY(Vec w,PetscScalar alpha,Vec x,Vec y)
1283: {

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

1298:   PetscLogEventBegin(VEC_WAXPY,x,y,w,0);
1299:    (*w->ops->waxpy)(w,alpha,x,y);
1300:   PetscLogEventEnd(VEC_WAXPY,x,y,w,0);
1301:   PetscObjectStateIncrease((PetscObject)w);
1302:   return(0);
1303: }

1307: /*@
1308:    VecPointwiseMult - Computes the componentwise multiplication w = x*y.

1310:    Collective on Vec

1312:    Input Parameters:
1313: .  x, y  - the vectors

1315:    Output Parameter:
1316: .  w - the result

1318:    Level: advanced

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

1322:    Concepts: vector^pointwise multiply

1324: .seealso: VecPointwiseDivide(), VecPointwiseMax(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
1325: @*/
1326: PetscErrorCode PETSCVEC_DLLEXPORT VecPointwiseMult(Vec w, Vec x,Vec y)
1327: {

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

1342:   PetscLogEventBegin(VEC_PointwiseMult,x,y,w,0);
1343:   (*w->ops->pointwisemult)(w,x,y);
1344:   PetscLogEventEnd(VEC_PointwiseMult,x,y,w,0);
1345:   PetscObjectStateIncrease((PetscObject)w);
1346:   return(0);
1347: }

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

1354:    Collective on Vec

1356:    Input Parameters:
1357: .  x, y  - the vectors

1359:    Output Parameter:
1360: .  w - the result

1362:    Level: advanced

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

1367:    Concepts: vector^pointwise multiply

1369: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
1370: @*/
1371: PetscErrorCode PETSCVEC_DLLEXPORT VecPointwiseMax(Vec w,Vec x,Vec y)
1372: {

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

1387:   (*w->ops->pointwisemax)(w,x,y);
1388:   PetscObjectStateIncrease((PetscObject)w);
1389:   return(0);
1390: }


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

1398:    Collective on Vec

1400:    Input Parameters:
1401: .  x, y  - the vectors

1403:    Output Parameter:
1404: .  w - the result

1406:    Level: advanced

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

1411:    Concepts: vector^pointwise multiply

1413: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
1414: @*/
1415: PetscErrorCode PETSCVEC_DLLEXPORT VecPointwiseMin(Vec w,Vec x,Vec y)
1416: {

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

1431:   (*w->ops->pointwisemin)(w,x,y);
1432:   PetscObjectStateIncrease((PetscObject)w);
1433:   return(0);
1434: }

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

1441:    Collective on Vec

1443:    Input Parameters:
1444: .  x, y  - the vectors

1446:    Output Parameter:
1447: .  w - the result

1449:    Level: advanced

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

1453:    Concepts: vector^pointwise multiply

1455: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMin(), VecPointwiseMax(), VecMaxPointwiseDivide()
1456: @*/
1457: PetscErrorCode PETSCVEC_DLLEXPORT VecPointwiseMaxAbs(Vec w,Vec x,Vec y)
1458: {

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

1473:   (*w->ops->pointwisemaxabs)(w,x,y);
1474:   PetscObjectStateIncrease((PetscObject)w);
1475:   return(0);
1476: }

1480: /*@
1481:    VecPointwiseDivide - Computes the componentwise division w = x/y.

1483:    Collective on Vec

1485:    Input Parameters:
1486: .  x, y  - the vectors

1488:    Output Parameter:
1489: .  w - the result

1491:    Level: advanced

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

1495:    Concepts: vector^pointwise divide

1497: .seealso: VecPointwiseMult(), VecPointwiseMax(), VecPointwiseMin(), VecPointwiseMaxAbs(), VecMaxPointwiseDivide()
1498: @*/
1499: PetscErrorCode PETSCVEC_DLLEXPORT VecPointwiseDivide(Vec w,Vec x,Vec y)
1500: {

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

1515:   (*w->ops->pointwisedivide)(w,x,y);
1516:   PetscObjectStateIncrease((PetscObject)w);
1517:   return(0);
1518: }

1522: /*@
1523:    VecMaxPointwiseDivide - Computes the maximum of the componentwise division w = abs(x/y).

1525:    Collective on Vec

1527:    Input Parameters:
1528: .  x, y  - the vectors

1530:    Output Parameter:
1531: .  max - the result

1533:    Level: advanced

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

1537: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMax(), VecPointwiseMin(), VecPointwiseMaxAbs()
1538: @*/
1539: PetscErrorCode PETSCVEC_DLLEXPORT VecMaxPointwiseDivide(Vec x,Vec y,PetscReal *max)
1540: {

1550:   if (x->N != y->N) SETERRQ(PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths");
1551:   if (x->n != y->n) SETERRQ(PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");

1553:   (*x->ops->maxpointwisedivide)(x,y,max);
1554:   return(0);
1555: }

1559: /*@C
1560:    VecDuplicate - Creates a new vector of the same type as an existing vector.

1562:    Collective on Vec

1564:    Input Parameters:
1565: .  v - a vector to mimic

1567:    Output Parameter:
1568: .  newv - location to put new vector

1570:    Notes:
1571:    VecDuplicate() does not copy the vector, but rather allocates storage
1572:    for the new vector.  Use VecCopy() to copy a vector.

1574:    Use VecDestroy() to free the space. Use VecDuplicateVecs() to get several
1575:    vectors. 

1577:    Level: beginner

1579: .seealso: VecDestroy(), VecDuplicateVecs(), VecCreate(), VecCopy()
1580: @*/
1581: PetscErrorCode PETSCVEC_DLLEXPORT VecDuplicate(Vec x,Vec *newv)
1582: {

1589:   (*x->ops->duplicate)(x,newv);
1590:   PetscObjectStateIncrease((PetscObject)*newv);
1591:   return(0);
1592: }

1596: /*@C
1597:    VecDestroy - Destroys a vector.

1599:    Collective on Vec

1601:    Input Parameters:
1602: .  v  - the vector

1604:    Level: beginner

1606: .seealso: VecDuplicate(), VecDestroyVecs()
1607: @*/
1608: PetscErrorCode PETSCVEC_DLLEXPORT VecDestroy(Vec v)
1609: {

1614:   if (--v->refct > 0) return(0);
1615:   /* destroy the internal part */
1616:   if (v->ops->destroy) {
1617:     (*v->ops->destroy)(v);
1618:   }
1619:   /* destroy the external/common part */
1620:   if (v->mapping) {
1621:     ISLocalToGlobalMappingDestroy(v->mapping);
1622:   }
1623:   if (v->bmapping) {
1624:     ISLocalToGlobalMappingDestroy(v->bmapping);
1625:   }
1626:   if (v->map) {
1627:     PetscMapDestroy(v->map);
1628:   }
1629:   PetscHeaderDestroy(v);
1630:   return(0);
1631: }

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

1638:    Collective on Vec

1640:    Input Parameters:
1641: +  m - the number of vectors to obtain
1642: -  v - a vector to mimic

1644:    Output Parameter:
1645: .  V - location to put pointer to array of vectors

1647:    Notes:
1648:    Use VecDestroyVecs() to free the space. Use VecDuplicate() to form a single
1649:    vector.

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

1656:    Level: intermediate

1658: .seealso:  VecDestroyVecs(), VecDuplicate(), VecCreate(), VecDuplicateVecsF90()
1659: @*/
1660: PetscErrorCode PETSCVEC_DLLEXPORT VecDuplicateVecs(Vec v,PetscInt m,Vec *V[])
1661: {

1668:   (*v->ops->duplicatevecs)(v, m,V);
1669:   return(0);
1670: }

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

1677:    Collective on Vec

1679:    Input Parameters:
1680: +  vv - pointer to array of vector pointers
1681: -  m - the number of vectors previously obtained

1683:    Fortran Note:
1684:    The Fortran interface is slightly different from that given below.
1685:    See the Fortran chapter of the users manual and 
1686:    petsc/src/vec/examples for details.

1688:    Level: intermediate

1690: .seealso: VecDuplicateVecs(), VecDestroyVecsf90()
1691: @*/
1692: PetscErrorCode PETSCVEC_DLLEXPORT VecDestroyVecs(Vec vv[],PetscInt m)
1693: {

1700:   (*(*vv)->ops->destroyvecs)(vv,m);
1701:   return(0);
1702: }

1706: /*@
1707:    VecSetValues - Inserts or adds values into certain locations of a vector. 

1709:    Not Collective

1711:    Input Parameters:
1712: +  x - vector to insert in
1713: .  ni - number of elements to add
1714: .  ix - indices where to add
1715: .  y - array of values
1716: -  iora - either INSERT_VALUES or ADD_VALUES, where
1717:    ADD_VALUES adds values to any existing entries, and
1718:    INSERT_VALUES replaces existing entries with new values

1720:    Notes: 
1721:    VecSetValues() sets x[ix[i]] = y[i], for i=0,...,ni-1.

1723:    Calls to VecSetValues() with the INSERT_VALUES and ADD_VALUES 
1724:    options cannot be mixed without intervening calls to the assembly
1725:    routines.

1727:    These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd() 
1728:    MUST be called after all calls to VecSetValues() have been completed.

1730:    VecSetValues() uses 0-based indices in Fortran as well as in C.

1732:    Negative indices may be passed in ix, these rows are 
1733:    simply ignored. This allows easily inserting element load matrices
1734:    with homogeneous Dirchlet boundary conditions that you don't want represented
1735:    in the vector.

1737:    Level: beginner

1739:    Concepts: vector^setting values

1741: .seealso:  VecAssemblyBegin(), VecAssemblyEnd(), VecSetValuesLocal(),
1742:            VecSetValue(), VecSetValuesBlocked(), InsertMode, INSERT_VALUES, ADD_VALUES, VecGetValues()
1743: @*/
1744: PetscErrorCode PETSCVEC_DLLEXPORT VecSetValues(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
1745: {

1753:   PetscLogEventBegin(VEC_SetValues,x,0,0,0);
1754:   (*x->ops->setvalues)(x,ni,ix,y,iora);
1755:   PetscLogEventEnd(VEC_SetValues,x,0,0,0);
1756:   PetscObjectStateIncrease((PetscObject)x);
1757:   return(0);
1758: }

1762: /*@
1763:    VecGetValues - Gets values from certain locations of a vector. Currently 
1764:           can only get values on the same processor

1766:     Collective on Vec
1767:  
1768:    Input Parameters:
1769: +  x - vector to get values from
1770: .  ni - number of elements to get
1771: -  ix - indices where to get them from (in global 1d numbering)

1773:    Output Parameter:
1774: .   y - array of values

1776:    Notes: 
1777:    The user provides the allocated array y; it is NOT allocated in this routine

1779:    VecGetValues() gets y[i] = x[ix[i]], for i=0,...,ni-1.

1781:    VecAssemblyBegin() and VecAssemblyEnd()  MUST be called before calling this

1783:    VecGetValues() uses 0-based indices in Fortran as well as in C.

1785:    Level: beginner

1787:    Concepts: vector^getting values

1789: .seealso:  VecAssemblyBegin(), VecAssemblyEnd(), VecGetValuesLocal(),
1790:            VecGetValue(), VecGetValuesBlocked(), InsertMode, INSERT_VALUES, ADD_VALUES, VecSetValues()
1791: @*/
1792: PetscErrorCode PETSCVEC_DLLEXPORT VecGetValues(Vec x,PetscInt ni,const PetscInt ix[],PetscScalar y[])
1793: {

1801:   (*x->ops->getvalues)(x,ni,ix,y);
1802:   return(0);
1803: }

1807: /*@
1808:    VecSetValuesBlocked - Inserts or adds blocks of values into certain locations of a vector. 

1810:    Not Collective

1812:    Input Parameters:
1813: +  x - vector to insert in
1814: .  ni - number of blocks to add
1815: .  ix - indices where to add in block count, rather than element count
1816: .  y - array of values
1817: -  iora - either INSERT_VALUES or ADD_VALUES, where
1818:    ADD_VALUES adds values to any existing entries, and
1819:    INSERT_VALUES replaces existing entries with new values

1821:    Notes: 
1822:    VecSetValuesBlocked() sets x[bs*ix[i]+j] = y[bs*i+j], 
1823:    for j=0,...,bs, for i=0,...,ni-1. where bs was set with VecSetBlockSize().

1825:    Calls to VecSetValuesBlocked() with the INSERT_VALUES and ADD_VALUES 
1826:    options cannot be mixed without intervening calls to the assembly
1827:    routines.

1829:    These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd() 
1830:    MUST be called after all calls to VecSetValuesBlocked() have been completed.

1832:    VecSetValuesBlocked() uses 0-based indices in Fortran as well as in C.

1834:    Negative indices may be passed in ix, these rows are 
1835:    simply ignored. This allows easily inserting element load matrices
1836:    with homogeneous Dirchlet boundary conditions that you don't want represented
1837:    in the vector.

1839:    Level: intermediate

1841:    Concepts: vector^setting values blocked

1843: .seealso:  VecAssemblyBegin(), VecAssemblyEnd(), VecSetValuesBlockedLocal(),
1844:            VecSetValues()
1845: @*/
1846: PetscErrorCode PETSCVEC_DLLEXPORT VecSetValuesBlocked(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
1847: {

1855:   PetscLogEventBegin(VEC_SetValues,x,0,0,0);
1856:   (*x->ops->setvaluesblocked)(x,ni,ix,y,iora);
1857:   PetscLogEventEnd(VEC_SetValues,x,0,0,0);
1858:   PetscObjectStateIncrease((PetscObject)x);
1859:   return(0);
1860: }

1864: /*@
1865:    VecSetLocalToGlobalMapping - Sets a local numbering to global numbering used
1866:    by the routine VecSetValuesLocal() to allow users to insert vector entries
1867:    using a local (per-processor) numbering.

1869:    Collective on Vec

1871:    Input Parameters:
1872: +  x - vector
1873: -  mapping - mapping created with ISLocalToGlobalMappingCreate() or ISLocalToGlobalMappingCreateIS()

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

1878:    Level: intermediate

1880:    Concepts: vector^setting values with local numbering

1882: seealso:  VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetValuesLocal(),
1883:            VecSetLocalToGlobalMappingBlocked(), VecSetValuesBlockedLocal()
1884: @*/
1885: PetscErrorCode PETSCVEC_DLLEXPORT VecSetLocalToGlobalMapping(Vec x,ISLocalToGlobalMapping mapping)
1886: {


1893:   if (x->mapping) {
1894:     SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Mapping already set for vector");
1895:   }

1897:   if (x->ops->setlocaltoglobalmapping) {
1898:     (*x->ops->setlocaltoglobalmapping)(x,mapping);
1899:   } else {
1900:     x->mapping = mapping;
1901:     PetscObjectReference((PetscObject)mapping);
1902:   }
1903:   return(0);
1904: }

1908: /*@
1909:    VecSetLocalToGlobalMappingBlock - Sets a local numbering to global numbering used
1910:    by the routine VecSetValuesBlockedLocal() to allow users to insert vector entries
1911:    using a local (per-processor) numbering.

1913:    Collective on Vec

1915:    Input Parameters:
1916: +  x - vector
1917: -  mapping - mapping created with ISLocalToGlobalMappingCreate() or ISLocalToGlobalMappingCreateIS()

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

1922:    Level: intermediate

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

1926: .seealso:  VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetValuesLocal(),
1927:            VecSetLocalToGlobalMapping(), VecSetValuesBlockedLocal()
1928: @*/
1929: PetscErrorCode PETSCVEC_DLLEXPORT VecSetLocalToGlobalMappingBlock(Vec x,ISLocalToGlobalMapping mapping)
1930: {


1937:   if (x->bmapping) {
1938:     SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Mapping already set for vector");
1939:   }
1940:   x->bmapping = mapping;
1941:   PetscObjectReference((PetscObject)mapping);
1942:   return(0);
1943: }

1947: /*@
1948:    VecSetValuesLocal - Inserts or adds values into certain locations of a vector,
1949:    using a local ordering of the nodes. 

1951:    Not Collective

1953:    Input Parameters:
1954: +  x - vector to insert in
1955: .  ni - number of elements to add
1956: .  ix - indices where to add
1957: .  y - array of values
1958: -  iora - either INSERT_VALUES or ADD_VALUES, where
1959:    ADD_VALUES adds values to any existing entries, and
1960:    INSERT_VALUES replaces existing entries with new values

1962:    Level: intermediate

1964:    Notes: 
1965:    VecSetValuesLocal() sets x[ix[i]] = y[i], for i=0,...,ni-1.

1967:    Calls to VecSetValues() with the INSERT_VALUES and ADD_VALUES 
1968:    options cannot be mixed without intervening calls to the assembly
1969:    routines.

1971:    These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd() 
1972:    MUST be called after all calls to VecSetValuesLocal() have been completed.

1974:    VecSetValuesLocal() uses 0-based indices in Fortran as well as in C.

1976:    Concepts: vector^setting values with local numbering

1978: .seealso:  VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetLocalToGlobalMapping(),
1979:            VecSetValuesBlockedLocal()
1980: @*/
1981: PetscErrorCode PETSCVEC_DLLEXPORT VecSetValuesLocal(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
1982: {
1984:   PetscInt       lixp[128],*lix = lixp;


1992:   PetscLogEventBegin(VEC_SetValues,x,0,0,0);
1993:   if (!x->ops->setvalueslocal) {
1994:     if (!x->mapping) {
1995:       SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Local to global never set with VecSetLocalToGlobalMapping()");
1996:     }
1997:     if (ni > 128) {
1998:       PetscMalloc(ni*sizeof(PetscInt),&lix);
1999:     }
2000:     ISLocalToGlobalMappingApply(x->mapping,ni,(PetscInt*)ix,lix);
2001:     (*x->ops->setvalues)(x,ni,lix,y,iora);
2002:     if (ni > 128) {
2003:       PetscFree(lix);
2004:     }
2005:   } else {
2006:     (*x->ops->setvalueslocal)(x,ni,ix,y,iora);
2007:   }
2008:   PetscLogEventEnd(VEC_SetValues,x,0,0,0);
2009:   PetscObjectStateIncrease((PetscObject)x);
2010:   return(0);
2011: }

2015: /*@
2016:    VecSetValuesBlockedLocal - Inserts or adds values into certain locations of a vector,
2017:    using a local ordering of the nodes. 

2019:    Not Collective

2021:    Input Parameters:
2022: +  x - vector to insert in
2023: .  ni - number of blocks to add
2024: .  ix - indices where to add in block count, not element count
2025: .  y - array of values
2026: -  iora - either INSERT_VALUES or ADD_VALUES, where
2027:    ADD_VALUES adds values to any existing entries, and
2028:    INSERT_VALUES replaces existing entries with new values

2030:    Level: intermediate

2032:    Notes: 
2033:    VecSetValuesBlockedLocal() sets x[bs*ix[i]+j] = y[bs*i+j], 
2034:    for j=0,..bs-1, for i=0,...,ni-1, where bs has been set with VecSetBlockSize().

2036:    Calls to VecSetValuesBlockedLocal() with the INSERT_VALUES and ADD_VALUES 
2037:    options cannot be mixed without intervening calls to the assembly
2038:    routines.

2040:    These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd() 
2041:    MUST be called after all calls to VecSetValuesBlockedLocal() have been completed.

2043:    VecSetValuesBlockedLocal() uses 0-based indices in Fortran as well as in C.


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

2048: .seealso:  VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetValuesBlocked(), 
2049:            VecSetLocalToGlobalMappingBlocked()
2050: @*/
2051: PetscErrorCode PETSCVEC_DLLEXPORT VecSetValuesBlockedLocal(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
2052: {
2054:   PetscInt       lixp[128],*lix = lixp;

2061:   if (!x->bmapping) {
2062:     SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Local to global never set with VecSetLocalToGlobalMappingBlocked()");
2063:   }
2064:   if (ni > 128) {
2065:     PetscMalloc(ni*sizeof(PetscInt),&lix);
2066:   }

2068:   PetscLogEventBegin(VEC_SetValues,x,0,0,0);
2069:   ISLocalToGlobalMappingApply(x->bmapping,ni,(PetscInt*)ix,lix);
2070:   (*x->ops->setvaluesblocked)(x,ni,lix,y,iora);
2071:   PetscLogEventEnd(VEC_SetValues,x,0,0,0);
2072:   if (ni > 128) {
2073:     PetscFree(lix);
2074:   }
2075:   PetscObjectStateIncrease((PetscObject)x);
2076:   return(0);
2077: }

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

2085:    Collective on Vec

2087:    Input Parameter:
2088: .  vec - the vector

2090:    Level: beginner

2092:    Concepts: assembly^vectors

2094: .seealso: VecAssemblyEnd(), VecSetValues()
2095: @*/
2096: PetscErrorCode PETSCVEC_DLLEXPORT VecAssemblyBegin(Vec vec)
2097: {
2099:   PetscTruth     flg;


2105:   PetscOptionsHasName(vec->prefix,"-vec_view_stash",&flg);
2106:   if (flg) {
2107:     VecStashView(vec,PETSC_VIEWER_STDOUT_(vec->comm));
2108:   }

2110:   PetscLogEventBegin(VEC_AssemblyBegin,vec,0,0,0);
2111:   if (vec->ops->assemblybegin) {
2112:     (*vec->ops->assemblybegin)(vec);
2113:   }
2114:   PetscLogEventEnd(VEC_AssemblyBegin,vec,0,0,0);
2115:   PetscObjectStateIncrease((PetscObject)vec);
2116:   return(0);
2117: }

2121: /*@
2122:    VecAssemblyEnd - Completes assembling the vector.  This routine should
2123:    be called after VecAssemblyBegin().

2125:    Collective on Vec

2127:    Input Parameter:
2128: .  vec - the vector

2130:    Options Database Keys:
2131: +  -vec_view - Prints vector in ASCII format
2132: .  -vec_view_matlab - Prints vector in ASCII Matlab format to stdout
2133: .  -vec_view_matlab_file - Prints vector in Matlab format to matlaboutput.mat
2134: .  -vec_view_draw - Activates vector viewing using drawing tools
2135: .  -display <name> - Sets display name (default is host)
2136: .  -draw_pause <sec> - Sets number of seconds to pause after display
2137: -  -vec_view_socket - Activates vector viewing using a socket
2138:  
2139:    Level: beginner

2141: .seealso: VecAssemblyBegin(), VecSetValues()
2142: @*/
2143: PetscErrorCode PETSCVEC_DLLEXPORT VecAssemblyEnd(Vec vec)
2144: {
2146:   PetscTruth     flg;

2150:   PetscLogEventBegin(VEC_AssemblyEnd,vec,0,0,0);
2152:   if (vec->ops->assemblyend) {
2153:     (*vec->ops->assemblyend)(vec);
2154:   }
2155:   PetscLogEventEnd(VEC_AssemblyEnd,vec,0,0,0);
2156:   PetscOptionsBegin(vec->comm,vec->prefix,"Vector Options","Vec");
2157:     PetscOptionsName("-vec_view","Print vector to stdout","VecView",&flg);
2158:     if (flg) {
2159:       VecView(vec,PETSC_VIEWER_STDOUT_(vec->comm));
2160:     }
2161:     PetscOptionsName("-vec_view_matlab","Print vector to stdout in a format Matlab can read","VecView",&flg);
2162:     if (flg) {
2163:       PetscViewerPushFormat(PETSC_VIEWER_STDOUT_(vec->comm),PETSC_VIEWER_ASCII_MATLAB);
2164:       VecView(vec,PETSC_VIEWER_STDOUT_(vec->comm));
2165:       PetscViewerPopFormat(PETSC_VIEWER_STDOUT_(vec->comm));
2166:     }
2167: #if defined(PETSC_HAVE_MATLAB)
2168:     PetscOptionsName("-vec_view_matlab_file","Print vector to matlaboutput.mat format Matlab can read","VecView",&flg);
2169:     if (flg) {
2170:       VecView(vec,PETSC_VIEWER_MATLAB_(vec->comm));
2171:     }
2172: #endif
2173: #if defined(PETSC_USE_SOCKET_VIEWER)
2174:     PetscOptionsName("-vec_view_socket","Send vector to socket (can be read from matlab)","VecView",&flg);
2175:     if (flg) {
2176:       VecView(vec,PETSC_VIEWER_SOCKET_(vec->comm));
2177:       PetscViewerFlush(PETSC_VIEWER_SOCKET_(vec->comm));
2178:     }
2179: #endif
2180:     PetscOptionsName("-vec_view_binary","Save vector to file in binary format","VecView",&flg);
2181:     if (flg) {
2182:       VecView(vec,PETSC_VIEWER_BINARY_(vec->comm));
2183:       PetscViewerFlush(PETSC_VIEWER_BINARY_(vec->comm));
2184:     }
2185:   PetscOptionsEnd();
2186:   /* These invoke PetscDrawGetDraw which invokes PetscOptionsBegin/End, */
2187:   /* hence they should not be inside the above PetscOptionsBegin/End block. */
2188:   PetscOptionsHasName(vec->prefix,"-vec_view_draw",&flg);
2189:   if (flg) {
2190:     VecView(vec,PETSC_VIEWER_DRAW_(vec->comm));
2191:     PetscViewerFlush(PETSC_VIEWER_DRAW_(vec->comm));
2192:   }
2193:   PetscOptionsHasName(vec->prefix,"-vec_view_draw_lg",&flg);
2194:   if (flg) {
2195:     PetscViewerSetFormat(PETSC_VIEWER_DRAW_(vec->comm),PETSC_VIEWER_DRAW_LG);
2196:     VecView(vec,PETSC_VIEWER_DRAW_(vec->comm));
2197:     PetscViewerFlush(PETSC_VIEWER_DRAW_(vec->comm));
2198:   }
2199:   return(0);
2200: }


2205: /*@C
2206:    VecMTDot - Computes indefinite vector multiple dot products. 
2207:    That is, it does NOT use the complex conjugate.

2209:    Collective on Vec

2211:    Input Parameters:
2212: +  nv - number of vectors
2213: .  x - one vector
2214: -  y - array of vectors.  Note that vectors are pointers

2216:    Output Parameter:
2217: .  val - array of the dot products

2219:    Notes for Users of Complex Numbers:
2220:    For complex vectors, VecMTDot() computes the indefinite form
2221: $      val = (x,y) = y^T x,
2222:    where y^T denotes the transpose of y.

2224:    Use VecMDot() for the inner product
2225: $      val = (x,y) = y^H x,
2226:    where y^H denotes the conjugate transpose of y.

2228:    Level: intermediate

2230:    Concepts: inner product^multiple
2231:    Concepts: vector^multiple inner products

2233: .seealso: VecMDot(), VecTDot()
2234: @*/
2235: PetscErrorCode PETSCVEC_DLLEXPORT VecMTDot(PetscInt nv,Vec x,const Vec y[],PetscScalar *val)
2236: {

2247:   if (x->N != (*y)->N) SETERRQ(PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths");
2248:   if (x->n != (*y)->n) SETERRQ(PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");

2250:   PetscLogEventBegin(VEC_MTDot,x,*y,0,0);
2251:   (*x->ops->mtdot)(nv,x,y,val);
2252:   PetscLogEventEnd(VEC_MTDot,x,*y,0,0);
2253:   return(0);
2254: }

2258: /*@C
2259:    VecMDot - Computes vector multiple dot products. 

2261:    Collective on Vec

2263:    Input Parameters:
2264: +  nv - number of vectors
2265: .  x - one vector
2266: -  y - array of vectors. 

2268:    Output Parameter:
2269: .  val - array of the dot products

2271:    Notes for Users of Complex Numbers:
2272:    For complex vectors, VecMDot() computes 
2273: $     val = (x,y) = y^H x,
2274:    where y^H denotes the conjugate transpose of y.

2276:    Use VecMTDot() for the indefinite form
2277: $     val = (x,y) = y^T x,
2278:    where y^T denotes the transpose of y.

2280:    Level: intermediate

2282:    Concepts: inner product^multiple
2283:    Concepts: vector^multiple inner products

2285: .seealso: VecMTDot(), VecDot()
2286: @*/
2287: PetscErrorCode PETSCVEC_DLLEXPORT VecMDot(PetscInt nv,Vec x,const Vec y[],PetscScalar *val)
2288: {

2299:   if (x->N != (*y)->N) SETERRQ(PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths");
2300:   if (x->n != (*y)->n) SETERRQ(PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");

2302:   PetscLogEventBarrierBegin(VEC_MDotBarrier,x,*y,0,0,x->comm);
2303:   (*x->ops->mdot)(nv,x,y,val);
2304:   PetscLogEventBarrierEnd(VEC_MDotBarrier,x,*y,0,0,x->comm);
2305:   return(0);
2306: }

2310: /*@C
2311:    VecMAXPY - Computes y = y + sum alpha[j] x[j]

2313:    Collective on Vec

2315:    Input Parameters:
2316: +  nv - number of scalars and x-vectors
2317: .  alpha - array of scalars
2318: .  y - one vector
2319: -  x - array of vectors

2321:    Level: intermediate

2323:    Concepts: BLAS

2325: .seealso: VecAXPY(), VecWAXPY(), VecAYPX()
2326: @*/
2327: PetscErrorCode PETSCVEC_DLLEXPORT VecMAXPY(Vec y,PetscInt nv,const PetscScalar alpha[],Vec *x)
2328: {

2339:   if (y->N != (*x)->N) SETERRQ(PETSC_ERR_ARG_INCOMP,"Incompatible vector global lengths");
2340:   if (y->n != (*x)->n) SETERRQ(PETSC_ERR_ARG_INCOMP,"Incompatible vector local lengths");

2342:   PetscLogEventBegin(VEC_MAXPY,*x,y,0,0);
2343:   (*y->ops->maxpy)(y,nv,alpha,x);
2344:   PetscLogEventEnd(VEC_MAXPY,*x,y,0,0);
2345:   PetscObjectStateIncrease((PetscObject)y);
2346:   return(0);
2347: }

2349: /*MC
2350:    VecGetArray - Returns a pointer to a contiguous array that contains this 
2351:    processor's portion of the vector data. For the standard PETSc
2352:    vectors, VecGetArray() returns a pointer to the local data array and
2353:    does not use any copies. If the underlying vector data is not stored
2354:    in a contiquous array this routine will copy the data to a contiquous
2355:    array and return a pointer to that. You MUST call VecRestoreArray() 
2356:    when you no longer need access to the array.

2358:    Synopsis:
2359:    PetscErrorCode VecGetArray(Vec x,PetscScalar *a[])

2361:    Not Collective

2363:    Input Parameter:
2364: .  x - the vector

2366:    Output Parameter:
2367: .  a - location to put pointer to the array

2369:    Fortran Note:
2370:    This routine is used differently from Fortran 77
2371: $    Vec         x
2372: $    PetscScalar x_array(1)
2373: $    PetscOffset i_x
2374: $    PetscErrorCode ierr
2375: $       call VecGetArray(x,x_array,i_x,ierr)
2376: $
2377: $   Access first local entry in vector with
2378: $      value = x_array(i_x + 1)
2379: $
2380: $      ...... other code
2381: $       call VecRestoreArray(x,x_array,i_x,ierr)
2382:    For Fortran 90 see VecGetArrayF90()

2384:    See the Fortran chapter of the users manual and 
2385:    petsc/src/snes/examples/tutorials/ex5f.F for details.

2387:    Level: beginner

2389:    Concepts: vector^accessing local values

2391: .seealso: VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(), VecGetArray2d()
2392: M*/
2395: PetscErrorCode VecGetArray_Private(Vec x,PetscScalar *a[])
2396: {

2403:   (*x->ops->getarray)(x,a);
2404:   return(0);
2405: }


2410: /*@C
2411:    VecGetArrays - Returns a pointer to the arrays in a set of vectors
2412:    that were created by a call to VecDuplicateVecs().  You MUST call
2413:    VecRestoreArrays() when you no longer need access to the array.

2415:    Not Collective

2417:    Input Parameter:
2418: +  x - the vectors
2419: -  n - the number of vectors

2421:    Output Parameter:
2422: .  a - location to put pointer to the array

2424:    Fortran Note:
2425:    This routine is not supported in Fortran.

2427:    Level: intermediate

2429: .seealso: VecGetArray(), VecRestoreArrays()
2430: @*/
2431: PetscErrorCode PETSCVEC_DLLEXPORT VecGetArrays(const Vec x[],PetscInt n,PetscScalar **a[])
2432: {
2434:   PetscInt       i;
2435:   PetscScalar    **q;

2441:   if (n <= 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Must get at least one array n = %D",n);
2442:   PetscMalloc(n*sizeof(PetscScalar*),&q);
2443:   for (i=0; i<n; ++i) {
2444:     VecGetArray(x[i],&q[i]);
2445:   }
2446:   *a = q;
2447:   return(0);
2448: }

2452: /*@C
2453:    VecRestoreArrays - Restores a group of vectors after VecGetArrays()
2454:    has been called.

2456:    Not Collective

2458:    Input Parameters:
2459: +  x - the vector
2460: .  n - the number of vectors
2461: -  a - location of pointer to arrays obtained from VecGetArrays()

2463:    Notes:
2464:    For regular PETSc vectors this routine does not involve any copies. For
2465:    any special vectors that do not store local vector data in a contiguous
2466:    array, this routine will copy the data back into the underlying 
2467:    vector data structure from the arrays obtained with VecGetArrays().

2469:    Fortran Note:
2470:    This routine is not supported in Fortran.

2472:    Level: intermediate

2474: .seealso: VecGetArrays(), VecRestoreArray()
2475: @*/
2476: PetscErrorCode PETSCVEC_DLLEXPORT VecRestoreArrays(const Vec x[],PetscInt n,PetscScalar **a[])
2477: {
2479:   PetscInt       i;
2480:   PetscScalar    **q = *a;


2487:   for(i=0;i<n;++i) {
2488:     VecRestoreArray(x[i],&q[i]);
2489:  }
2490:   PetscFree(q);
2491:   return(0);
2492: }

2494: /*MC
2495:    VecRestoreArray - Restores a vector after VecGetArray() has been called.

2497:    Not Collective

2499:    Synopsis:
2500:    PetscErrorCode VecRestoreArray(Vec x,PetscScalar *a[])

2502:    Input Parameters:
2503: +  x - the vector
2504: -  a - location of pointer to array obtained from VecGetArray()

2506:    Level: beginner

2508:    Notes:
2509:    For regular PETSc vectors this routine does not involve any copies. For
2510:    any special vectors that do not store local vector data in a contiguous
2511:    array, this routine will copy the data back into the underlying 
2512:    vector data structure from the array obtained with VecGetArray().

2514:    This routine actually zeros out the a pointer. This is to prevent accidental
2515:    us of the array after it has been restored. If you pass null for a it will 
2516:    not zero the array pointer a.

2518:    Fortran Note:
2519:    This routine is used differently from Fortran 77
2520: $    Vec         x
2521: $    PetscScalar x_array(1)
2522: $    PetscOffset i_x
2523: $    PetscErrorCode ierr
2524: $       call VecGetArray(x,x_array,i_x,ierr)
2525: $
2526: $   Access first local entry in vector with
2527: $      value = x_array(i_x + 1)
2528: $
2529: $      ...... other code
2530: $       call VecRestoreArray(x,x_array,i_x,ierr)

2532:    See the Fortran chapter of the users manual and 
2533:    petsc/src/snes/examples/tutorials/ex5f.F for details.
2534:    For Fortran 90 see VecRestoreArrayF90()

2536: .seealso: VecGetArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(), VecRestoreArray2d()
2537: M*/
2540: PetscErrorCode VecRestoreArray_Private(Vec x,PetscScalar *a[])
2541: {

2548: #if defined(PETSC_USE_DEBUG)
2549:   CHKMEMQ;
2550: #endif
2551:   if (x->ops->restorearray) {
2552:     (*x->ops->restorearray)(x,a);
2553:   }
2554:   PetscObjectStateIncrease((PetscObject)x);
2555:   return(0);
2556: }

2558: #undef  __FUNCT__
2560: /*@
2561:   VecViewFromOptions - This function visualizes the vector based upon user options.

2563:   Collective on Vec

2565:   Input Parameters:
2566: . vec   - The vector
2567: . title - The title

2569:   Level: intermediate

2571: .keywords: Vec, view, options, database
2572: .seealso: VecSetFromOptions(), VecView()
2573: @*/
2574: PetscErrorCode PETSCVEC_DLLEXPORT VecViewFromOptions(Vec vec, char *title)
2575: {
2576:   PetscViewer    viewer;
2577:   PetscDraw      draw;
2578:   PetscTruth     opt;
2579:   char           *titleStr;
2580:   char           typeName[1024];
2581:   char           fileName[PETSC_MAX_PATH_LEN];
2582:   size_t         len;

2586:   PetscOptionsHasName(vec->prefix, "-vec_view", &opt);
2587:   if (opt) {
2588:     PetscOptionsGetString(vec->prefix, "-vec_view", typeName, 1024, &opt);
2589:     PetscStrlen(typeName, &len);
2590:     if (len > 0) {
2591:       PetscViewerCreate(vec->comm, &viewer);
2592:       PetscViewerSetType(viewer, typeName);
2593:       PetscOptionsGetString(vec->prefix, "-vec_view_file", fileName, 1024, &opt);
2594:       if (opt) {
2595:         PetscViewerSetFilename(viewer, fileName);
2596:       } else {
2597:         PetscViewerSetFilename(viewer, vec->name);
2598:       }
2599:       VecView(vec, viewer);
2600:       PetscViewerFlush(viewer);
2601:       PetscViewerDestroy(viewer);
2602:     } else {
2603:       VecView(vec, PETSC_VIEWER_STDOUT_(vec->comm));
2604:     }
2605:   }
2606:   PetscOptionsHasName(vec->prefix, "-vec_view_draw", &opt);
2607:   if (opt) {
2608:     PetscViewerDrawOpen(vec->comm, 0, 0, 0, 0, 300, 300, &viewer);
2609:     PetscViewerDrawGetDraw(viewer, 0, &draw);
2610:     if (title) {
2611:       titleStr = title;
2612:     } else {
2613:       PetscObjectName((PetscObject) vec);
2614:       titleStr = vec->name;
2615:     }
2616:     PetscDrawSetTitle(draw, titleStr);
2617:     VecView(vec, viewer);
2618:     PetscViewerFlush(viewer);
2619:     PetscDrawPause(draw);
2620:     PetscViewerDestroy(viewer);
2621:   }
2622:   return(0);
2623: }

2627: /*@C
2628:    VecView - Views a vector object. 

2630:    Collective on Vec

2632:    Input Parameters:
2633: +  v - the vector
2634: -  viewer - an optional visualization context

2636:    Notes:
2637:    The available visualization contexts include
2638: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
2639: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
2640:          output where only the first processor opens
2641:          the file.  All other processors send their 
2642:          data to the first processor to print. 

2644:    You can change the format the vector is printed using the 
2645:    option PetscViewerSetFormat().

2647:    The user can open alternative visualization contexts with
2648: +    PetscViewerASCIIOpen() - Outputs vector to a specified file
2649: .    PetscViewerBinaryOpen() - Outputs vector in binary to a
2650:          specified file; corresponding input uses VecLoad()
2651: .    PetscViewerDrawOpen() - Outputs vector to an X window display
2652: -    PetscViewerSocketOpen() - Outputs vector to Socket viewer

2654:    The user can call PetscViewerSetFormat() to specify the output
2655:    format of ASCII printed objects (when using PETSC_VIEWER_STDOUT_SELF,
2656:    PETSC_VIEWER_STDOUT_WORLD and PetscViewerASCIIOpen).  Available formats include
2657: +    PETSC_VIEWER_ASCII_DEFAULT - default, prints vector contents
2658: .    PETSC_VIEWER_ASCII_MATLAB - prints vector contents in Matlab format
2659: .    PETSC_VIEWER_ASCII_INDEX - prints vector contents, including indices of vector elements
2660: -    PETSC_VIEWER_ASCII_COMMON - prints vector contents, using a 
2661:          format common among all vector types

2663:    Level: beginner

2665:    Concepts: vector^printing
2666:    Concepts: vector^saving to disk

2668: .seealso: PetscViewerASCIIOpen(), PetscViewerDrawOpen(), PetscDrawLGCreate(),
2669:           PetscViewerSocketOpen(), PetscViewerBinaryOpen(), VecLoad(), PetscViewerCreate(),
2670:           PetscRealView(), PetscScalarView(), PetscIntView()
2671: @*/
2672: PetscErrorCode PETSCVEC_DLLEXPORT VecView(Vec vec,PetscViewer viewer)
2673: {
2674:   PetscErrorCode    ierr;
2675:   PetscViewerFormat format;

2680:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(vec->comm);
2683:   if (vec->stash.n || vec->bstash.n) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call VecAssemblyBegin/End() before viewing this vector");

2685:   /*
2686:      Check if default viewer has been overridden, but user request it anyways
2687:   */
2688:   PetscViewerGetFormat(viewer,&format);
2689:   if (vec->ops->viewnative && format == PETSC_VIEWER_NATIVE) {
2690:     PetscViewerPopFormat(viewer);
2691:     (*vec->ops->viewnative)(vec,viewer);
2692:     PetscViewerPushFormat(viewer,PETSC_VIEWER_NATIVE);
2693:   } else {
2694:     (*vec->ops->view)(vec,viewer);
2695:   }
2696:   return(0);
2697: }

2701: /*@
2702:    VecGetSize - Returns the global number of elements of the vector.

2704:    Not Collective

2706:    Input Parameter:
2707: .  x - the vector

2709:    Output Parameters:
2710: .  size - the global length of the vector

2712:    Level: beginner

2714:    Concepts: vector^local size

2716: .seealso: VecGetLocalSize()
2717: @*/
2718: PetscErrorCode PETSCVEC_DLLEXPORT VecGetSize(Vec x,PetscInt *size)
2719: {

2726:   (*x->ops->getsize)(x,size);
2727:   return(0);
2728: }

2732: /*@
2733:    VecGetLocalSize - Returns the number of elements of the vector stored 
2734:    in local memory. This routine may be implementation dependent, so use 
2735:    with care.

2737:    Not Collective

2739:    Input Parameter:
2740: .  x - the vector

2742:    Output Parameter:
2743: .  size - the length of the local piece of the vector

2745:    Level: beginner

2747:    Concepts: vector^size

2749: .seealso: VecGetSize()
2750: @*/
2751: PetscErrorCode PETSCVEC_DLLEXPORT VecGetLocalSize(Vec x,PetscInt *size)
2752: {

2759:   (*x->ops->getlocalsize)(x,size);
2760:   return(0);
2761: }

2765: /*@C
2766:    VecGetOwnershipRange - Returns the range of indices owned by 
2767:    this processor, assuming that the vectors are laid out with the
2768:    first n1 elements on the first processor, next n2 elements on the
2769:    second, etc.  For certain parallel layouts this range may not be 
2770:    well defined. 

2772:    Not Collective

2774:    Input Parameter:
2775: .  x - the vector

2777:    Output Parameters:
2778: +  low - the first local element, pass in PETSC_NULL if not interested
2779: -  high - one more than the last local element, pass in PETSC_NULL if not interested

2781:    Note:
2782:    The high argument is one more than the last element stored locally.

2784:    Fortran: PETSC_NULL_INTEGER should be used instead of PETSC_NULL

2786:    Level: beginner

2788:    Concepts: ownership^of vectors
2789:    Concepts: vector^ownership of elements

2791: @*/
2792: PetscErrorCode PETSCVEC_DLLEXPORT VecGetOwnershipRange(Vec x,PetscInt *low,PetscInt *high)
2793: {

2801:   PetscMapGetLocalRange(x->map,low,high);
2802:   return(0);
2803: }

2807: /*@C
2808:    VecGetPetscMap - Returns the map associated with the vector

2810:    Not Collective

2812:    Input Parameter:
2813: .  x - the vector

2815:    Output Parameters:
2816: .  map - the map

2818:    Level: developer

2820: @*/
2821: PetscErrorCode PETSCVEC_DLLEXPORT VecGetPetscMap(Vec x,PetscMap *map)
2822: {
2827:   *map = x->map;
2828:   return(0);
2829: }

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

2836:    Collective on Vec

2838:    Input Parameter:
2839: +  x - the vector
2840: -  op - the option

2842:    Supported Options:
2843: +     VEC_IGNORE_OFF_PROC_ENTRIES, which causes VecSetValues() to ignore 
2844:       entries destined to be stored on a seperate processor. This can be used
2845:       to eliminate the global reduction in the VecAssemblyXXXX() if you know 
2846:       that you have only used VecSetValues() to set local elements
2847: -   VEC_TREAT_OFF_PROC_ENTRIES restores the treatment of off processor entries.

2849:    Level: intermediate

2851: @*/
2852: PetscErrorCode PETSCVEC_DLLEXPORT VecSetOption(Vec x,VecOption op)
2853: {

2859:   if (x->ops->setoption) {
2860:     (*x->ops->setoption)(x,op);
2861:   }
2862:   return(0);
2863: }

2867: /* Default routines for obtaining and releasing; */
2868: /* may be used by any implementation */
2869: PetscErrorCode VecDuplicateVecs_Default(Vec w,PetscInt m,Vec *V[])
2870: {
2872:   PetscInt       i;

2877:   if (m <= 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"m must be > 0: m = %D",m);
2878:   PetscMalloc(m*sizeof(Vec*),V);
2879:   for (i=0; i<m; i++) {VecDuplicate(w,*V+i);}
2880:   return(0);
2881: }

2885: PetscErrorCode VecDestroyVecs_Default(Vec v[], PetscInt m)
2886: {
2888:   PetscInt       i;

2892:   if (m <= 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"m must be > 0: m = %D",m);
2893:   for (i=0; i<m; i++) {VecDestroy(v[i]);}
2894:   PetscFree(v);
2895:   return(0);
2896: }

2900: /*@
2901:    VecPlaceArray - Allows one to replace the array in a vector with an
2902:    array provided by the user. This is useful to avoid copying an array
2903:    into a vector.

2905:    Not Collective

2907:    Input Parameters:
2908: +  vec - the vector
2909: -  array - the array

2911:    Notes:
2912:    You can return to the original array with a call to VecResetArray()

2914:    Level: developer

2916: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecResetArray()

2918: @*/
2919: PetscErrorCode PETSCVEC_DLLEXPORT VecPlaceArray(Vec vec,const PetscScalar array[])
2920: {

2927:   if (vec->ops->placearray) {
2928:     (*vec->ops->placearray)(vec,array);
2929:   } else {
2930:     SETERRQ(PETSC_ERR_SUP,"Cannot place array in this type of vector");
2931:   }
2932:   PetscObjectStateIncrease((PetscObject)vec);
2933:   return(0);
2934: }

2938: /*@
2939:    VecResetArray - Resets a vector to use its default memory. Call this 
2940:    after the use of VecPlaceArray().

2942:    Not Collective

2944:    Input Parameters:
2945: .  vec - the vector

2947:    Level: developer

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

2951: @*/
2952: PetscErrorCode PETSCVEC_DLLEXPORT VecResetArray(Vec vec)
2953: {

2959:   if (vec->ops->resetarray) {
2960:     (*vec->ops->resetarray)(vec);
2961:   } else {
2962:     SETERRQ(PETSC_ERR_SUP,"Cannot reset array in this type of vector");
2963:   }
2964:   PetscObjectStateIncrease((PetscObject)vec);
2965:   return(0);
2966: }

2970: /*@C
2971:    VecReplaceArray - Allows one to replace the array in a vector with an
2972:    array provided by the user. This is useful to avoid copying an array
2973:    into a vector.

2975:    Not Collective

2977:    Input Parameters:
2978: +  vec - the vector
2979: -  array - the array

2981:    Notes:
2982:    This permanently replaces the array and frees the memory associated
2983:    with the old array.

2985:    The memory passed in MUST be obtained with PetscMalloc() and CANNOT be
2986:    freed by the user. It will be freed when the vector is destroy. 

2988:    Not supported from Fortran

2990:    Level: developer

2992: .seealso: VecGetArray(), VecRestoreArray(), VecPlaceArray(), VecResetArray()

2994: @*/
2995: PetscErrorCode PETSCVEC_DLLEXPORT VecReplaceArray(Vec vec,const PetscScalar array[])
2996: {

3002:   if (vec->ops->replacearray) {
3003:     (*vec->ops->replacearray)(vec,array);
3004:  } else {
3005:     SETERRQ(PETSC_ERR_SUP,"Cannot replace array in this type of vector");
3006:   }
3007:   PetscObjectStateIncrease((PetscObject)vec);
3008:   return(0);
3009: }

3011: /*MC
3012:     VecDuplicateVecsF90 - Creates several vectors of the same type as an existing vector
3013:     and makes them accessible via a Fortran90 pointer.

3015:     Synopsis:
3016:     VecDuplicateVecsF90(Vec x,int n,{Vec, pointer :: y(:)},integer ierr)

3018:     Collective on Vec

3020:     Input Parameters:
3021: +   x - a vector to mimic
3022: -   n - the number of vectors to obtain

3024:     Output Parameters:
3025: +   y - Fortran90 pointer to the array of vectors
3026: -   ierr - error code

3028:     Example of Usage: 
3029: .vb
3030:     Vec x
3031:     Vec, pointer :: y(:)
3032:     ....
3033:     call VecDuplicateVecsF90(x,2,y,ierr)
3034:     call VecSet(y(2),alpha,ierr)
3035:     call VecSet(y(2),alpha,ierr)
3036:     ....
3037:     call VecDestroyVecsF90(y,2,ierr)
3038: .ve

3040:     Notes:
3041:     Not yet supported for all F90 compilers

3043:     Use VecDestroyVecsF90() to free the space.

3045:     Level: beginner

3047: .seealso:  VecDestroyVecsF90(), VecDuplicateVecs()

3049: M*/

3051: /*MC
3052:     VecRestoreArrayF90 - Restores a vector to a usable state after a call to
3053:     VecGetArrayF90().

3055:     Synopsis:
3056:     VecRestoreArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)

3058:     Not collective

3060:     Input Parameters:
3061: +   x - vector
3062: -   xx_v - the Fortran90 pointer to the array

3064:     Output Parameter:
3065: .   ierr - error code

3067:     Example of Usage: 
3068: .vb
3069:     PetscScalar, pointer :: xx_v(:)
3070:     ....
3071:     call VecGetArrayF90(x,xx_v,ierr)
3072:     a = xx_v(3)
3073:     call VecRestoreArrayF90(x,xx_v,ierr)
3074: .ve
3075:    
3076:     Notes:
3077:     Not yet supported for all F90 compilers

3079:     Level: beginner

3081: .seealso:  VecGetArrayF90(), VecGetArray(), VecRestoreArray()

3083: M*/

3085: /*MC
3086:     VecDestroyVecsF90 - Frees a block of vectors obtained with VecDuplicateVecsF90().

3088:     Synopsis:
3089:     VecDestroyVecsF90({Vec, pointer :: x(:)},integer n,integer ierr)

3091:     Input Parameters:
3092: +   x - pointer to array of vector pointers
3093: -   n - the number of vectors previously obtained

3095:     Output Parameter:
3096: .   ierr - error code

3098:     Notes:
3099:     Not yet supported for all F90 compilers

3101:     Level: beginner

3103: .seealso:  VecDestroyVecs(), VecDuplicateVecsF90()

3105: M*/

3107: /*MC
3108:     VecGetArrayF90 - Accesses a vector array from Fortran90. For default PETSc
3109:     vectors, VecGetArrayF90() returns a pointer to the local data array. Otherwise,
3110:     this routine is implementation dependent. You MUST call VecRestoreArrayF90() 
3111:     when you no longer need access to the array.

3113:     Synopsis:
3114:     VecGetArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)

3116:     Not Collective 

3118:     Input Parameter:
3119: .   x - vector

3121:     Output Parameters:
3122: +   xx_v - the Fortran90 pointer to the array
3123: -   ierr - error code

3125:     Example of Usage: 
3126: .vb
3127:     PetscScalar, pointer :: xx_v(:)
3128:     ....
3129:     call VecGetArrayF90(x,xx_v,ierr)
3130:     a = xx_v(3)
3131:     call VecRestoreArrayF90(x,xx_v,ierr)
3132: .ve

3134:     Notes:
3135:     Not yet supported for all F90 compilers

3137:     Level: beginner

3139: .seealso:  VecRestoreArrayF90(), VecGetArray(), VecRestoreArray()

3141: M*/

3145: /*@C 
3146:   VecLoadIntoVector - Loads a vector that has been stored in binary format
3147:   with VecView().

3149:   Collective on PetscViewer 

3151:   Input Parameters:
3152: + viewer - binary file viewer, obtained from PetscViewerBinaryOpen()
3153: - vec - vector to contain files values (must be of correct length)

3155:   Level: intermediate

3157:   Notes:
3158:   The input file must contain the full global vector, as
3159:   written by the routine VecView().

3161:   Use VecLoad() to create the vector as the values are read in

3163:   Notes for advanced users:
3164:   Most users should not need to know the details of the binary storage
3165:   format, since VecLoad() and VecView() completely hide these details.
3166:   But for anyone who's interested, the standard binary matrix storage
3167:   format is
3168: .vb
3169:      int    VEC_FILE_COOKIE
3170:      int    number of rows
3171:      PetscScalar *values of all nonzeros
3172: .ve

3174:    Note for Cray users, the int's stored in the binary file are 32 bit
3175: integers; not 64 as they are represented in the memory, so if you
3176: write your own routines to read/write these binary files from the Cray
3177: you need to adjust the integer sizes that you read in, see
3178: PetscBinaryRead() and PetscBinaryWrite() to see how this may be
3179: done.

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

3187:    Concepts: vector^loading from file

3189: .seealso: PetscViewerBinaryOpen(), VecView(), MatLoad(), VecLoad() 
3190: @*/
3191: PetscErrorCode PETSCVEC_DLLEXPORT VecLoadIntoVector(PetscViewer viewer,Vec vec)
3192: {

3199:   if (!vec->ops->loadintovector) {
3200:     SETERRQ(PETSC_ERR_SUP,"Vector does not support load");
3201:   }
3202:   (*vec->ops->loadintovector)(viewer,vec);
3203:   PetscObjectStateIncrease((PetscObject)vec);
3204:   return(0);
3205: }

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

3212:    Collective on Vec

3214:    Input Parameter:
3215: .  v - the vector 

3217:    Output Parameter:
3218: .  v - the vector reciprocal

3220:    Level: intermediate

3222:    Concepts: vector^reciprocal

3224: @*/
3225: PetscErrorCode PETSCVEC_DLLEXPORT VecReciprocal(Vec vec)
3226: {

3232:   if (vec->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
3233:   if (!vec->ops->reciprocal) {
3234:     SETERRQ(PETSC_ERR_SUP,"Vector does not support reciprocal operation");
3235:   }
3236:   (*vec->ops->reciprocal)(vec);
3237:   PetscObjectStateIncrease((PetscObject)vec);
3238:   return(0);
3239: }

3243: PetscErrorCode PETSCVEC_DLLEXPORT VecSetOperation(Vec vec,VecOperation op, void (*f)(void))
3244: {
3247:   /* save the native version of the viewer */
3248:   if (op == VECOP_VIEW && !vec->ops->viewnative) {
3249:     vec->ops->viewnative = vec->ops->view;
3250:   }
3251:   (((void(**)(void))vec->ops)[(int)op]) = f;
3252:   return(0);
3253: }

3257: /*@
3258:    VecStashSetInitialSize - sets the sizes of the vec-stash, that is
3259:    used during the assembly process to store values that belong to 
3260:    other processors.

3262:    Collective on Vec

3264:    Input Parameters:
3265: +  vec   - the vector
3266: .  size  - the initial size of the stash.
3267: -  bsize - the initial size of the block-stash(if used).

3269:    Options Database Keys:
3270: +   -vecstash_initial_size <size> or <size0,size1,...sizep-1>
3271: -   -vecstash_block_initial_size <bsize> or <bsize0,bsize1,...bsizep-1>

3273:    Level: intermediate

3275:    Notes: 
3276:      The block-stash is used for values set with VecSetValuesBlocked() while
3277:      the stash is used for values set with VecSetValues()

3279:      Run with the option -log_info and look for output of the form
3280:      VecAssemblyBegin_MPIXXX:Stash has MM entries, uses nn mallocs.
3281:      to determine the appropriate value, MM, to use for size and 
3282:      VecAssemblyBegin_MPIXXX:Block-Stash has BMM entries, uses nn mallocs.
3283:      to determine the value, BMM to use for bsize

3285:    Concepts: vector^stash
3286:    Concepts: stash^vector

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

3290: @*/
3291: PetscErrorCode PETSCVEC_DLLEXPORT VecStashSetInitialSize(Vec vec,PetscInt size,PetscInt bsize)
3292: {

3297:   VecStashSetInitialSize_Private(&vec->stash,size);
3298:   VecStashSetInitialSize_Private(&vec->bstash,bsize);
3299:   return(0);
3300: }

3304: /*@
3305:    VecStashView - Prints the entries in the vector stash and block stash.

3307:    Collective on Vec

3309:    Input Parameters:
3310: +  vec   - the vector
3311: -  viewer - the viewer

3313:    Level: advanced

3315:    Concepts: vector^stash
3316:    Concepts: stash^vector

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

3320: @*/
3321: PetscErrorCode PETSCVEC_DLLEXPORT VecStashView(Vec v,PetscViewer viewer)
3322: {
3324:   PetscMPIInt    rank;
3325:   PetscInt       i,j;
3326:   PetscTruth     match;
3327:   VecStash       *s;
3328:   PetscScalar    val;


3335:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&match);
3336:   if (!match) SETERRQ1(PETSC_ERR_SUP,"Stash viewer only works with ASCII viewer not %s\n",((PetscObject)v)->type_name);
3337:   PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
3338:   MPI_Comm_rank(v->comm,&rank);
3339:   s = &v->bstash;

3341:   /* print block stash */
3342:   PetscViewerASCIISynchronizedPrintf(viewer,"[%d]Vector Block stash size %D block size %D\n",rank,s->n,s->bs);
3343:   for (i=0; i<s->n; i++) {
3344:     PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Element %D ",rank,s->idx[i]);
3345:     for (j=0; j<s->bs; j++) {
3346:       val = s->array[i*s->bs+j];
3347: #if defined(PETSC_USE_COMPLEX)
3348:       PetscViewerASCIISynchronizedPrintf(viewer,"(%18.16e %18.16e) ",PetscRealPart(val),PetscImaginaryPart(val));
3349: #else
3350:       PetscViewerASCIISynchronizedPrintf(viewer,"%18.16e ",val);
3351: #endif
3352:     }
3353:     PetscViewerASCIISynchronizedPrintf(viewer,"\n");
3354:   }
3355:   PetscViewerFlush(viewer);

3357:   s = &v->stash;

3359:   /* print basic stash */
3360:   PetscViewerASCIISynchronizedPrintf(viewer,"[%d]Vector stash size %D\n",rank,s->n);
3361:   for (i=0; i<s->n; i++) {
3362:     val = s->array[i];
3363: #if defined(PETSC_USE_COMPLEX)
3364:       PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Element %D (%18.16e %18.16e) ",rank,s->idx[i],PetscRealPart(val),PetscImaginaryPart(val));
3365: #else
3366:     PetscViewerASCIISynchronizedPrintf(viewer,"[%d] Element %D %18.16e\n",rank,s->idx[i],val);
3367: #endif
3368:   }
3369:   PetscViewerFlush(viewer);

3371:   PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
3372:   return(0);
3373: }

3377: /*@C
3378:    VecGetArray2d - Returns a pointer to a 2d contiguous array that contains this 
3379:    processor's portion of the vector data.  You MUST call VecRestoreArray2d() 
3380:    when you no longer need access to the array.

3382:    Not Collective

3384:    Input Parameter:
3385: +  x - the vector
3386: .  m - first dimension of two dimensional array
3387: .  n - second dimension of two dimensional array
3388: .  mstart - first index you will use in first coordinate direction (often 0)
3389: -  nstart - first index in the second coordinate direction (often 0)

3391:    Output Parameter:
3392: .  a - location to put pointer to the array

3394:    Level: beginner

3396:   Notes:
3397:    For a vector obtained from DACreateLocalVector() mstart and nstart are likely
3398:    obtained from the corner indices obtained from DAGetGhostCorners() while for
3399:    DACreateGlobalVector() they are the corner indices from DAGetCorners(). In both cases
3400:    the arguments from DAGet[Ghost}Corners() are reversed in the call to VecGetArray2d().
3401:    
3402:    For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.

3404:    Concepts: vector^accessing local values as 2d array

3406: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3407:           VecRestoreArray2d(), DAVecGetarray(), DAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3408:           VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3409: @*/
3410: PetscErrorCode PETSCVEC_DLLEXPORT VecGetArray2d(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3411: {
3413:   PetscInt       i,N;
3414:   PetscScalar    *aa;

3420:   VecGetLocalSize(x,&N);
3421:   if (m*n != N) SETERRQ3(PETSC_ERR_ARG_INCOMP,"Local array size %D does not match 2d array dimensions %D by %D",N,m,n);
3422:   VecGetArray(x,&aa);

3424:   PetscMalloc(m*sizeof(PetscScalar*),a);
3425:   for (i=0; i<m; i++) (*a)[i] = aa + i*n - nstart;
3426:   *a -= mstart;
3427:   return(0);
3428: }

3432: /*@C
3433:    VecRestoreArray2d - Restores a vector after VecGetArray2d() has been called.

3435:    Not Collective

3437:    Input Parameters:
3438: +  x - the vector
3439: .  m - first dimension of two dimensional array
3440: .  n - second dimension of the two dimensional array
3441: .  mstart - first index you will use in first coordinate direction (often 0)
3442: .  nstart - first index in the second coordinate direction (often 0)
3443: -  a - location of pointer to array obtained from VecGetArray2d()

3445:    Level: beginner

3447:    Notes:
3448:    For regular PETSc vectors this routine does not involve any copies. For
3449:    any special vectors that do not store local vector data in a contiguous
3450:    array, this routine will copy the data back into the underlying 
3451:    vector data structure from the array obtained with VecGetArray().

3453:    This routine actually zeros out the a pointer. 

3455: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3456:           VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DAVecGetArray(), DAVecRestoreArray()
3457:           VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3458: @*/
3459: PetscErrorCode PETSCVEC_DLLEXPORT VecRestoreArray2d(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3460: {
3462:   void           *dummy;

3468:   dummy = (void*)(*a + mstart);
3469:   PetscFree(dummy);
3470:   VecRestoreArray(x,PETSC_NULL);
3471:   return(0);
3472: }

3476: /*@C
3477:    VecGetArray1d - Returns a pointer to a 1d contiguous array that contains this 
3478:    processor's portion of the vector data.  You MUST call VecRestoreArray1d() 
3479:    when you no longer need access to the array.

3481:    Not Collective

3483:    Input Parameter:
3484: +  x - the vector
3485: .  m - first dimension of two dimensional array
3486: -  mstart - first index you will use in first coordinate direction (often 0)

3488:    Output Parameter:
3489: .  a - location to put pointer to the array

3491:    Level: beginner

3493:   Notes:
3494:    For a vector obtained from DACreateLocalVector() mstart are likely
3495:    obtained from the corner indices obtained from DAGetGhostCorners() while for
3496:    DACreateGlobalVector() they are the corner indices from DAGetCorners(). 
3497:    
3498:    For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.

3500: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3501:           VecRestoreArray2d(), DAVecGetArray(), DAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3502:           VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3503: @*/
3504: PetscErrorCode PETSCVEC_DLLEXPORT VecGetArray1d(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
3505: {
3507:   PetscInt       N;

3513:   VecGetLocalSize(x,&N);
3514:   if (m != N) SETERRQ2(PETSC_ERR_ARG_OUTOFRANGE,"Local array size %D does not match 1d array dimensions %D",N,m);
3515:   VecGetArray(x,a);
3516:   *a  -= mstart;
3517:   return(0);
3518: }

3522: /*@C
3523:    VecRestoreArray1d - Restores a vector after VecGetArray1d() has been called.

3525:    Not Collective

3527:    Input Parameters:
3528: +  x - the vector
3529: .  m - first dimension of two dimensional array
3530: .  mstart - first index you will use in first coordinate direction (often 0)
3531: -  a - location of pointer to array obtained from VecGetArray21()

3533:    Level: beginner

3535:    Notes:
3536:    For regular PETSc vectors this routine does not involve any copies. For
3537:    any special vectors that do not store local vector data in a contiguous
3538:    array, this routine will copy the data back into the underlying 
3539:    vector data structure from the array obtained with VecGetArray1d().

3541:    This routine actually zeros out the a pointer. 

3543:    Concepts: vector^accessing local values as 1d array

3545: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3546:           VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DAVecGetArray(), DAVecRestoreArray()
3547:           VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
3548: @*/
3549: PetscErrorCode PETSCVEC_DLLEXPORT VecRestoreArray1d(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
3550: {

3556:   VecRestoreArray(x,PETSC_NULL);
3557:   return(0);
3558: }

3562: /*@
3563:    VecConjugate - Conjugates a vector.

3565:    Collective on Vec

3567:    Input Parameters:
3568: .  x - the vector

3570:    Level: intermediate

3572:    Concepts: vector^conjugate

3574: @*/
3575: PetscErrorCode PETSCVEC_DLLEXPORT VecConjugate(Vec x)
3576: {
3577: #ifdef PETSC_USE_COMPLEX

3583:   if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
3584:   (*x->ops->conjugate)(x);
3585:   /* we need to copy norms here */
3586:   PetscObjectStateIncrease((PetscObject)x);
3587:   return(0);
3588: #else
3589:   return(0);
3590: #endif
3591: }

3595: /*@C
3596:    VecGetArray3d - Returns a pointer to a 3d contiguous array that contains this 
3597:    processor's portion of the vector data.  You MUST call VecRestoreArray3d() 
3598:    when you no longer need access to the array.

3600:    Not Collective

3602:    Input Parameter:
3603: +  x - the vector
3604: .  m - first dimension of three dimensional array
3605: .  n - second dimension of three dimensional array
3606: .  p - third dimension of three dimensional array
3607: .  mstart - first index you will use in first coordinate direction (often 0)
3608: .  nstart - first index in the second coordinate direction (often 0)
3609: -  pstart - first index in the third coordinate direction (often 0)

3611:    Output Parameter:
3612: .  a - location to put pointer to the array

3614:    Level: beginner

3616:   Notes:
3617:    For a vector obtained from DACreateLocalVector() mstart, nstart, and pstart are likely
3618:    obtained from the corner indices obtained from DAGetGhostCorners() while for
3619:    DACreateGlobalVector() they are the corner indices from DAGetCorners(). In both cases
3620:    the arguments from DAGet[Ghost}Corners() are reversed in the call to VecGetArray3d().
3621:    
3622:    For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.

3624:    Concepts: vector^accessing local values as 3d array

3626: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3627:           VecRestoreArray2d(), DAVecGetarray(), DAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3628:           VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3629: @*/
3630: PetscErrorCode PETSCVEC_DLLEXPORT VecGetArray3d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3631: {
3633:   PetscInt       i,N,j;
3634:   PetscScalar    *aa,**b;

3640:   VecGetLocalSize(x,&N);
3641:   if (m*n*p != N) SETERRQ4(PETSC_ERR_ARG_INCOMP,"Local array size %D does not match 3d array dimensions %D by %D by %D",N,m,n,p);
3642:   VecGetArray(x,&aa);

3644:   PetscMalloc(m*sizeof(PetscScalar**)+m*n*sizeof(PetscScalar*),a);
3645:   b    = (PetscScalar **)((*a) + m);
3646:   for (i=0; i<m; i++)   (*a)[i] = b + i*n - nstart;
3647:   for (i=0; i<m; i++) {
3648:     for (j=0; j<n; j++) {
3649:       b[i*n+j] = aa + i*n*p + j*p - pstart;
3650:     }
3651:   }
3652:   *a -= mstart;
3653:   return(0);
3654: }

3658: /*@C
3659:    VecRestoreArray3d - Restores a vector after VecGetArray3d() has been called.

3661:    Not Collective

3663:    Input Parameters:
3664: +  x - the vector
3665: .  m - first dimension of three dimensional array
3666: .  n - second dimension of the three dimensional array
3667: .  p - third dimension of the three dimensional array
3668: .  mstart - first index you will use in first coordinate direction (often 0)
3669: .  nstart - first index in the second coordinate direction (often 0)
3670: .  pstart - first index in the third coordinate direction (often 0)
3671: -  a - location of pointer to array obtained from VecGetArray3d()

3673:    Level: beginner

3675:    Notes:
3676:    For regular PETSc vectors this routine does not involve any copies. For
3677:    any special vectors that do not store local vector data in a contiguous
3678:    array, this routine will copy the data back into the underlying 
3679:    vector data structure from the array obtained with VecGetArray().

3681:    This routine actually zeros out the a pointer. 

3683: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3684:           VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DAVecGetArray(), DAVecRestoreArray()
3685:           VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3686: @*/
3687: PetscErrorCode PETSCVEC_DLLEXPORT VecRestoreArray3d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3688: {
3690:   void           *dummy;

3696:   dummy = (void*)(*a + mstart);
3697:   PetscFree(dummy);
3698:   VecRestoreArray(x,PETSC_NULL);
3699:   return(0);
3700: }

3704: /*@C
3705:    VecGetArray4d - Returns a pointer to a 4d contiguous array that contains this 
3706:    processor's portion of the vector data.  You MUST call VecRestoreArray4d() 
3707:    when you no longer need access to the array.

3709:    Not Collective

3711:    Input Parameter:
3712: +  x - the vector
3713: .  m - first dimension of four dimensional array
3714: .  n - second dimension of four dimensional array
3715: .  p - third dimension of four dimensional array
3716: .  q - fourth dimension of four dimensional array
3717: .  mstart - first index you will use in first coordinate direction (often 0)
3718: .  nstart - first index in the second coordinate direction (often 0)
3719: .  pstart - first index in the third coordinate direction (often 0)
3720: -  qstart - first index in the fourth coordinate direction (often 0)

3722:    Output Parameter:
3723: .  a - location to put pointer to the array

3725:    Level: beginner

3727:   Notes:
3728:    For a vector obtained from DACreateLocalVector() mstart, nstart, and pstart are likely
3729:    obtained from the corner indices obtained from DAGetGhostCorners() while for
3730:    DACreateGlobalVector() they are the corner indices from DAGetCorners(). In both cases
3731:    the arguments from DAGet[Ghost}Corners() are reversed in the call to VecGetArray3d().
3732:    
3733:    For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.

3735:    Concepts: vector^accessing local values as 3d array

3737: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3738:           VecRestoreArray2d(), DAVecGetarray(), DAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3739:           VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3740: @*/
3741: PetscErrorCode PETSCVEC_DLLEXPORT VecGetArray4d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3742: {
3744:   PetscInt       i,N,j,k;
3745:   PetscScalar    *aa,***b,**c;

3751:   VecGetLocalSize(x,&N);
3752:   if (m*n*p*q != N) SETERRQ5(PETSC_ERR_ARG_INCOMP,"Local array size %D does not match 4d array dimensions %D by %D by %D by %D",N,m,n,p,q);
3753:   VecGetArray(x,&aa);

3755:   PetscMalloc(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p*sizeof(PetscScalar*),a);
3756:   b    = (PetscScalar ***)((*a) + m);
3757:   c    = (PetscScalar **)(b + m*n);
3758:   for (i=0; i<m; i++)   (*a)[i] = b + i*n - nstart;
3759:   for (i=0; i<m; i++) {
3760:     for (j=0; j<n; j++) {
3761:       b[i*n+j] = c + i*n*p + j*p - pstart;
3762:     }
3763:   }
3764:   for (i=0; i<m; i++) {
3765:     for (j=0; j<n; j++) {
3766:       for (k=0; k<p; k++) {
3767:         c[i*n*p+j*p+k] = aa + i*n*p*q + j*p*q + k*q - qstart;
3768:       }
3769:     }
3770:   }
3771:   *a -= mstart;
3772:   return(0);
3773: }

3777: /*@C
3778:    VecRestoreArray4d - Restores a vector after VecGetArray3d() has been called.

3780:    Not Collective

3782:    Input Parameters:
3783: +  x - the vector
3784: .  m - first dimension of four dimensional array
3785: .  n - second dimension of the four dimensional array
3786: .  p - third dimension of the four dimensional array
3787: .  q - fourth dimension of the four dimensional array
3788: .  mstart - first index you will use in first coordinate direction (often 0)
3789: .  nstart - first index in the second coordinate direction (often 0)
3790: .  pstart - first index in the third coordinate direction (often 0)
3791: .  qstart - first index in the fourth coordinate direction (often 0)
3792: -  a - location of pointer to array obtained from VecGetArray4d()

3794:    Level: beginner

3796:    Notes:
3797:    For regular PETSc vectors this routine does not involve any copies. For
3798:    any special vectors that do not store local vector data in a contiguous
3799:    array, this routine will copy the data back into the underlying 
3800:    vector data structure from the array obtained with VecGetArray().

3802:    This routine actually zeros out the a pointer. 

3804: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3805:           VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DAVecGetArray(), DAVecRestoreArray()
3806:           VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3807: @*/
3808: PetscErrorCode PETSCVEC_DLLEXPORT VecRestoreArray4d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3809: {
3811:   void           *dummy;

3817:   dummy = (void*)(*a + mstart);
3818:   PetscFree(dummy);
3819:   VecRestoreArray(x,PETSC_NULL);
3820:   return(0);
3821: }

3823: EXTERN PetscErrorCode VecStashGetInfo_Private(VecStash*,PetscInt*,PetscInt*);
3826: /*@ 
3827:    VecStashGetInfo - Gets how many values are currently in the vector stash, i.e. need
3828:        to be communicated to other processors during the VecAssemblyBegin/End() process

3830:     Not collective

3832:    Input Parameter:
3833: .   vec - the vector

3835:    Output Parameters:
3836: +   nstash   - the size of the stash
3837: .   reallocs - the number of additional mallocs incurred.
3838: .   bnstash   - the size of the block stash
3839: -   breallocs - the number of additional mallocs incurred.in the block stash
3840:  
3841:    Level: advanced

3843: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), Vec, VecStashSetInitialSize(), VecStashView()
3844:   
3845: @*/
3846: PetscErrorCode PETSCVEC_DLLEXPORT VecStashGetInfo(Vec vec,PetscInt *nstash,PetscInt *reallocs,PetscInt *bnstash,PetscInt *brealloc)
3847: {
3850:   VecStashGetInfo_Private(&vec->stash,nstash,reallocs);
3851:   VecStashGetInfo_Private(&vec->bstash,nstash,reallocs);
3852:   return(0);
3853: }