Actual source code: pbvec.c
petsc-3.4.5 2014-06-29
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,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(*v,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 = 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: };
171: /*
172: VecCreate_MPI_Private - Basic create routine called by VecCreate_MPI() (i.e. VecCreateMPI()),
173: VecCreateMPIWithArray(), VecCreate_Shared() (i.e. VecCreateShared()), VecCreateGhost(),
174: VecDuplicate_MPI(), VecCreateGhostWithArray(), VecDuplicate_MPI(), and VecDuplicate_Shared()
176: If alloc is true and array is NULL then this routine allocates the space, otherwise
177: no space is allocated.
178: */
179: PetscErrorCode VecCreate_MPI_Private(Vec v,PetscBool alloc,PetscInt nghost,const PetscScalar array[])
180: {
181: Vec_MPI *s;
185: PetscNewLog(v,Vec_MPI,&s);
186: v->data = (void*)s;
187: PetscMemcpy(v->ops,&DvOps,sizeof(DvOps));
188: s->nghost = nghost;
189: v->petscnative = PETSC_TRUE;
191: PetscLayoutSetUp(v->map);
193: s->array = (PetscScalar*)array;
194: s->array_allocated = 0;
195: if (alloc && !array) {
196: PetscInt n = v->map->n+nghost;
197: PetscMalloc(n*sizeof(PetscScalar),&s->array);
198: PetscLogObjectMemory(v,n*sizeof(PetscScalar));
199: PetscMemzero(s->array,v->map->n*sizeof(PetscScalar));
200: s->array_allocated = s->array;
201: }
203: /* By default parallel vectors do not have local representation */
204: s->localrep = 0;
205: s->localupdate = 0;
207: v->stash.insertmode = NOT_SET_VALUES;
208: /* create the stashes. The block-size for bstash is set later when
209: VecSetValuesBlocked is called.
210: */
211: VecStashCreate_Private(PetscObjectComm((PetscObject)v),1,&v->stash);
212: VecStashCreate_Private(PetscObjectComm((PetscObject)v),v->map->bs,&v->bstash);
214: #if defined(PETSC_HAVE_MATLAB_ENGINE)
215: PetscObjectComposeFunction((PetscObject)v,"PetscMatlabEnginePut_C",VecMatlabEnginePut_Default);
216: PetscObjectComposeFunction((PetscObject)v,"PetscMatlabEngineGet_C",VecMatlabEngineGet_Default);
217: #endif
218: PetscObjectChangeTypeName((PetscObject)v,VECMPI);
219: return(0);
220: }
222: /*MC
223: VECMPI - VECMPI = "mpi" - The basic parallel vector
225: Options Database Keys:
226: . -vec_type mpi - sets the vector type to VECMPI during a call to VecSetFromOptions()
228: Level: beginner
230: .seealso: VecCreate(), VecSetType(), VecSetFromOptions(), VecCreateMpiWithArray(), VECMPI, VecType, VecCreateMPI(), VecCreateMpi()
231: M*/
235: PETSC_EXTERN PetscErrorCode VecCreate_MPI(Vec vv)
236: {
240: VecCreate_MPI_Private(vv,PETSC_TRUE,0,0);
241: return(0);
242: }
244: /*MC
245: VECSTANDARD = "standard" - A VECSEQ on one process and VECMPI on more than one process
247: Options Database Keys:
248: . -vec_type standard - sets a vector type to standard on calls to VecSetFromOptions()
250: Level: beginner
252: .seealso: VecCreateSeq(), VecCreateMPI()
253: M*/
257: PETSC_EXTERN PetscErrorCode VecCreate_Standard(Vec v)
258: {
260: PetscMPIInt size;
263: MPI_Comm_size(PetscObjectComm((PetscObject)v),&size);
264: if (size == 1) {
265: VecSetType(v,VECSEQ);
266: } else {
267: VecSetType(v,VECMPI);
268: }
269: return(0);
270: }
274: /*@C
275: VecCreateMPIWithArray - Creates a parallel, array-style vector,
276: where the user provides the array space to store the vector values.
278: Collective on MPI_Comm
280: Input Parameters:
281: + comm - the MPI communicator to use
282: . bs - block size, same meaning as VecSetBlockSize()
283: . n - local vector length, cannot be PETSC_DECIDE
284: . N - global vector length (or PETSC_DECIDE to have calculated)
285: - array - the user provided array to store the vector values
287: Output Parameter:
288: . vv - the vector
290: Notes:
291: Use VecDuplicate() or VecDuplicateVecs() to form additional vectors of the
292: same type as an existing vector.
294: If the user-provided array is NULL, then VecPlaceArray() can be used
295: at a later stage to SET the array for storing the vector values.
297: PETSc does NOT free the array when the vector is destroyed via VecDestroy().
298: The user should not free the array until the vector is destroyed.
300: Level: intermediate
302: Concepts: vectors^creating with array
304: .seealso: VecCreateSeqWithArray(), VecCreate(), VecDuplicate(), VecDuplicateVecs(), VecCreateGhost(),
305: VecCreateMPI(), VecCreateGhostWithArray(), VecPlaceArray()
307: @*/
308: PetscErrorCode VecCreateMPIWithArray(MPI_Comm comm,PetscInt bs,PetscInt n,PetscInt N,const PetscScalar array[],Vec *vv)
309: {
313: if (n == PETSC_DECIDE) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must set local size of vector");
314: PetscSplitOwnership(comm,&n,&N);
315: VecCreate(comm,vv);
316: VecSetSizes(*vv,n,N);
317: VecSetBlockSize(*vv,bs);
318: VecCreate_MPI_Private(*vv,PETSC_FALSE,0,array);
319: return(0);
320: }
324: /*@C
325: VecCreateGhostWithArray - Creates a parallel vector with ghost padding on each processor;
326: the caller allocates the array space.
328: Collective on MPI_Comm
330: Input Parameters:
331: + comm - the MPI communicator to use
332: . n - local vector length
333: . N - global vector length (or PETSC_DECIDE to have calculated if n is given)
334: . nghost - number of local ghost points
335: . ghosts - global indices of ghost points (or NULL if not needed), these do not need to be in increasing order (sorted)
336: - array - the space to store the vector values (as long as n + nghost)
338: Output Parameter:
339: . vv - the global vector representation (without ghost points as part of vector)
341: Notes:
342: Use VecGhostGetLocalForm() to access the local, ghosted representation
343: of the vector.
345: This also automatically sets the ISLocalToGlobalMapping() for this vector.
347: Level: advanced
349: Concepts: vectors^creating with array
350: Concepts: vectors^ghosted
352: .seealso: VecCreate(), VecGhostGetLocalForm(), VecGhostRestoreLocalForm(),
353: VecCreateGhost(), VecCreateSeqWithArray(), VecCreateMPIWithArray(),
354: VecCreateGhostBlock(), VecCreateGhostBlockWithArray(), VecMPISetGhost()
356: @*/
357: PetscErrorCode VecCreateGhostWithArray(MPI_Comm comm,PetscInt n,PetscInt N,PetscInt nghost,const PetscInt ghosts[],const PetscScalar array[],Vec *vv)
358: {
359: PetscErrorCode ierr;
360: Vec_MPI *w;
361: PetscScalar *larray;
362: IS from,to;
363: ISLocalToGlobalMapping ltog;
364: PetscInt rstart,i,*indices;
367: *vv = 0;
369: if (n == PETSC_DECIDE) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must set local size");
370: if (nghost == PETSC_DECIDE) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must set local ghost size");
371: if (nghost < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Ghost length must be >= 0");
372: PetscSplitOwnership(comm,&n,&N);
373: /* Create global representation */
374: VecCreate(comm,vv);
375: VecSetSizes(*vv,n,N);
376: VecCreate_MPI_Private(*vv,PETSC_TRUE,nghost,array);
377: w = (Vec_MPI*)(*vv)->data;
378: /* Create local representation */
379: VecGetArray(*vv,&larray);
380: VecCreateSeqWithArray(PETSC_COMM_SELF,1,n+nghost,larray,&w->localrep);
381: PetscLogObjectParent(*vv,w->localrep);
382: VecRestoreArray(*vv,&larray);
384: /*
385: Create scatter context for scattering (updating) ghost values
386: */
387: ISCreateGeneral(comm,nghost,ghosts,PETSC_COPY_VALUES,&from);
388: ISCreateStride(PETSC_COMM_SELF,nghost,n,1,&to);
389: VecScatterCreate(*vv,from,w->localrep,to,&w->localupdate);
390: PetscLogObjectParent(*vv,w->localupdate);
391: ISDestroy(&to);
392: ISDestroy(&from);
394: /* set local to global mapping for ghosted vector */
395: PetscMalloc((n+nghost)*sizeof(PetscInt),&indices);
396: VecGetOwnershipRange(*vv,&rstart,NULL);
397: for (i=0; i<n; i++) {
398: indices[i] = rstart + i;
399: }
400: for (i=0; i<nghost; i++) {
401: indices[n+i] = ghosts[i];
402: }
403: ISLocalToGlobalMappingCreate(comm,n+nghost,indices,PETSC_OWN_POINTER,<og);
404: VecSetLocalToGlobalMapping(*vv,ltog);
405: ISLocalToGlobalMappingDestroy(<og);
406: return(0);
407: }
411: /*@
412: VecCreateGhost - Creates a parallel vector with ghost padding on each processor.
414: Collective on MPI_Comm
416: Input Parameters:
417: + comm - the MPI communicator to use
418: . n - local vector length
419: . N - global vector length (or PETSC_DECIDE to have calculated if n is given)
420: . nghost - number of local ghost points
421: - ghosts - global indices of ghost points, these do not need to be in increasing order (sorted)
423: Output Parameter:
424: . vv - the global vector representation (without ghost points as part of vector)
426: Notes:
427: Use VecGhostGetLocalForm() to access the local, ghosted representation
428: of the vector.
430: This also automatically sets the ISLocalToGlobalMapping() for this vector.
432: Level: advanced
434: Concepts: vectors^ghosted
436: .seealso: VecCreateSeq(), VecCreate(), VecDuplicate(), VecDuplicateVecs(), VecCreateMPI(),
437: VecGhostGetLocalForm(), VecGhostRestoreLocalForm(), VecGhostUpdateBegin(),
438: VecCreateGhostWithArray(), VecCreateMPIWithArray(), VecGhostUpdateEnd(),
439: VecCreateGhostBlock(), VecCreateGhostBlockWithArray(), VecMPISetGhost()
441: @*/
442: PetscErrorCode VecCreateGhost(MPI_Comm comm,PetscInt n,PetscInt N,PetscInt nghost,const PetscInt ghosts[],Vec *vv)
443: {
447: VecCreateGhostWithArray(comm,n,N,nghost,ghosts,0,vv);
448: return(0);
449: }
453: /*@
454: VecMPISetGhost - Sets the ghost points for an MPI ghost vector
456: Collective on Vec
458: Input Parameters:
459: + vv - the MPI vector
460: . nghost - number of local ghost points
461: - ghosts - global indices of ghost points, these do not need to be in increasing order (sorted)
464: Notes:
465: Use VecGhostGetLocalForm() to access the local, ghosted representation
466: of the vector.
468: This also automatically sets the ISLocalToGlobalMapping() for this vector.
470: You must call this AFTER you have set the type of the vector (with VecSetType()) and the size (with VecSetSizes()).
472: Level: advanced
474: Concepts: vectors^ghosted
476: .seealso: VecCreateSeq(), VecCreate(), VecDuplicate(), VecDuplicateVecs(), VecCreateMPI(),
477: VecGhostGetLocalForm(), VecGhostRestoreLocalForm(), VecGhostUpdateBegin(),
478: VecCreateGhostWithArray(), VecCreateMPIWithArray(), VecGhostUpdateEnd(),
479: VecCreateGhostBlock(), VecCreateGhostBlockWithArray()
481: @*/
482: PetscErrorCode VecMPISetGhost(Vec vv,PetscInt nghost,const PetscInt ghosts[])
483: {
485: PetscBool flg;
488: PetscObjectTypeCompare((PetscObject)vv,VECMPI,&flg);
489: /* if already fully existant VECMPI then basically destroy it and rebuild with ghosting */
490: if (flg) {
491: PetscInt n,N;
492: Vec_MPI *w;
493: PetscScalar *larray;
494: IS from,to;
495: ISLocalToGlobalMapping ltog;
496: PetscInt rstart,i,*indices;
497: MPI_Comm comm;
499: PetscObjectGetComm((PetscObject)vv,&comm);
500: n = vv->map->n;
501: N = vv->map->N;
502: (*vv->ops->destroy)(vv);
503: VecSetSizes(vv,n,N);
504: VecCreate_MPI_Private(vv,PETSC_TRUE,nghost,NULL);
505: w = (Vec_MPI*)(vv)->data;
506: /* Create local representation */
507: VecGetArray(vv,&larray);
508: VecCreateSeqWithArray(PETSC_COMM_SELF,1,n+nghost,larray,&w->localrep);
509: PetscLogObjectParent(vv,w->localrep);
510: VecRestoreArray(vv,&larray);
512: /*
513: Create scatter context for scattering (updating) ghost values
514: */
515: ISCreateGeneral(comm,nghost,ghosts,PETSC_COPY_VALUES,&from);
516: ISCreateStride(PETSC_COMM_SELF,nghost,n,1,&to);
517: VecScatterCreate(vv,from,w->localrep,to,&w->localupdate);
518: PetscLogObjectParent(vv,w->localupdate);
519: ISDestroy(&to);
520: ISDestroy(&from);
522: /* set local to global mapping for ghosted vector */
523: PetscMalloc((n+nghost)*sizeof(PetscInt),&indices);
524: VecGetOwnershipRange(vv,&rstart,NULL);
526: for (i=0; i<n; i++) indices[i] = rstart + i;
527: for (i=0; i<nghost; i++) indices[n+i] = ghosts[i];
529: ISLocalToGlobalMappingCreate(comm,n+nghost,indices,PETSC_OWN_POINTER,<og);
530: VecSetLocalToGlobalMapping(vv,ltog);
531: ISLocalToGlobalMappingDestroy(<og);
532: } else if (vv->ops->create == VecCreate_MPI) SETERRQ(PetscObjectComm((PetscObject)vv),PETSC_ERR_ARG_WRONGSTATE,"Must set local or global size before setting ghosting");
533: else if (!((PetscObject)vv)->type_name) SETERRQ(PetscObjectComm((PetscObject)vv),PETSC_ERR_ARG_WRONGSTATE,"Must set type to VECMPI before ghosting");
534: return(0);
535: }
538: /* ------------------------------------------------------------------------------------------*/
541: /*@C
542: VecCreateGhostBlockWithArray - Creates a parallel vector with ghost padding on each processor;
543: the caller allocates the array space. Indices in the ghost region are based on blocks.
545: Collective on MPI_Comm
547: Input Parameters:
548: + comm - the MPI communicator to use
549: . bs - block size
550: . n - local vector length
551: . N - global vector length (or PETSC_DECIDE to have calculated if n is given)
552: . nghost - number of local ghost blocks
553: . 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)
554: - array - the space to store the vector values (as long as n + nghost*bs)
556: Output Parameter:
557: . vv - the global vector representation (without ghost points as part of vector)
559: Notes:
560: Use VecGhostGetLocalForm() to access the local, ghosted representation
561: of the vector.
563: n is the local vector size (total local size not the number of blocks) while nghost
564: is the number of blocks in the ghost portion, i.e. the number of elements in the ghost
565: portion is bs*nghost
567: Level: advanced
569: Concepts: vectors^creating ghosted
570: Concepts: vectors^creating with array
572: .seealso: VecCreate(), VecGhostGetLocalForm(), VecGhostRestoreLocalForm(),
573: VecCreateGhost(), VecCreateSeqWithArray(), VecCreateMPIWithArray(),
574: VecCreateGhostWithArray(), VecCreateGhostBlock()
576: @*/
577: PetscErrorCode VecCreateGhostBlockWithArray(MPI_Comm comm,PetscInt bs,PetscInt n,PetscInt N,PetscInt nghost,const PetscInt ghosts[],const PetscScalar array[],Vec *vv)
578: {
579: PetscErrorCode ierr;
580: Vec_MPI *w;
581: PetscScalar *larray;
582: IS from,to;
583: ISLocalToGlobalMapping ltog;
584: PetscInt rstart,i,nb,*indices;
587: *vv = 0;
589: if (n == PETSC_DECIDE) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must set local size");
590: if (nghost == PETSC_DECIDE) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must set local ghost size");
591: if (nghost < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Ghost length must be >= 0");
592: if (n % bs) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Local size must be a multiple of block size");
593: PetscSplitOwnership(comm,&n,&N);
594: /* Create global representation */
595: VecCreate(comm,vv);
596: VecSetSizes(*vv,n,N);
597: VecSetBlockSize(*vv,bs);
598: VecCreate_MPI_Private(*vv,PETSC_TRUE,nghost*bs,array);
599: w = (Vec_MPI*)(*vv)->data;
600: /* Create local representation */
601: VecGetArray(*vv,&larray);
602: VecCreateSeqWithArray(PETSC_COMM_SELF,bs,n+bs*nghost,larray,&w->localrep);
603: PetscLogObjectParent(*vv,w->localrep);
604: VecRestoreArray(*vv,&larray);
606: /*
607: Create scatter context for scattering (updating) ghost values
608: */
609: ISCreateBlock(comm,bs,nghost,ghosts,PETSC_COPY_VALUES,&from);
610: ISCreateStride(PETSC_COMM_SELF,bs*nghost,n,1,&to);
611: VecScatterCreate(*vv,from,w->localrep,to,&w->localupdate);
612: PetscLogObjectParent(*vv,w->localupdate);
613: ISDestroy(&to);
614: ISDestroy(&from);
616: /* set local to global mapping for ghosted vector */
617: nb = n/bs;
618: PetscMalloc((nb+nghost)*sizeof(PetscInt),&indices);
619: VecGetOwnershipRange(*vv,&rstart,NULL);
621: for (i=0; i<nb; i++) indices[i] = rstart + i*bs;
622: for (i=0; i<nghost; i++) indices[nb+i] = ghosts[i];
624: ISLocalToGlobalMappingCreate(comm,nb+nghost,indices,PETSC_OWN_POINTER,<og);
625: VecSetLocalToGlobalMappingBlock(*vv,ltog);
626: ISLocalToGlobalMappingDestroy(<og);
627: return(0);
628: }
632: /*@
633: VecCreateGhostBlock - Creates a parallel vector with ghost padding on each processor.
634: The indicing of the ghost points is done with blocks.
636: Collective on MPI_Comm
638: Input Parameters:
639: + comm - the MPI communicator to use
640: . bs - the block size
641: . n - local vector length
642: . N - global vector length (or PETSC_DECIDE to have calculated if n is given)
643: . nghost - number of local ghost blocks
644: - ghosts - global indices of ghost blocks, counts are by block, not by individual index, these do not need to be in increasing order (sorted)
646: Output Parameter:
647: . vv - the global vector representation (without ghost points as part of vector)
649: Notes:
650: Use VecGhostGetLocalForm() to access the local, ghosted representation
651: of the vector.
653: n is the local vector size (total local size not the number of blocks) while nghost
654: is the number of blocks in the ghost portion, i.e. the number of elements in the ghost
655: portion is bs*nghost
657: Level: advanced
659: Concepts: vectors^ghosted
661: .seealso: VecCreateSeq(), VecCreate(), VecDuplicate(), VecDuplicateVecs(), VecCreateMPI(),
662: VecGhostGetLocalForm(), VecGhostRestoreLocalForm(),
663: VecCreateGhostWithArray(), VecCreateMPIWithArray(), VecCreateGhostBlockWithArray()
665: @*/
666: PetscErrorCode VecCreateGhostBlock(MPI_Comm comm,PetscInt bs,PetscInt n,PetscInt N,PetscInt nghost,const PetscInt ghosts[],Vec *vv)
667: {
671: VecCreateGhostBlockWithArray(comm,bs,n,N,nghost,ghosts,0,vv);
672: return(0);
673: }