Actual source code: pbvec.c
petsc-3.6.4 2016-04-12
2: /*
3: This file contains routines for Parallel vector operations.
4: */
5: #include <../src/vec/vec/impls/mpi/pvecimpl.h> /*I "petscvec.h" I*/
9: static PetscErrorCode VecDot_MPI(Vec xin,Vec yin,PetscScalar *z)
10: {
11: PetscScalar sum,work;
15: VecDot_Seq(xin,yin,&work);
16: MPI_Allreduce(&work,&sum,1,MPIU_SCALAR,MPIU_SUM,PetscObjectComm((PetscObject)xin));
17: *z = sum;
18: return(0);
19: }
23: static PetscErrorCode VecTDot_MPI(Vec xin,Vec yin,PetscScalar *z)
24: {
25: PetscScalar sum,work;
29: VecTDot_Seq(xin,yin,&work);
30: MPI_Allreduce(&work,&sum,1,MPIU_SCALAR,MPIU_SUM,PetscObjectComm((PetscObject)xin));
31: *z = sum;
32: return(0);
33: }
35: extern PetscErrorCode VecView_MPI_Draw(Vec,PetscViewer);
39: static PetscErrorCode VecPlaceArray_MPI(Vec vin,const PetscScalar *a)
40: {
42: Vec_MPI *v = (Vec_MPI*)vin->data;
45: if (v->unplacedarray) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"VecPlaceArray() was already called on this vector, without a call to VecResetArray()");
46: v->unplacedarray = v->array; /* save previous array so reset can bring it back */
47: v->array = (PetscScalar*)a;
48: if (v->localrep) {
49: VecPlaceArray(v->localrep,a);
50: }
51: return(0);
52: }
54: extern PetscErrorCode VecGetValues_MPI(Vec,PetscInt,const PetscInt[],PetscScalar[]);
58: static PetscErrorCode VecDuplicate_MPI(Vec win,Vec *v)
59: {
61: Vec_MPI *vw,*w = (Vec_MPI*)win->data;
62: PetscScalar *array;
65: VecCreate(PetscObjectComm((PetscObject)win),v);
66: PetscLayoutReference(win->map,&(*v)->map);
68: VecCreate_MPI_Private(*v,PETSC_TRUE,w->nghost,0);
69: vw = (Vec_MPI*)(*v)->data;
70: PetscMemcpy((*v)->ops,win->ops,sizeof(struct _VecOps));
72: /* save local representation of the parallel vector (and scatter) if it exists */
73: if (w->localrep) {
74: VecGetArray(*v,&array);
75: VecCreateSeqWithArray(PETSC_COMM_SELF,PetscAbs(win->map->bs),win->map->n+w->nghost,array,&vw->localrep);
76: PetscMemcpy(vw->localrep->ops,w->localrep->ops,sizeof(struct _VecOps));
77: VecRestoreArray(*v,&array);
78: PetscLogObjectParent((PetscObject)*v,(PetscObject)vw->localrep);
80: vw->localupdate = w->localupdate;
81: if (vw->localupdate) {
82: PetscObjectReference((PetscObject)vw->localupdate);
83: }
84: }
86: /* New vector should inherit stashing property of parent */
87: (*v)->stash.donotstash = win->stash.donotstash;
88: (*v)->stash.ignorenegidx = win->stash.ignorenegidx;
90: PetscObjectListDuplicate(((PetscObject)win)->olist,&((PetscObject)(*v))->olist);
91: PetscFunctionListDuplicate(((PetscObject)win)->qlist,&((PetscObject)(*v))->qlist);
93: (*v)->map->bs = PetscAbs(win->map->bs);
94: (*v)->bstash.bs = win->bstash.bs;
95: return(0);
96: }
98: extern PetscErrorCode VecSetOption_MPI(Vec,VecOption,PetscBool);
99: extern PetscErrorCode VecResetArray_MPI(Vec);
101: static struct _VecOps DvOps = { VecDuplicate_MPI, /* 1 */
102: VecDuplicateVecs_Default,
103: VecDestroyVecs_Default,
104: VecDot_MPI,
105: VecMDot_MPI,
106: VecNorm_MPI,
107: VecTDot_MPI,
108: VecMTDot_MPI,
109: VecScale_Seq,
110: VecCopy_Seq, /* 10 */
111: VecSet_Seq,
112: VecSwap_Seq,
113: VecAXPY_Seq,
114: VecAXPBY_Seq,
115: VecMAXPY_Seq,
116: VecAYPX_Seq,
117: VecWAXPY_Seq,
118: VecAXPBYPCZ_Seq,
119: VecPointwiseMult_Seq,
120: VecPointwiseDivide_Seq,
121: VecSetValues_MPI, /* 20 */
122: VecAssemblyBegin_MPI,
123: VecAssemblyEnd_MPI,
124: 0,
125: VecGetSize_MPI,
126: VecGetSize_Seq,
127: 0,
128: VecMax_MPI,
129: VecMin_MPI,
130: VecSetRandom_Seq,
131: VecSetOption_MPI,
132: VecSetValuesBlocked_MPI,
133: VecDestroy_MPI,
134: VecView_MPI,
135: VecPlaceArray_MPI,
136: VecReplaceArray_Seq,
137: VecDot_Seq,
138: VecTDot_Seq,
139: VecNorm_Seq,
140: VecMDot_Seq,
141: VecMTDot_Seq,
142: VecLoad_Default,
143: VecReciprocal_Default,
144: VecConjugate_Seq,
145: 0,
146: 0,
147: VecResetArray_MPI,
148: 0,
149: VecMaxPointwiseDivide_Seq,
150: VecPointwiseMax_Seq,
151: VecPointwiseMaxAbs_Seq,
152: VecPointwiseMin_Seq,
153: VecGetValues_MPI,
154: 0,
155: 0,
156: 0,
157: 0,
158: 0,
159: 0,
160: VecStrideGather_Default,
161: VecStrideScatter_Default,
162: 0,
163: 0,
164: 0,
165: 0,
166: 0,
167: VecStrideSubSetGather_Default,
168: VecStrideSubSetScatter_Default,
169: 0,
170: 0
171: };
175: /*
176: VecCreate_MPI_Private - Basic create routine called by VecCreate_MPI() (i.e. VecCreateMPI()),
177: VecCreateMPIWithArray(), VecCreate_Shared() (i.e. VecCreateShared()), VecCreateGhost(),
178: VecDuplicate_MPI(), VecCreateGhostWithArray(), VecDuplicate_MPI(), and VecDuplicate_Shared()
180: If alloc is true and array is NULL then this routine allocates the space, otherwise
181: no space is allocated.
182: */
183: PetscErrorCode VecCreate_MPI_Private(Vec v,PetscBool alloc,PetscInt nghost,const PetscScalar array[])
184: {
185: Vec_MPI *s;
189: PetscNewLog(v,&s);
190: v->data = (void*)s;
191: PetscMemcpy(v->ops,&DvOps,sizeof(DvOps));
192: s->nghost = nghost;
193: v->petscnative = PETSC_TRUE;
195: PetscLayoutSetUp(v->map);
197: s->array = (PetscScalar*)array;
198: s->array_allocated = 0;
199: if (alloc && !array) {
200: PetscInt n = v->map->n+nghost;
201: PetscMalloc1(n,&s->array);
202: PetscLogObjectMemory((PetscObject)v,n*sizeof(PetscScalar));
203: PetscMemzero(s->array,n*sizeof(PetscScalar));
204: s->array_allocated = s->array;
205: }
207: /* By default parallel vectors do not have local representation */
208: s->localrep = 0;
209: s->localupdate = 0;
211: v->stash.insertmode = NOT_SET_VALUES;
212: /* create the stashes. The block-size for bstash is set later when
213: VecSetValuesBlocked is called.
214: */
215: VecStashCreate_Private(PetscObjectComm((PetscObject)v),1,&v->stash);
216: VecStashCreate_Private(PetscObjectComm((PetscObject)v),PetscAbs(v->map->bs),&v->bstash);
218: #if defined(PETSC_HAVE_MATLAB_ENGINE)
219: PetscObjectComposeFunction((PetscObject)v,"PetscMatlabEnginePut_C",VecMatlabEnginePut_Default);
220: PetscObjectComposeFunction((PetscObject)v,"PetscMatlabEngineGet_C",VecMatlabEngineGet_Default);
221: #endif
222: PetscObjectChangeTypeName((PetscObject)v,VECMPI);
223: return(0);
224: }
226: /*MC
227: VECMPI - VECMPI = "mpi" - The basic parallel vector
229: Options Database Keys:
230: . -vec_type mpi - sets the vector type to VECMPI during a call to VecSetFromOptions()
232: Level: beginner
234: .seealso: VecCreate(), VecSetType(), VecSetFromOptions(), VecCreateMpiWithArray(), VECMPI, VecType, VecCreateMPI(), VecCreateMpi()
235: M*/
239: PETSC_EXTERN PetscErrorCode VecCreate_MPI(Vec vv)
240: {
244: VecCreate_MPI_Private(vv,PETSC_TRUE,0,0);
245: return(0);
246: }
248: /*MC
249: VECSTANDARD = "standard" - A VECSEQ on one process and VECMPI on more than one process
251: Options Database Keys:
252: . -vec_type standard - sets a vector type to standard on calls to VecSetFromOptions()
254: Level: beginner
256: .seealso: VecCreateSeq(), VecCreateMPI()
257: M*/
261: PETSC_EXTERN PetscErrorCode VecCreate_Standard(Vec v)
262: {
264: PetscMPIInt size;
267: MPI_Comm_size(PetscObjectComm((PetscObject)v),&size);
268: if (size == 1) {
269: VecSetType(v,VECSEQ);
270: } else {
271: VecSetType(v,VECMPI);
272: }
273: return(0);
274: }
278: /*@C
279: VecCreateMPIWithArray - Creates a parallel, array-style vector,
280: where the user provides the array space to store the vector values.
282: Collective on MPI_Comm
284: Input Parameters:
285: + comm - the MPI communicator to use
286: . bs - block size, same meaning as VecSetBlockSize()
287: . n - local vector length, cannot be PETSC_DECIDE
288: . N - global vector length (or PETSC_DECIDE to have calculated)
289: - array - the user provided array to store the vector values
291: Output Parameter:
292: . vv - the vector
294: Notes:
295: Use VecDuplicate() or VecDuplicateVecs() to form additional vectors of the
296: same type as an existing vector.
298: If the user-provided array is NULL, then VecPlaceArray() can be used
299: at a later stage to SET the array for storing the vector values.
301: PETSc does NOT free the array when the vector is destroyed via VecDestroy().
302: The user should not free the array until the vector is destroyed.
304: Level: intermediate
306: Concepts: vectors^creating with array
308: .seealso: VecCreateSeqWithArray(), VecCreate(), VecDuplicate(), VecDuplicateVecs(), VecCreateGhost(),
309: VecCreateMPI(), VecCreateGhostWithArray(), VecPlaceArray()
311: @*/
312: PetscErrorCode VecCreateMPIWithArray(MPI_Comm comm,PetscInt bs,PetscInt n,PetscInt N,const PetscScalar array[],Vec *vv)
313: {
317: if (n == PETSC_DECIDE) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must set local size of vector");
318: PetscSplitOwnership(comm,&n,&N);
319: VecCreate(comm,vv);
320: VecSetSizes(*vv,n,N);
321: VecSetBlockSize(*vv,bs);
322: VecCreate_MPI_Private(*vv,PETSC_FALSE,0,array);
323: return(0);
324: }
328: /*@C
329: VecCreateGhostWithArray - Creates a parallel vector with ghost padding on each processor;
330: the caller allocates the array space.
332: Collective on MPI_Comm
334: Input Parameters:
335: + comm - the MPI communicator to use
336: . n - local vector length
337: . N - global vector length (or PETSC_DECIDE to have calculated if n is given)
338: . nghost - number of local ghost points
339: . ghosts - global indices of ghost points (or NULL if not needed), these do not need to be in increasing order (sorted)
340: - array - the space to store the vector values (as long as n + nghost)
342: Output Parameter:
343: . vv - the global vector representation (without ghost points as part of vector)
345: Notes:
346: Use VecGhostGetLocalForm() to access the local, ghosted representation
347: of the vector.
349: This also automatically sets the ISLocalToGlobalMapping() for this vector.
351: Level: advanced
353: Concepts: vectors^creating with array
354: Concepts: vectors^ghosted
356: .seealso: VecCreate(), VecGhostGetLocalForm(), VecGhostRestoreLocalForm(),
357: VecCreateGhost(), VecCreateSeqWithArray(), VecCreateMPIWithArray(),
358: VecCreateGhostBlock(), VecCreateGhostBlockWithArray(), VecMPISetGhost()
360: @*/
361: PetscErrorCode VecCreateGhostWithArray(MPI_Comm comm,PetscInt n,PetscInt N,PetscInt nghost,const PetscInt ghosts[],const PetscScalar array[],Vec *vv)
362: {
363: PetscErrorCode ierr;
364: Vec_MPI *w;
365: PetscScalar *larray;
366: IS from,to;
367: ISLocalToGlobalMapping ltog;
368: PetscInt rstart,i,*indices;
371: *vv = 0;
373: if (n == PETSC_DECIDE) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must set local size");
374: if (nghost == PETSC_DECIDE) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must set local ghost size");
375: if (nghost < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Ghost length must be >= 0");
376: PetscSplitOwnership(comm,&n,&N);
377: /* Create global representation */
378: VecCreate(comm,vv);
379: VecSetSizes(*vv,n,N);
380: VecCreate_MPI_Private(*vv,PETSC_TRUE,nghost,array);
381: w = (Vec_MPI*)(*vv)->data;
382: /* Create local representation */
383: VecGetArray(*vv,&larray);
384: VecCreateSeqWithArray(PETSC_COMM_SELF,1,n+nghost,larray,&w->localrep);
385: PetscLogObjectParent((PetscObject)*vv,(PetscObject)w->localrep);
386: VecRestoreArray(*vv,&larray);
388: /*
389: Create scatter context for scattering (updating) ghost values
390: */
391: ISCreateGeneral(comm,nghost,ghosts,PETSC_COPY_VALUES,&from);
392: ISCreateStride(PETSC_COMM_SELF,nghost,n,1,&to);
393: VecScatterCreate(*vv,from,w->localrep,to,&w->localupdate);
394: PetscLogObjectParent((PetscObject)*vv,(PetscObject)w->localupdate);
395: ISDestroy(&to);
396: ISDestroy(&from);
398: /* set local to global mapping for ghosted vector */
399: PetscMalloc1(n+nghost,&indices);
400: VecGetOwnershipRange(*vv,&rstart,NULL);
401: for (i=0; i<n; i++) {
402: indices[i] = rstart + i;
403: }
404: for (i=0; i<nghost; i++) {
405: indices[n+i] = ghosts[i];
406: }
407: ISLocalToGlobalMappingCreate(comm,1,n+nghost,indices,PETSC_OWN_POINTER,<og);
408: VecSetLocalToGlobalMapping(*vv,ltog);
409: ISLocalToGlobalMappingDestroy(<og);
410: return(0);
411: }
415: /*@
416: VecCreateGhost - Creates a parallel vector with ghost padding on each processor.
418: Collective on MPI_Comm
420: Input Parameters:
421: + comm - the MPI communicator to use
422: . n - local vector length
423: . N - global vector length (or PETSC_DECIDE to have calculated if n is given)
424: . nghost - number of local ghost points
425: - ghosts - global indices of ghost points, these do not need to be in increasing order (sorted)
427: Output Parameter:
428: . vv - the global vector representation (without ghost points as part of vector)
430: Notes:
431: Use VecGhostGetLocalForm() to access the local, ghosted representation
432: of the vector.
434: This also automatically sets the ISLocalToGlobalMapping() for this vector.
436: Level: advanced
438: Concepts: vectors^ghosted
440: .seealso: VecCreateSeq(), VecCreate(), VecDuplicate(), VecDuplicateVecs(), VecCreateMPI(),
441: VecGhostGetLocalForm(), VecGhostRestoreLocalForm(), VecGhostUpdateBegin(),
442: VecCreateGhostWithArray(), VecCreateMPIWithArray(), VecGhostUpdateEnd(),
443: VecCreateGhostBlock(), VecCreateGhostBlockWithArray(), VecMPISetGhost()
445: @*/
446: PetscErrorCode VecCreateGhost(MPI_Comm comm,PetscInt n,PetscInt N,PetscInt nghost,const PetscInt ghosts[],Vec *vv)
447: {
451: VecCreateGhostWithArray(comm,n,N,nghost,ghosts,0,vv);
452: return(0);
453: }
457: /*@
458: VecMPISetGhost - Sets the ghost points for an MPI ghost vector
460: Collective on Vec
462: Input Parameters:
463: + vv - the MPI vector
464: . nghost - number of local ghost points
465: - ghosts - global indices of ghost points, these do not need to be in increasing order (sorted)
468: Notes:
469: Use VecGhostGetLocalForm() to access the local, ghosted representation
470: of the vector.
472: This also automatically sets the ISLocalToGlobalMapping() for this vector.
474: You must call this AFTER you have set the type of the vector (with VecSetType()) and the size (with VecSetSizes()).
476: Level: advanced
478: Concepts: vectors^ghosted
480: .seealso: VecCreateSeq(), VecCreate(), VecDuplicate(), VecDuplicateVecs(), VecCreateMPI(),
481: VecGhostGetLocalForm(), VecGhostRestoreLocalForm(), VecGhostUpdateBegin(),
482: VecCreateGhostWithArray(), VecCreateMPIWithArray(), VecGhostUpdateEnd(),
483: VecCreateGhostBlock(), VecCreateGhostBlockWithArray()
485: @*/
486: PetscErrorCode VecMPISetGhost(Vec vv,PetscInt nghost,const PetscInt ghosts[])
487: {
489: PetscBool flg;
492: PetscObjectTypeCompare((PetscObject)vv,VECMPI,&flg);
493: /* if already fully existant VECMPI then basically destroy it and rebuild with ghosting */
494: if (flg) {
495: PetscInt n,N;
496: Vec_MPI *w;
497: PetscScalar *larray;
498: IS from,to;
499: ISLocalToGlobalMapping ltog;
500: PetscInt rstart,i,*indices;
501: MPI_Comm comm;
503: PetscObjectGetComm((PetscObject)vv,&comm);
504: n = vv->map->n;
505: N = vv->map->N;
506: (*vv->ops->destroy)(vv);
507: VecSetSizes(vv,n,N);
508: VecCreate_MPI_Private(vv,PETSC_TRUE,nghost,NULL);
509: w = (Vec_MPI*)(vv)->data;
510: /* Create local representation */
511: VecGetArray(vv,&larray);
512: VecCreateSeqWithArray(PETSC_COMM_SELF,1,n+nghost,larray,&w->localrep);
513: PetscLogObjectParent((PetscObject)vv,(PetscObject)w->localrep);
514: VecRestoreArray(vv,&larray);
516: /*
517: Create scatter context for scattering (updating) ghost values
518: */
519: ISCreateGeneral(comm,nghost,ghosts,PETSC_COPY_VALUES,&from);
520: ISCreateStride(PETSC_COMM_SELF,nghost,n,1,&to);
521: VecScatterCreate(vv,from,w->localrep,to,&w->localupdate);
522: PetscLogObjectParent((PetscObject)vv,(PetscObject)w->localupdate);
523: ISDestroy(&to);
524: ISDestroy(&from);
526: /* set local to global mapping for ghosted vector */
527: PetscMalloc1(n+nghost,&indices);
528: VecGetOwnershipRange(vv,&rstart,NULL);
530: for (i=0; i<n; i++) indices[i] = rstart + i;
531: for (i=0; i<nghost; i++) indices[n+i] = ghosts[i];
533: ISLocalToGlobalMappingCreate(comm,1,n+nghost,indices,PETSC_OWN_POINTER,<og);
534: VecSetLocalToGlobalMapping(vv,ltog);
535: ISLocalToGlobalMappingDestroy(<og);
536: } else if (vv->ops->create == VecCreate_MPI) SETERRQ(PetscObjectComm((PetscObject)vv),PETSC_ERR_ARG_WRONGSTATE,"Must set local or global size before setting ghosting");
537: else if (!((PetscObject)vv)->type_name) SETERRQ(PetscObjectComm((PetscObject)vv),PETSC_ERR_ARG_WRONGSTATE,"Must set type to VECMPI before ghosting");
538: return(0);
539: }
542: /* ------------------------------------------------------------------------------------------*/
545: /*@C
546: VecCreateGhostBlockWithArray - Creates a parallel vector with ghost padding on each processor;
547: the caller allocates the array space. Indices in the ghost region are based on blocks.
549: Collective on MPI_Comm
551: Input Parameters:
552: + comm - the MPI communicator to use
553: . bs - block size
554: . n - local vector length
555: . N - global vector length (or PETSC_DECIDE to have calculated if n is given)
556: . nghost - number of local ghost blocks
557: . ghosts - global indices of ghost blocks (or NULL if not needed), counts are by block not by index, these do not need to be in increasing order (sorted)
558: - array - the space to store the vector values (as long as n + nghost*bs)
560: Output Parameter:
561: . vv - the global vector representation (without ghost points as part of vector)
563: Notes:
564: Use VecGhostGetLocalForm() to access the local, ghosted representation
565: of the vector.
567: n is the local vector size (total local size not the number of blocks) while nghost
568: is the number of blocks in the ghost portion, i.e. the number of elements in the ghost
569: portion is bs*nghost
571: Level: advanced
573: Concepts: vectors^creating ghosted
574: Concepts: vectors^creating with array
576: .seealso: VecCreate(), VecGhostGetLocalForm(), VecGhostRestoreLocalForm(),
577: VecCreateGhost(), VecCreateSeqWithArray(), VecCreateMPIWithArray(),
578: VecCreateGhostWithArray(), VecCreateGhostBlock()
580: @*/
581: PetscErrorCode VecCreateGhostBlockWithArray(MPI_Comm comm,PetscInt bs,PetscInt n,PetscInt N,PetscInt nghost,const PetscInt ghosts[],const PetscScalar array[],Vec *vv)
582: {
583: PetscErrorCode ierr;
584: Vec_MPI *w;
585: PetscScalar *larray;
586: IS from,to;
587: ISLocalToGlobalMapping ltog;
588: PetscInt rstart,i,nb,*indices;
591: *vv = 0;
593: if (n == PETSC_DECIDE) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must set local size");
594: if (nghost == PETSC_DECIDE) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must set local ghost size");
595: if (nghost < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Ghost length must be >= 0");
596: if (n % bs) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Local size must be a multiple of block size");
597: PetscSplitOwnership(comm,&n,&N);
598: /* Create global representation */
599: VecCreate(comm,vv);
600: VecSetSizes(*vv,n,N);
601: VecSetBlockSize(*vv,bs);
602: VecCreate_MPI_Private(*vv,PETSC_TRUE,nghost*bs,array);
603: w = (Vec_MPI*)(*vv)->data;
604: /* Create local representation */
605: VecGetArray(*vv,&larray);
606: VecCreateSeqWithArray(PETSC_COMM_SELF,bs,n+bs*nghost,larray,&w->localrep);
607: PetscLogObjectParent((PetscObject)*vv,(PetscObject)w->localrep);
608: VecRestoreArray(*vv,&larray);
610: /*
611: Create scatter context for scattering (updating) ghost values
612: */
613: ISCreateBlock(comm,bs,nghost,ghosts,PETSC_COPY_VALUES,&from);
614: ISCreateStride(PETSC_COMM_SELF,bs*nghost,n,1,&to);
615: VecScatterCreate(*vv,from,w->localrep,to,&w->localupdate);
616: PetscLogObjectParent((PetscObject)*vv,(PetscObject)w->localupdate);
617: ISDestroy(&to);
618: ISDestroy(&from);
620: /* set local to global mapping for ghosted vector */
621: nb = n/bs;
622: PetscMalloc1(nb+nghost,&indices);
623: VecGetOwnershipRange(*vv,&rstart,NULL);
625: for (i=0; i<nb; i++) indices[i] = rstart + i;
626: for (i=0; i<nghost; i++) indices[nb+i] = ghosts[i];
628: ISLocalToGlobalMappingCreate(comm,bs,nb+nghost,indices,PETSC_OWN_POINTER,<og);
629: VecSetLocalToGlobalMapping(*vv,ltog);
630: ISLocalToGlobalMappingDestroy(<og);
631: return(0);
632: }
636: /*@
637: VecCreateGhostBlock - Creates a parallel vector with ghost padding on each processor.
638: The indicing of the ghost points is done with blocks.
640: Collective on MPI_Comm
642: Input Parameters:
643: + comm - the MPI communicator to use
644: . bs - the block size
645: . n - local vector length
646: . N - global vector length (or PETSC_DECIDE to have calculated if n is given)
647: . nghost - number of local ghost blocks
648: - ghosts - global indices of ghost blocks, counts are by block, not by individual index, these do not need to be in increasing order (sorted)
650: Output Parameter:
651: . vv - the global vector representation (without ghost points as part of vector)
653: Notes:
654: Use VecGhostGetLocalForm() to access the local, ghosted representation
655: of the vector.
657: n is the local vector size (total local size not the number of blocks) while nghost
658: is the number of blocks in the ghost portion, i.e. the number of elements in the ghost
659: portion is bs*nghost
661: Level: advanced
663: Concepts: vectors^ghosted
665: .seealso: VecCreateSeq(), VecCreate(), VecDuplicate(), VecDuplicateVecs(), VecCreateMPI(),
666: VecGhostGetLocalForm(), VecGhostRestoreLocalForm(),
667: VecCreateGhostWithArray(), VecCreateMPIWithArray(), VecCreateGhostBlockWithArray()
669: @*/
670: PetscErrorCode VecCreateGhostBlock(MPI_Comm comm,PetscInt bs,PetscInt n,PetscInt N,PetscInt nghost,const PetscInt ghosts[],Vec *vv)
671: {
675: VecCreateGhostBlockWithArray(comm,bs,n,N,nghost,ghosts,0,vv);
676: return(0);
677: }