Actual source code: vscat.c
petsc-3.4.5 2014-06-29
2: /*
3: Code for creating scatters between vectors. This file
4: includes the code for scattering between sequential vectors and
5: some special cases for parallel scatters.
6: */
8: #include <petsc-private/vecimpl.h> /*I "petscvec.h" I*/
10: /* Logging support */
11: PetscClassId VEC_SCATTER_CLASSID;
13: #if defined(PETSC_USE_DEBUG)
14: /*
15: Checks if any indices are less than zero and generates an error
16: */
19: static PetscErrorCode VecScatterCheckIndices_Private(PetscInt nmax,PetscInt n,const PetscInt *idx)
20: {
21: PetscInt i;
24: for (i=0; i<n; i++) {
25: if (idx[i] < 0) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Negative index %D at %D location",idx[i],i);
26: if (idx[i] >= nmax) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Index %D at %D location greater than max %D",idx[i],i,nmax);
27: }
28: return(0);
29: }
30: #endif
32: /*
33: This is special scatter code for when the entire parallel vector is copied to each processor.
35: This code was written by Cameron Cooper, Occidental College, Fall 1995,
36: will working at ANL as a SERS student.
37: */
40: PetscErrorCode VecScatterBegin_MPI_ToAll(VecScatter ctx,Vec x,Vec y,InsertMode addv,ScatterMode mode)
41: {
43: PetscInt yy_n,xx_n;
44: PetscScalar *xv,*yv;
47: VecGetLocalSize(y,&yy_n);
48: VecGetLocalSize(x,&xx_n);
49: VecGetArray(y,&yv);
50: if (x != y) {VecGetArray(x,&xv);}
52: if (mode & SCATTER_REVERSE) {
53: PetscScalar *xvt,*xvt2;
54: VecScatter_MPI_ToAll *scat = (VecScatter_MPI_ToAll*)ctx->todata;
55: PetscInt i;
56: PetscMPIInt *disply = scat->displx;
58: if (addv == INSERT_VALUES) {
59: PetscInt rstart,rend;
60: /*
61: copy the correct part of the local vector into the local storage of
62: the MPI one. Note: this operation only makes sense if all the local
63: vectors have the same values
64: */
65: VecGetOwnershipRange(y,&rstart,&rend);
66: PetscMemcpy(yv,xv+rstart,yy_n*sizeof(PetscScalar));
67: } else {
68: MPI_Comm comm;
69: PetscMPIInt rank;
70: PetscObjectGetComm((PetscObject)y,&comm);
71: MPI_Comm_rank(comm,&rank);
72: if (scat->work1) xvt = scat->work1;
73: else {
74: PetscMalloc(xx_n*sizeof(PetscScalar),&xvt);
75: scat->work1 = xvt;
76: }
77: if (!rank) { /* I am the zeroth processor, values are accumulated here */
78: if (scat->work2) xvt2 = scat->work2;
79: else {
80: PetscMalloc(xx_n*sizeof(PetscScalar),&xvt2);
81: scat->work2 = xvt2;
82: }
83: MPI_Gatherv(yv,yy_n,MPIU_SCALAR,xvt2,scat->count,disply,MPIU_SCALAR,0,PetscObjectComm((PetscObject)ctx));
84: MPI_Reduce(xv,xvt,xx_n,MPIU_SCALAR,MPIU_SUM,0,PetscObjectComm((PetscObject)ctx));
85: if (addv == ADD_VALUES) {
86: for (i=0; i<xx_n; i++) xvt[i] += xvt2[i];
87: #if !defined(PETSC_USE_COMPLEX)
88: } else if (addv == MAX_VALUES) {
89: for (i=0; i<xx_n; i++) xvt[i] = PetscMax(xvt[i],xvt2[i]);
90: #endif
91: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Wrong insert option");
92: MPI_Scatterv(xvt,scat->count,disply,MPIU_SCALAR,yv,yy_n,MPIU_SCALAR,0,PetscObjectComm((PetscObject)ctx));
93: } else {
94: MPI_Gatherv(yv,yy_n,MPIU_SCALAR,0, 0,0,MPIU_SCALAR,0,PetscObjectComm((PetscObject)ctx));
95: MPI_Reduce(xv,xvt,xx_n,MPIU_SCALAR,MPIU_SUM,0,PetscObjectComm((PetscObject)ctx));
96: MPI_Scatterv(0,scat->count,disply,MPIU_SCALAR,yv,yy_n,MPIU_SCALAR,0,PetscObjectComm((PetscObject)ctx));
97: }
98: }
99: } else {
100: PetscScalar *yvt;
101: VecScatter_MPI_ToAll *scat = (VecScatter_MPI_ToAll*)ctx->todata;
102: PetscInt i;
103: PetscMPIInt *displx = scat->displx;
105: if (addv == INSERT_VALUES) {
106: MPI_Allgatherv(xv,xx_n,MPIU_SCALAR,yv,scat->count,displx,MPIU_SCALAR,PetscObjectComm((PetscObject)ctx));
107: } else {
108: if (scat->work1) yvt = scat->work1;
109: else {
110: PetscMalloc(yy_n*sizeof(PetscScalar),&yvt);
111: scat->work1 = yvt;
112: }
113: MPI_Allgatherv(xv,xx_n,MPIU_SCALAR,yvt,scat->count,displx,MPIU_SCALAR,PetscObjectComm((PetscObject)ctx));
114: if (addv == ADD_VALUES) {
115: for (i=0; i<yy_n; i++) yv[i] += yvt[i];
116: #if !defined(PETSC_USE_COMPLEX)
117: } else if (addv == MAX_VALUES) {
118: for (i=0; i<yy_n; i++) yv[i] = PetscMax(yv[i],yvt[i]);
119: #endif
120: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Wrong insert option");
121: }
122: }
123: VecRestoreArray(y,&yv);
124: if (x != y) {VecRestoreArray(x,&xv);}
125: return(0);
126: }
130: PetscErrorCode VecScatterView_MPI_ToAll(VecScatter in,PetscViewer viewer)
131: {
133: PetscBool isascii;
136: PetscObjectTypeCompare((PetscObject)in,PETSCVIEWERASCII,&isascii);
137: if (isascii) {
138: PetscViewerASCIIPrintf(viewer,"Entire parallel vector is copied to each process\n");
139: }
140: return(0);
141: }
143: /*
144: This is special scatter code for when the entire parallel vector is copied to processor 0.
146: */
149: PetscErrorCode VecScatterBegin_MPI_ToOne(VecScatter ctx,Vec x,Vec y,InsertMode addv,ScatterMode mode)
150: {
152: PetscMPIInt rank;
153: PetscInt yy_n,xx_n;
154: PetscScalar *xv,*yv;
155: MPI_Comm comm;
158: VecGetLocalSize(y,&yy_n);
159: VecGetLocalSize(x,&xx_n);
160: VecGetArray(x,&xv);
161: VecGetArray(y,&yv);
163: PetscObjectGetComm((PetscObject)x,&comm);
164: MPI_Comm_rank(comm,&rank);
166: /* -------- Reverse scatter; spread from processor 0 to other processors */
167: if (mode & SCATTER_REVERSE) {
168: PetscScalar *yvt;
169: VecScatter_MPI_ToAll *scat = (VecScatter_MPI_ToAll*)ctx->todata;
170: PetscInt i;
171: PetscMPIInt *disply = scat->displx;
173: if (addv == INSERT_VALUES) {
174: MPI_Scatterv(xv,scat->count,disply,MPIU_SCALAR,yv,yy_n,MPIU_SCALAR,0,PetscObjectComm((PetscObject)ctx));
175: } else {
176: if (scat->work2) yvt = scat->work2;
177: else {
178: PetscMalloc(xx_n*sizeof(PetscScalar),&yvt);
179: scat->work2 = yvt;
180: }
181: MPI_Scatterv(xv,scat->count,disply,MPIU_SCALAR,yvt,yy_n,MPIU_SCALAR,0,PetscObjectComm((PetscObject)ctx));
182: if (addv == ADD_VALUES) {
183: for (i=0; i<yy_n; i++) yv[i] += yvt[i];
184: #if !defined(PETSC_USE_COMPLEX)
185: } else if (addv == MAX_VALUES) {
186: for (i=0; i<yy_n; i++) yv[i] = PetscMax(yv[i],yvt[i]);
187: #endif
188: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Wrong insert option");
189: }
190: /* --------- Forward scatter; gather all values onto processor 0 */
191: } else {
192: PetscScalar *yvt = 0;
193: VecScatter_MPI_ToAll *scat = (VecScatter_MPI_ToAll*)ctx->todata;
194: PetscInt i;
195: PetscMPIInt *displx = scat->displx;
197: if (addv == INSERT_VALUES) {
198: MPI_Gatherv(xv,xx_n,MPIU_SCALAR,yv,scat->count,displx,MPIU_SCALAR,0,PetscObjectComm((PetscObject)ctx));
199: } else {
200: if (!rank) {
201: if (scat->work1) yvt = scat->work1;
202: else {
203: PetscMalloc(yy_n*sizeof(PetscScalar),&yvt);
204: scat->work1 = yvt;
205: }
206: }
207: MPI_Gatherv(xv,xx_n,MPIU_SCALAR,yvt,scat->count,displx,MPIU_SCALAR,0,PetscObjectComm((PetscObject)ctx));
208: if (!rank) {
209: if (addv == ADD_VALUES) {
210: for (i=0; i<yy_n; i++) yv[i] += yvt[i];
211: #if !defined(PETSC_USE_COMPLEX)
212: } else if (addv == MAX_VALUES) {
213: for (i=0; i<yy_n; i++) yv[i] = PetscMax(yv[i],yvt[i]);
214: #endif
215: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Wrong insert option");
216: }
217: }
218: }
219: VecRestoreArray(x,&xv);
220: VecRestoreArray(y,&yv);
221: return(0);
222: }
224: /*
225: The follow to are used for both VecScatterBegin_MPI_ToAll() and VecScatterBegin_MPI_ToOne()
226: */
229: PetscErrorCode VecScatterDestroy_MPI_ToAll(VecScatter ctx)
230: {
231: VecScatter_MPI_ToAll *scat = (VecScatter_MPI_ToAll*)ctx->todata;
232: PetscErrorCode ierr;
235: PetscFree(scat->work1);
236: PetscFree(scat->work2);
237: PetscFree3(ctx->todata,scat->count,scat->displx);
238: return(0);
239: }
243: PetscErrorCode VecScatterDestroy_SGToSG(VecScatter ctx)
244: {
248: PetscFree2(((VecScatter_Seq_General*)ctx->todata)->vslots,((VecScatter_Seq_General*)ctx->fromdata)->vslots);
249: PetscFree2(ctx->todata,ctx->fromdata);
250: return(0);
251: }
255: PetscErrorCode VecScatterDestroy_SGToSS(VecScatter ctx)
256: {
260: PetscFree(((VecScatter_Seq_General*)ctx->fromdata)->vslots);
261: PetscFree2(ctx->todata,ctx->fromdata);
262: return(0);
263: }
267: PetscErrorCode VecScatterDestroy_SSToSG(VecScatter ctx)
268: {
272: PetscFree(((VecScatter_Seq_General*)ctx->todata)->vslots);
273: PetscFree2(ctx->todata,ctx->fromdata);
274: return(0);
275: }
279: PetscErrorCode VecScatterDestroy_SSToSS(VecScatter ctx)
280: {
284: PetscFree2(ctx->todata,ctx->fromdata);
285: return(0);
286: }
288: /* -------------------------------------------------------------------------*/
291: PetscErrorCode VecScatterCopy_MPI_ToAll(VecScatter in,VecScatter out)
292: {
293: VecScatter_MPI_ToAll *in_to = (VecScatter_MPI_ToAll*)in->todata,*sto;
294: PetscErrorCode ierr;
295: PetscMPIInt size,*count,*displx;
296: PetscInt i;
299: out->begin = in->begin;
300: out->end = in->end;
301: out->copy = in->copy;
302: out->destroy = in->destroy;
303: out->view = in->view;
305: MPI_Comm_size(PetscObjectComm((PetscObject)out),&size);
306: PetscMalloc3(1,VecScatter_MPI_ToAll,&sto,size,PetscMPIInt,&count,size,PetscMPIInt,&displx);
307: sto->type = in_to->type;
308: sto->count = count;
309: sto->displx = displx;
310: for (i=0; i<size; i++) {
311: sto->count[i] = in_to->count[i];
312: sto->displx[i] = in_to->displx[i];
313: }
314: sto->work1 = 0;
315: sto->work2 = 0;
316: out->todata = (void*)sto;
317: out->fromdata = (void*)0;
318: return(0);
319: }
321: /* --------------------------------------------------------------------------------------*/
322: /*
323: Scatter: sequential general to sequential general
324: */
327: PetscErrorCode VecScatterBegin_SGToSG(VecScatter ctx,Vec x,Vec y,InsertMode addv,ScatterMode mode)
328: {
329: VecScatter_Seq_General *gen_to = (VecScatter_Seq_General*)ctx->todata;
330: VecScatter_Seq_General *gen_from = (VecScatter_Seq_General*)ctx->fromdata;
331: PetscErrorCode ierr;
332: PetscInt i,n = gen_from->n,*fslots,*tslots;
333: PetscScalar *xv,*yv;
336: VecGetArray(x,&xv);
337: if (x != y) {
338: VecGetArray(y,&yv);
339: } else yv = xv;
341: if (mode & SCATTER_REVERSE) {
342: gen_to = (VecScatter_Seq_General*)ctx->fromdata;
343: gen_from = (VecScatter_Seq_General*)ctx->todata;
344: }
345: fslots = gen_from->vslots;
346: tslots = gen_to->vslots;
348: if (addv == INSERT_VALUES) {
349: for (i=0; i<n; i++) yv[tslots[i]] = xv[fslots[i]];
350: } else if (addv == ADD_VALUES) {
351: for (i=0; i<n; i++) yv[tslots[i]] += xv[fslots[i]];
352: #if !defined(PETSC_USE_COMPLEX)
353: } else if (addv == MAX_VALUES) {
354: for (i=0; i<n; i++) yv[tslots[i]] = PetscMax(yv[tslots[i]],xv[fslots[i]]);
355: #endif
356: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Wrong insert option");
357: VecRestoreArray(x,&xv);
358: if (x != y) {VecRestoreArray(y,&yv);}
359: return(0);
360: }
362: /*
363: Scatter: sequential general to sequential stride 1
364: */
367: PetscErrorCode VecScatterBegin_SGToSS_Stride1(VecScatter ctx,Vec x,Vec y,InsertMode addv,ScatterMode mode)
368: {
369: VecScatter_Seq_Stride *gen_to = (VecScatter_Seq_Stride*)ctx->todata;
370: VecScatter_Seq_General *gen_from = (VecScatter_Seq_General*)ctx->fromdata;
371: PetscInt i,n = gen_from->n,*fslots = gen_from->vslots;
372: PetscErrorCode ierr;
373: PetscInt first = gen_to->first;
374: PetscScalar *xv,*yv;
377: VecGetArray(x,&xv);
378: if (x != y) {
379: VecGetArray(y,&yv);
380: } else yv = xv;
381: if (mode & SCATTER_REVERSE) {
382: xv += first;
383: if (addv == INSERT_VALUES) {
384: for (i=0; i<n; i++) yv[fslots[i]] = xv[i];
385: } else if (addv == ADD_VALUES) {
386: for (i=0; i<n; i++) yv[fslots[i]] += xv[i];
387: #if !defined(PETSC_USE_COMPLEX)
388: } else if (addv == MAX_VALUES) {
389: for (i=0; i<n; i++) yv[fslots[i]] = PetscMax(yv[fslots[i]],xv[i]);
390: #endif
391: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Wrong insert option");
392: } else {
393: yv += first;
394: if (addv == INSERT_VALUES) {
395: for (i=0; i<n; i++) yv[i] = xv[fslots[i]];
396: } else if (addv == ADD_VALUES) {
397: for (i=0; i<n; i++) yv[i] += xv[fslots[i]];
398: #if !defined(PETSC_USE_COMPLEX)
399: } else if (addv == MAX_VALUES) {
400: for (i=0; i<n; i++) yv[i] = PetscMax(yv[i],xv[fslots[i]]);
401: #endif
402: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Wrong insert option");
403: }
404: VecRestoreArray(x,&xv);
405: if (x != y) {VecRestoreArray(y,&yv);}
406: return(0);
407: }
409: /*
410: Scatter: sequential general to sequential stride
411: */
414: PetscErrorCode VecScatterBegin_SGToSS(VecScatter ctx,Vec x,Vec y,InsertMode addv,ScatterMode mode)
415: {
416: VecScatter_Seq_Stride *gen_to = (VecScatter_Seq_Stride*)ctx->todata;
417: VecScatter_Seq_General *gen_from = (VecScatter_Seq_General*)ctx->fromdata;
418: PetscInt i,n = gen_from->n,*fslots = gen_from->vslots;
419: PetscErrorCode ierr;
420: PetscInt first = gen_to->first,step = gen_to->step;
421: PetscScalar *xv,*yv;
424: VecGetArray(x,&xv);
425: if (x != y) {
426: VecGetArray(y,&yv);
427: } else yv = xv;
429: if (mode & SCATTER_REVERSE) {
430: if (addv == INSERT_VALUES) {
431: for (i=0; i<n; i++) yv[fslots[i]] = xv[first + i*step];
432: } else if (addv == ADD_VALUES) {
433: for (i=0; i<n; i++) yv[fslots[i]] += xv[first + i*step];
434: #if !defined(PETSC_USE_COMPLEX)
435: } else if (addv == MAX_VALUES) {
436: for (i=0; i<n; i++) yv[fslots[i]] = PetscMax(yv[fslots[i]],xv[first + i*step]);
437: #endif
438: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Wrong insert option");
439: } else {
440: if (addv == INSERT_VALUES) {
441: for (i=0; i<n; i++) yv[first + i*step] = xv[fslots[i]];
442: } else if (addv == ADD_VALUES) {
443: for (i=0; i<n; i++) yv[first + i*step] += xv[fslots[i]];
444: #if !defined(PETSC_USE_COMPLEX)
445: } else if (addv == MAX_VALUES) {
446: for (i=0; i<n; i++) yv[first + i*step] = PetscMax(yv[first + i*step],xv[fslots[i]]);
447: #endif
448: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Wrong insert option");
449: }
450: VecRestoreArray(x,&xv);
451: if (x != y) {VecRestoreArray(y,&yv);}
452: return(0);
453: }
455: /*
456: Scatter: sequential stride 1 to sequential general
457: */
460: PetscErrorCode VecScatterBegin_SSToSG_Stride1(VecScatter ctx,Vec x,Vec y,InsertMode addv,ScatterMode mode)
461: {
462: VecScatter_Seq_Stride *gen_from = (VecScatter_Seq_Stride*)ctx->fromdata;
463: VecScatter_Seq_General *gen_to = (VecScatter_Seq_General*)ctx->todata;
464: PetscInt i,n = gen_from->n,*fslots = gen_to->vslots;
465: PetscErrorCode ierr;
466: PetscInt first = gen_from->first;
467: PetscScalar *xv,*yv;
470: VecGetArray(x,&xv);
471: if (x != y) {
472: VecGetArray(y,&yv);
473: } else yv = xv;
475: if (mode & SCATTER_REVERSE) {
476: yv += first;
477: if (addv == INSERT_VALUES) {
478: for (i=0; i<n; i++) yv[i] = xv[fslots[i]];
479: } else if (addv == ADD_VALUES) {
480: for (i=0; i<n; i++) yv[i] += xv[fslots[i]];
481: #if !defined(PETSC_USE_COMPLEX)
482: } else if (addv == MAX_VALUES) {
483: for (i=0; i<n; i++) yv[i] = PetscMax(yv[i],xv[fslots[i]]);
484: #endif
485: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Wrong insert option");
486: } else {
487: xv += first;
488: if (addv == INSERT_VALUES) {
489: for (i=0; i<n; i++) yv[fslots[i]] = xv[i];
490: } else if (addv == ADD_VALUES) {
491: for (i=0; i<n; i++) yv[fslots[i]] += xv[i];
492: #if !defined(PETSC_USE_COMPLEX)
493: } else if (addv == MAX_VALUES) {
494: for (i=0; i<n; i++) yv[fslots[i]] = PetscMax(yv[fslots[i]],xv[i]);
495: #endif
496: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Wrong insert option");
497: }
498: VecRestoreArray(x,&xv);
499: if (x != y) {VecRestoreArray(y,&yv);}
500: return(0);
501: }
505: /*
506: Scatter: sequential stride to sequential general
507: */
508: PetscErrorCode VecScatterBegin_SSToSG(VecScatter ctx,Vec x,Vec y,InsertMode addv,ScatterMode mode)
509: {
510: VecScatter_Seq_Stride *gen_from = (VecScatter_Seq_Stride*)ctx->fromdata;
511: VecScatter_Seq_General *gen_to = (VecScatter_Seq_General*)ctx->todata;
512: PetscInt i,n = gen_from->n,*fslots = gen_to->vslots;
513: PetscErrorCode ierr;
514: PetscInt first = gen_from->first,step = gen_from->step;
515: PetscScalar *xv,*yv;
518: VecGetArray(x,&xv);
519: if (x != y) {
520: VecGetArray(y,&yv);
521: } else yv = xv;
523: if (mode & SCATTER_REVERSE) {
524: if (addv == INSERT_VALUES) {
525: for (i=0; i<n; i++) yv[first + i*step] = xv[fslots[i]];
526: } else if (addv == ADD_VALUES) {
527: for (i=0; i<n; i++) yv[first + i*step] += xv[fslots[i]];
528: #if !defined(PETSC_USE_COMPLEX)
529: } else if (addv == MAX_VALUES) {
530: for (i=0; i<n; i++) yv[first + i*step] = PetscMax(yv[first + i*step],xv[fslots[i]]);
531: #endif
532: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Wrong insert option");
533: } else {
534: if (addv == INSERT_VALUES) {
535: for (i=0; i<n; i++) yv[fslots[i]] = xv[first + i*step];
536: } else if (addv == ADD_VALUES) {
537: for (i=0; i<n; i++) yv[fslots[i]] += xv[first + i*step];
538: #if !defined(PETSC_USE_COMPLEX)
539: } else if (addv == MAX_VALUES) {
540: for (i=0; i<n; i++) yv[fslots[i]] = PetscMax(yv[fslots[i]],xv[first + i*step]);
541: #endif
542: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Wrong insert option");
543: }
544: VecRestoreArray(x,&xv);
545: if (x != y) {VecRestoreArray(y,&yv);}
546: return(0);
547: }
551: PetscErrorCode VecScatterView_SSToSG(VecScatter in,PetscViewer viewer)
552: {
553: PetscErrorCode ierr;
554: VecScatter_Seq_Stride *in_from = (VecScatter_Seq_Stride*)in->todata;
555: VecScatter_Seq_General *in_to = (VecScatter_Seq_General*)in->fromdata;
556: PetscInt i;
557: PetscBool isascii;
560: PetscObjectTypeCompare((PetscObject)in,PETSCVIEWERASCII,&isascii);
561: if (isascii) {
562: PetscViewerASCIIPrintf(viewer,"Sequential stride to general scatter\n");
563: for (i=0; i<in_to->n; i++) {
564: PetscViewerASCIIPrintf(viewer,"%D to %D\n",in_from->first + in_from->step*i,in_to->vslots[i]);
565: }
566: }
567: return(0);
568: }
569: /*
570: Scatter: sequential stride to sequential stride
571: */
574: PetscErrorCode VecScatterBegin_SSToSS(VecScatter ctx,Vec x,Vec y,InsertMode addv,ScatterMode mode)
575: {
576: VecScatter_Seq_Stride *gen_to = (VecScatter_Seq_Stride*)ctx->todata;
577: VecScatter_Seq_Stride *gen_from = (VecScatter_Seq_Stride*)ctx->fromdata;
578: PetscInt i,n = gen_from->n,to_first = gen_to->first,to_step = gen_to->step;
579: PetscErrorCode ierr;
580: PetscInt from_first = gen_from->first,from_step = gen_from->step;
581: PetscScalar *xv,*yv;
584: VecGetArrayRead(x,(const PetscScalar **)&xv);
585: if (x != y) {
586: VecGetArray(y,&yv);
587: } else yv = xv;
589: if (mode & SCATTER_REVERSE) {
590: from_first = gen_to->first;
591: to_first = gen_from->first;
592: from_step = gen_to->step;
593: to_step = gen_from->step;
594: }
596: if (addv == INSERT_VALUES) {
597: if (to_step == 1 && from_step == 1) {
598: PetscMemcpy(yv+to_first,xv+from_first,n*sizeof(PetscScalar));
599: } else {
600: for (i=0; i<n; i++) yv[to_first + i*to_step] = xv[from_first+i*from_step];
601: }
602: } else if (addv == ADD_VALUES) {
603: if (to_step == 1 && from_step == 1) {
604: yv += to_first; xv += from_first;
605: for (i=0; i<n; i++) yv[i] += xv[i];
606: } else {
607: for (i=0; i<n; i++) yv[to_first + i*to_step] += xv[from_first+i*from_step];
608: }
609: #if !defined(PETSC_USE_COMPLEX)
610: } else if (addv == MAX_VALUES) {
611: if (to_step == 1 && from_step == 1) {
612: yv += to_first; xv += from_first;
613: for (i=0; i<n; i++) yv[i] = PetscMax(yv[i],xv[i]);
614: } else {
615: for (i=0; i<n; i++) yv[to_first + i*to_step] = PetscMax(yv[to_first + i*to_step],xv[from_first+i*from_step]);
616: }
617: #endif
618: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Wrong insert option");
619: VecRestoreArrayRead(x,(const PetscScalar**)&xv);
620: if (x != y) {VecRestoreArray(y,&yv);}
621: return(0);
622: }
624: /* --------------------------------------------------------------------------*/
629: PetscErrorCode VecScatterCopy_SGToSG(VecScatter in,VecScatter out)
630: {
631: PetscErrorCode ierr;
632: VecScatter_Seq_General *in_to = (VecScatter_Seq_General*)in->todata,*out_to = NULL;
633: VecScatter_Seq_General *in_from = (VecScatter_Seq_General*)in->fromdata,*out_from = NULL;
636: out->begin = in->begin;
637: out->end = in->end;
638: out->copy = in->copy;
639: out->destroy = in->destroy;
640: out->view = in->view;
642: PetscMalloc2(1,VecScatter_Seq_General,&out_to,1,VecScatter_Seq_General,&out_from);
643: PetscMalloc2(in_to->n,PetscInt,&out_to->vslots,in_from->n,PetscInt,&out_from->vslots);
644: out_to->n = in_to->n;
645: out_to->type = in_to->type;
646: out_to->nonmatching_computed = PETSC_FALSE;
647: out_to->slots_nonmatching = 0;
648: out_to->is_copy = PETSC_FALSE;
649: PetscMemcpy(out_to->vslots,in_to->vslots,(out_to->n)*sizeof(PetscInt));
651: out_from->n = in_from->n;
652: out_from->type = in_from->type;
653: out_from->nonmatching_computed = PETSC_FALSE;
654: out_from->slots_nonmatching = 0;
655: out_from->is_copy = PETSC_FALSE;
656: PetscMemcpy(out_from->vslots,in_from->vslots,(out_from->n)*sizeof(PetscInt));
658: out->todata = (void*)out_to;
659: out->fromdata = (void*)out_from;
660: return(0);
661: }
665: PetscErrorCode VecScatterView_SGToSG(VecScatter in,PetscViewer viewer)
666: {
667: PetscErrorCode ierr;
668: VecScatter_Seq_General *in_to = (VecScatter_Seq_General*)in->todata;
669: VecScatter_Seq_General *in_from = (VecScatter_Seq_General*)in->fromdata;
670: PetscInt i;
671: PetscBool isascii;
674: PetscObjectTypeCompare((PetscObject)in,PETSCVIEWERASCII,&isascii);
675: if (isascii) {
676: PetscViewerASCIIPrintf(viewer,"Sequential general scatter\n");
677: for (i=0; i<in_to->n; i++) {
678: PetscViewerASCIIPrintf(viewer,"%D to %D\n",in_from->vslots[i],in_to->vslots[i]);
679: }
680: }
681: return(0);
682: }
687: PetscErrorCode VecScatterCopy_SGToSS(VecScatter in,VecScatter out)
688: {
689: PetscErrorCode ierr;
690: VecScatter_Seq_Stride *in_to = (VecScatter_Seq_Stride*)in->todata,*out_to = NULL;
691: VecScatter_Seq_General *in_from = (VecScatter_Seq_General*)in->fromdata,*out_from = NULL;
694: out->begin = in->begin;
695: out->end = in->end;
696: out->copy = in->copy;
697: out->destroy = in->destroy;
698: out->view = in->view;
700: PetscMalloc2(1,VecScatter_Seq_Stride,&out_to,1,VecScatter_Seq_General,&out_from);
701: PetscMalloc(in_from->n*sizeof(PetscInt),&out_from->vslots);
702: out_to->n = in_to->n;
703: out_to->type = in_to->type;
704: out_to->first = in_to->first;
705: out_to->step = in_to->step;
706: out_to->type = in_to->type;
708: out_from->n = in_from->n;
709: out_from->type = in_from->type;
710: out_from->nonmatching_computed = PETSC_FALSE;
711: out_from->slots_nonmatching = 0;
712: out_from->is_copy = PETSC_FALSE;
713: PetscMemcpy(out_from->vslots,in_from->vslots,(out_from->n)*sizeof(PetscInt));
715: out->todata = (void*)out_to;
716: out->fromdata = (void*)out_from;
717: return(0);
718: }
722: PetscErrorCode VecScatterView_SGToSS(VecScatter in,PetscViewer viewer)
723: {
724: PetscErrorCode ierr;
725: VecScatter_Seq_Stride *in_to = (VecScatter_Seq_Stride*)in->todata;
726: VecScatter_Seq_General *in_from = (VecScatter_Seq_General*)in->fromdata;
727: PetscInt i;
728: PetscBool isascii;
731: PetscObjectTypeCompare((PetscObject)in,PETSCVIEWERASCII,&isascii);
732: if (isascii) {
733: PetscViewerASCIIPrintf(viewer,"Sequential general scatter to stride\n");
734: for (i=0; i<in_to->n; i++) {
735: PetscViewerASCIIPrintf(viewer,"%D to %D\n",in_from->vslots[i],in_to->first + in_to->step*i);
736: }
737: }
738: return(0);
739: }
741: /* --------------------------------------------------------------------------*/
742: /*
743: Scatter: parallel to sequential vector, sequential strides for both.
744: */
747: PetscErrorCode VecScatterCopy_SSToSS(VecScatter in,VecScatter out)
748: {
749: VecScatter_Seq_Stride *in_to = (VecScatter_Seq_Stride*)in->todata,*out_to = NULL;
750: VecScatter_Seq_Stride *in_from = (VecScatter_Seq_Stride*)in->fromdata,*out_from = NULL;
751: PetscErrorCode ierr;
754: out->begin = in->begin;
755: out->end = in->end;
756: out->copy = in->copy;
757: out->destroy = in->destroy;
758: out->view = in->view;
760: PetscMalloc2(1,VecScatter_Seq_Stride,&out_to,1,VecScatter_Seq_Stride,&out_from);
761: out_to->n = in_to->n;
762: out_to->type = in_to->type;
763: out_to->first = in_to->first;
764: out_to->step = in_to->step;
765: out_to->type = in_to->type;
766: out_from->n = in_from->n;
767: out_from->type = in_from->type;
768: out_from->first = in_from->first;
769: out_from->step = in_from->step;
770: out_from->type = in_from->type;
771: out->todata = (void*)out_to;
772: out->fromdata = (void*)out_from;
773: return(0);
774: }
778: PetscErrorCode VecScatterView_SSToSS(VecScatter in,PetscViewer viewer)
779: {
780: VecScatter_Seq_Stride *in_to = (VecScatter_Seq_Stride*)in->todata;
781: VecScatter_Seq_Stride *in_from = (VecScatter_Seq_Stride*)in->fromdata;
782: PetscErrorCode ierr;
783: PetscBool isascii;
786: PetscObjectTypeCompare((PetscObject)in,PETSCVIEWERASCII,&isascii);
787: if (isascii) {
788: PetscViewerASCIIPrintf(viewer,"Sequential stride count %D start %D step to start %D stride %D\n",in_to->n,in_to->first,in_to->step,in_from->first,in_from->step);
789: }
790: return(0);
791: }
794: extern PetscErrorCode VecScatterCreate_PtoS(PetscInt,const PetscInt*,PetscInt,const PetscInt*,Vec,Vec,PetscInt,VecScatter);
795: extern PetscErrorCode VecScatterCreate_PtoP(PetscInt,const PetscInt*,PetscInt,const PetscInt*,Vec,Vec,PetscInt,VecScatter);
796: extern PetscErrorCode VecScatterCreate_StoP(PetscInt,const PetscInt*,PetscInt,const PetscInt*,Vec,Vec,PetscInt,VecScatter);
798: /* =======================================================================*/
799: #define VEC_SEQ_ID 0
800: #define VEC_MPI_ID 1
801: #define IS_GENERAL_ID 0
802: #define IS_STRIDE_ID 1
803: #define IS_BLOCK_ID 2
805: /*
806: Blocksizes we have optimized scatters for
807: */
809: #define VecScatterOptimizedBS(mbs) ((2 <= mbs && mbs <= 8) || mbs == 12)
811: PetscErrorCode VecScatterCreateEmpty(MPI_Comm comm,VecScatter *newctx)
812: {
813: VecScatter ctx;
817: PetscHeaderCreate(ctx,_p_VecScatter,int,VEC_SCATTER_CLASSID,"VecScatter","VecScatter","Vec",comm,VecScatterDestroy,VecScatterView);
818: ctx->inuse = PETSC_FALSE;
819: ctx->beginandendtogether = PETSC_FALSE;
820: PetscOptionsGetBool(NULL,"-vecscatter_merge",&ctx->beginandendtogether,NULL);
821: if (ctx->beginandendtogether) {
822: PetscInfo(ctx,"Using combined (merged) vector scatter begin and end\n");
823: }
824: ctx->packtogether = PETSC_FALSE;
825: PetscOptionsGetBool(NULL,"-vecscatter_packtogether",&ctx->packtogether,NULL);
826: if (ctx->packtogether) {
827: PetscInfo(ctx,"Pack all messages before sending\n");
828: }
829: *newctx = ctx;
830: return(0);
831: }
835: /*@C
836: VecScatterCreate - Creates a vector scatter context.
838: Collective on Vec
840: Input Parameters:
841: + xin - a vector that defines the shape (parallel data layout of the vector)
842: of vectors from which we scatter
843: . yin - a vector that defines the shape (parallel data layout of the vector)
844: of vectors to which we scatter
845: . ix - the indices of xin to scatter (if NULL scatters all values)
846: - iy - the indices of yin to hold results (if NULL fills entire vector yin)
848: Output Parameter:
849: . newctx - location to store the new scatter context
851: Options Database Keys: (uses regular MPI_Sends by default)
852: + -vecscatter_view - Prints detail of communications
853: . -vecscatter_view ::ascii_info - Print less details about communication
854: . -vecscatter_ssend - Uses MPI_Ssend_init() instead of MPI_Send_init()
855: . -vecscatter_rsend - use ready receiver mode for MPI sends
856: . -vecscatter_merge - VecScatterBegin() handles all of the communication, VecScatterEnd() is a nop
857: eliminates the chance for overlap of computation and communication
858: . -vecscatter_sendfirst - Posts sends before receives
859: . -vecscatter_packtogether - Pack all messages before sending, receive all messages before unpacking
860: . -vecscatter_alltoall - Uses MPI all to all communication for scatter
861: . -vecscatter_window - Use MPI 2 window operations to move data
862: . -vecscatter_nopack - Avoid packing to work vector when possible (if used with -vecscatter_alltoall then will use MPI_Alltoallw()
863: - -vecscatter_reproduce - insure that the order of the communications are done the same for each scatter, this under certain circumstances
864: will make the results of scatters deterministic when otherwise they are not (it may be slower also).
866: $
867: $ --When packing is used--
868: $ MPI Datatypes (no packing) sendfirst merge packtogether persistent*
869: $ _nopack _sendfirst _merge _packtogether -vecscatter_
870: $ ----------------------------------------------------------------------------------------------------------------------------
871: $ Message passing Send p X X X always
872: $ Ssend p X X X always _ssend
873: $ Rsend p nonsense X X always _rsend
874: $ AlltoAll v or w X nonsense always X nonsense _alltoall
875: $ MPI_Win p nonsense p p nonsense _window
876: $
877: $ Since persistent sends and receives require a constant memory address they can only be used when data is packed into the work vector
878: $ because the in and out array may be different for each call to VecScatterBegin/End().
879: $
880: $ p indicates possible, but not implemented. X indicates implemented
881: $
883: Level: intermediate
885: Notes:
886: In calls to VecScatter() you can use different vectors than the xin and
887: yin you used above; BUT they must have the same parallel data layout, for example,
888: they could be obtained from VecDuplicate().
889: A VecScatter context CANNOT be used in two or more simultaneous scatters;
890: that is you cannot call a second VecScatterBegin() with the same scatter
891: context until the VecScatterEnd() has been called on the first VecScatterBegin().
892: In this case a separate VecScatter is needed for each concurrent scatter.
894: Currently the MPI_Send(), MPI_Ssend() and MPI_Rsend() all use PERSISTENT versions.
895: (this unfortunately requires that the same in and out arrays be used for each use, this
896: is why when not using MPI_alltoallw() we always need to pack the input into the work array before sending
897: and unpack upon receeving instead of using MPI datatypes to avoid the packing/unpacking).
899: Both ix and iy cannot be NULL at the same time.
901: Concepts: scatter^between vectors
902: Concepts: gather^between vectors
904: .seealso: VecScatterDestroy(), VecScatterCreateToAll(), VecScatterCreateToZero()
905: @*/
906: PetscErrorCode VecScatterCreate(Vec xin,IS ix,Vec yin,IS iy,VecScatter *newctx)
907: {
908: VecScatter ctx;
909: PetscErrorCode ierr;
910: PetscMPIInt size;
911: PetscInt xin_type = VEC_SEQ_ID,yin_type = VEC_SEQ_ID,*range;
912: PetscInt ix_type = IS_GENERAL_ID,iy_type = IS_GENERAL_ID;
913: MPI_Comm comm,ycomm;
914: PetscBool totalv,ixblock,iyblock,iystride,islocal,cando,flag;
915: IS tix = 0,tiy = 0;
916: PetscViewer viewer;
917: PetscViewerFormat format;
920: if (!ix && !iy) SETERRQ(PetscObjectComm((PetscObject)xin),PETSC_ERR_SUP,"Cannot pass default in for both input and output indices");
922: /*
923: Determine if the vectors are "parallel", ie. it shares a comm with other processors, or
924: sequential (it does not share a comm). The difference is that parallel vectors treat the
925: index set as providing indices in the global parallel numbering of the vector, with
926: sequential vectors treat the index set as providing indices in the local sequential
927: numbering
928: */
929: PetscObjectGetComm((PetscObject)xin,&comm);
930: MPI_Comm_size(comm,&size);
931: if (size > 1) xin_type = VEC_MPI_ID;
933: PetscObjectGetComm((PetscObject)yin,&ycomm);
934: MPI_Comm_size(ycomm,&size);
935: if (size > 1) {comm = ycomm; yin_type = VEC_MPI_ID;}
937: /* generate the Scatter context */
938: PetscHeaderCreate(ctx,_p_VecScatter,int,VEC_SCATTER_CLASSID,"VecScatter","VecScatter","Vec",comm,VecScatterDestroy,VecScatterView);
939: ctx->inuse = PETSC_FALSE;
941: ctx->beginandendtogether = PETSC_FALSE;
942: PetscOptionsGetBool(NULL,"-vecscatter_merge",&ctx->beginandendtogether,NULL);
943: if (ctx->beginandendtogether) {
944: PetscInfo(ctx,"Using combined (merged) vector scatter begin and end\n");
945: }
946: ctx->packtogether = PETSC_FALSE;
947: PetscOptionsGetBool(NULL,"-vecscatter_packtogether",&ctx->packtogether,NULL);
948: if (ctx->packtogether) {
949: PetscInfo(ctx,"Pack all messages before sending\n");
950: }
952: VecGetLocalSize(xin,&ctx->from_n);
953: VecGetLocalSize(yin,&ctx->to_n);
955: /*
956: if ix or iy is not included; assume just grabbing entire vector
957: */
958: if (!ix && xin_type == VEC_SEQ_ID) {
959: ISCreateStride(comm,ctx->from_n,0,1,&ix);
960: tix = ix;
961: } else if (!ix && xin_type == VEC_MPI_ID) {
962: if (yin_type == VEC_MPI_ID) {
963: PetscInt ntmp, low;
964: VecGetLocalSize(xin,&ntmp);
965: VecGetOwnershipRange(xin,&low,NULL);
966: ISCreateStride(comm,ntmp,low,1,&ix);
967: } else {
968: PetscInt Ntmp;
969: VecGetSize(xin,&Ntmp);
970: ISCreateStride(comm,Ntmp,0,1,&ix);
971: }
972: tix = ix;
973: } else if (!ix) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"ix not given, but not Seq or MPI vector");
975: if (!iy && yin_type == VEC_SEQ_ID) {
976: ISCreateStride(comm,ctx->to_n,0,1,&iy);
977: tiy = iy;
978: } else if (!iy && yin_type == VEC_MPI_ID) {
979: if (xin_type == VEC_MPI_ID) {
980: PetscInt ntmp, low;
981: VecGetLocalSize(yin,&ntmp);
982: VecGetOwnershipRange(yin,&low,NULL);
983: ISCreateStride(comm,ntmp,low,1,&iy);
984: } else {
985: PetscInt Ntmp;
986: VecGetSize(yin,&Ntmp);
987: ISCreateStride(comm,Ntmp,0,1,&iy);
988: }
989: tiy = iy;
990: } else if (!iy) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"iy not given, but not Seq or MPI vector");
992: /*
993: Determine types of index sets
994: */
995: PetscObjectTypeCompare((PetscObject)ix,ISBLOCK,&flag);
996: if (flag) ix_type = IS_BLOCK_ID;
997: PetscObjectTypeCompare((PetscObject)iy,ISBLOCK,&flag);
998: if (flag) iy_type = IS_BLOCK_ID;
999: PetscObjectTypeCompare((PetscObject)ix,ISSTRIDE,&flag);
1000: if (flag) ix_type = IS_STRIDE_ID;
1001: PetscObjectTypeCompare((PetscObject)iy,ISSTRIDE,&flag);
1002: if (flag) iy_type = IS_STRIDE_ID;
1004: /* ===========================================================================================================
1005: Check for special cases
1006: ==========================================================================================================*/
1007: /* ---------------------------------------------------------------------------*/
1008: if (xin_type == VEC_SEQ_ID && yin_type == VEC_SEQ_ID) {
1009: if (ix_type == IS_GENERAL_ID && iy_type == IS_GENERAL_ID) {
1010: PetscInt nx,ny;
1011: const PetscInt *idx,*idy;
1012: VecScatter_Seq_General *to = NULL,*from = NULL;
1014: ISGetLocalSize(ix,&nx);
1015: ISGetLocalSize(iy,&ny);
1016: if (nx != ny) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Local scatter sizes don't match");
1017: ISGetIndices(ix,&idx);
1018: ISGetIndices(iy,&idy);
1019: PetscMalloc2(1,VecScatter_Seq_General,&to,1,VecScatter_Seq_General,&from);
1020: PetscMalloc2(nx,PetscInt,&to->vslots,nx,PetscInt,&from->vslots);
1021: to->n = nx;
1022: #if defined(PETSC_USE_DEBUG)
1023: VecScatterCheckIndices_Private(ctx->to_n,ny,idy);
1024: #endif
1025: PetscMemcpy(to->vslots,idy,nx*sizeof(PetscInt));
1026: from->n = nx;
1027: #if defined(PETSC_USE_DEBUG)
1028: VecScatterCheckIndices_Private(ctx->from_n,nx,idx);
1029: #endif
1030: PetscMemcpy(from->vslots,idx,nx*sizeof(PetscInt));
1031: to->type = VEC_SCATTER_SEQ_GENERAL;
1032: from->type = VEC_SCATTER_SEQ_GENERAL;
1033: ctx->todata = (void*)to;
1034: ctx->fromdata = (void*)from;
1035: ctx->begin = VecScatterBegin_SGToSG;
1036: ctx->end = 0;
1037: ctx->destroy = VecScatterDestroy_SGToSG;
1038: ctx->copy = VecScatterCopy_SGToSG;
1039: ctx->view = VecScatterView_SGToSG;
1040: PetscInfo(xin,"Special case: sequential vector general scatter\n");
1041: goto functionend;
1042: } else if (ix_type == IS_STRIDE_ID && iy_type == IS_STRIDE_ID) {
1043: PetscInt nx,ny,to_first,to_step,from_first,from_step;
1044: VecScatter_Seq_Stride *from8 = NULL,*to8 = NULL;
1046: ISGetLocalSize(ix,&nx);
1047: ISGetLocalSize(iy,&ny);
1048: if (nx != ny) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Local scatter sizes don't match");
1049: ISStrideGetInfo(iy,&to_first,&to_step);
1050: ISStrideGetInfo(ix,&from_first,&from_step);
1051: PetscMalloc2(1,VecScatter_Seq_Stride,&to8,1,VecScatter_Seq_Stride,&from8);
1052: to8->n = nx;
1053: to8->first = to_first;
1054: to8->step = to_step;
1055: from8->n = nx;
1056: from8->first = from_first;
1057: from8->step = from_step;
1058: to8->type = VEC_SCATTER_SEQ_STRIDE;
1059: from8->type = VEC_SCATTER_SEQ_STRIDE;
1060: ctx->todata = (void*)to8;
1061: ctx->fromdata = (void*)from8;
1062: ctx->begin = VecScatterBegin_SSToSS;
1063: ctx->end = 0;
1064: ctx->destroy = VecScatterDestroy_SSToSS;
1065: ctx->copy = VecScatterCopy_SSToSS;
1066: ctx->view = VecScatterView_SSToSS;
1067: PetscInfo(xin,"Special case: sequential vector stride to stride\n");
1068: goto functionend;
1069: } else if (ix_type == IS_GENERAL_ID && iy_type == IS_STRIDE_ID) {
1070: PetscInt nx,ny,first,step;
1071: const PetscInt *idx;
1072: VecScatter_Seq_General *from9 = NULL;
1073: VecScatter_Seq_Stride *to9 = NULL;
1075: ISGetLocalSize(ix,&nx);
1076: ISGetIndices(ix,&idx);
1077: ISGetLocalSize(iy,&ny);
1078: ISStrideGetInfo(iy,&first,&step);
1079: if (nx != ny) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Local scatter sizes don't match");
1080: PetscMalloc2(1,VecScatter_Seq_Stride,&to9,1,VecScatter_Seq_General,&from9);
1081: PetscMalloc(nx*sizeof(PetscInt),&from9->vslots);
1082: to9->n = nx;
1083: to9->first = first;
1084: to9->step = step;
1085: from9->n = nx;
1086: #if defined(PETSC_USE_DEBUG)
1087: VecScatterCheckIndices_Private(ctx->from_n,nx,idx);
1088: #endif
1089: PetscMemcpy(from9->vslots,idx,nx*sizeof(PetscInt));
1090: ctx->todata = (void*)to9; ctx->fromdata = (void*)from9;
1091: if (step == 1) ctx->begin = VecScatterBegin_SGToSS_Stride1;
1092: else ctx->begin = VecScatterBegin_SGToSS;
1093: ctx->destroy = VecScatterDestroy_SGToSS;
1094: ctx->end = 0;
1095: ctx->copy = VecScatterCopy_SGToSS;
1096: ctx->view = VecScatterView_SGToSS;
1097: to9->type = VEC_SCATTER_SEQ_STRIDE;
1098: from9->type = VEC_SCATTER_SEQ_GENERAL;
1099: PetscInfo(xin,"Special case: sequential vector general to stride\n");
1100: goto functionend;
1101: } else if (ix_type == IS_STRIDE_ID && iy_type == IS_GENERAL_ID) {
1102: PetscInt nx,ny,first,step;
1103: const PetscInt *idy;
1104: VecScatter_Seq_General *to10 = NULL;
1105: VecScatter_Seq_Stride *from10 = NULL;
1107: ISGetLocalSize(ix,&nx);
1108: ISGetIndices(iy,&idy);
1109: ISGetLocalSize(iy,&ny);
1110: ISStrideGetInfo(ix,&first,&step);
1111: if (nx != ny) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Local scatter sizes don't match");
1112: PetscMalloc2(1,VecScatter_Seq_General,&to10,1,VecScatter_Seq_Stride,&from10);
1113: PetscMalloc(nx*sizeof(PetscInt),&to10->vslots);
1114: from10->n = nx;
1115: from10->first = first;
1116: from10->step = step;
1117: to10->n = nx;
1118: #if defined(PETSC_USE_DEBUG)
1119: VecScatterCheckIndices_Private(ctx->to_n,ny,idy);
1120: #endif
1121: PetscMemcpy(to10->vslots,idy,nx*sizeof(PetscInt));
1122: ctx->todata = (void*)to10;
1123: ctx->fromdata = (void*)from10;
1124: if (step == 1) ctx->begin = VecScatterBegin_SSToSG_Stride1;
1125: else ctx->begin = VecScatterBegin_SSToSG;
1126: ctx->destroy = VecScatterDestroy_SSToSG;
1127: ctx->end = 0;
1128: ctx->copy = 0;
1129: ctx->view = VecScatterView_SSToSG;
1130: to10->type = VEC_SCATTER_SEQ_GENERAL;
1131: from10->type = VEC_SCATTER_SEQ_STRIDE;
1132: PetscInfo(xin,"Special case: sequential vector stride to general\n");
1133: goto functionend;
1134: } else {
1135: PetscInt nx,ny;
1136: const PetscInt *idx,*idy;
1137: VecScatter_Seq_General *to11 = NULL,*from11 = NULL;
1138: PetscBool idnx,idny;
1140: ISGetLocalSize(ix,&nx);
1141: ISGetLocalSize(iy,&ny);
1142: if (nx != ny) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Local scatter sizes don't match, in %D out %D",nx,ny);
1144: ISIdentity(ix,&idnx);
1145: ISIdentity(iy,&idny);
1146: if (idnx && idny) {
1147: VecScatter_Seq_Stride *to13 = NULL,*from13 = NULL;
1148: PetscMalloc2(1,VecScatter_Seq_Stride,&to13,1,VecScatter_Seq_Stride,&from13);
1149: to13->n = nx;
1150: to13->first = 0;
1151: to13->step = 1;
1152: from13->n = nx;
1153: from13->first = 0;
1154: from13->step = 1;
1155: to13->type = VEC_SCATTER_SEQ_STRIDE;
1156: from13->type = VEC_SCATTER_SEQ_STRIDE;
1157: ctx->todata = (void*)to13;
1158: ctx->fromdata = (void*)from13;
1159: ctx->begin = VecScatterBegin_SSToSS;
1160: ctx->end = 0;
1161: ctx->destroy = VecScatterDestroy_SSToSS;
1162: ctx->copy = VecScatterCopy_SSToSS;
1163: ctx->view = VecScatterView_SSToSS;
1164: PetscInfo(xin,"Special case: sequential copy\n");
1165: goto functionend;
1166: }
1168: ISGetIndices(iy,&idy);
1169: ISGetIndices(ix,&idx);
1170: PetscMalloc2(1,VecScatter_Seq_General,&to11,1,VecScatter_Seq_General,&from11);
1171: PetscMalloc2(nx,PetscInt,&to11->vslots,nx,PetscInt,&from11->vslots);
1172: to11->n = nx;
1173: #if defined(PETSC_USE_DEBUG)
1174: VecScatterCheckIndices_Private(ctx->to_n,ny,idy);
1175: #endif
1176: PetscMemcpy(to11->vslots,idy,nx*sizeof(PetscInt));
1177: from11->n = nx;
1178: #if defined(PETSC_USE_DEBUG)
1179: VecScatterCheckIndices_Private(ctx->from_n,nx,idx);
1180: #endif
1181: PetscMemcpy(from11->vslots,idx,nx*sizeof(PetscInt));
1182: to11->type = VEC_SCATTER_SEQ_GENERAL;
1183: from11->type = VEC_SCATTER_SEQ_GENERAL;
1184: ctx->todata = (void*)to11;
1185: ctx->fromdata = (void*)from11;
1186: ctx->begin = VecScatterBegin_SGToSG;
1187: ctx->end = 0;
1188: ctx->destroy = VecScatterDestroy_SGToSG;
1189: ctx->copy = VecScatterCopy_SGToSG;
1190: ctx->view = VecScatterView_SGToSG;
1191: ISRestoreIndices(ix,&idx);
1192: ISRestoreIndices(iy,&idy);
1193: PetscInfo(xin,"Sequential vector scatter with block indices\n");
1194: goto functionend;
1195: }
1196: }
1197: /* ---------------------------------------------------------------------------*/
1198: if (xin_type == VEC_MPI_ID && yin_type == VEC_SEQ_ID) {
1200: /* ===========================================================================================================
1201: Check for special cases
1202: ==========================================================================================================*/
1203: islocal = PETSC_FALSE;
1204: /* special case extracting (subset of) local portion */
1205: if (ix_type == IS_STRIDE_ID && iy_type == IS_STRIDE_ID) {
1206: PetscInt nx,ny,to_first,to_step,from_first,from_step;
1207: PetscInt start,end,min,max;
1208: VecScatter_Seq_Stride *from12 = NULL,*to12 = NULL;
1210: VecGetOwnershipRange(xin,&start,&end);
1211: ISGetLocalSize(ix,&nx);
1212: ISStrideGetInfo(ix,&from_first,&from_step);
1213: ISGetLocalSize(iy,&ny);
1214: ISStrideGetInfo(iy,&to_first,&to_step);
1215: if (nx != ny) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Local scatter sizes don't match");
1216: ISGetMinMax(ix,&min,&max);
1217: if (min >= start && max < end) islocal = PETSC_TRUE;
1218: else islocal = PETSC_FALSE;
1219: MPI_Allreduce(&islocal,&cando,1,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)xin));
1220: if (cando) {
1221: PetscMalloc2(1,VecScatter_Seq_Stride,&to12,1,VecScatter_Seq_Stride,&from12);
1222: to12->n = nx;
1223: to12->first = to_first;
1224: to12->step = to_step;
1225: from12->n = nx;
1226: from12->first = from_first-start;
1227: from12->step = from_step;
1228: to12->type = VEC_SCATTER_SEQ_STRIDE;
1229: from12->type = VEC_SCATTER_SEQ_STRIDE;
1230: ctx->todata = (void*)to12;
1231: ctx->fromdata = (void*)from12;
1232: ctx->begin = VecScatterBegin_SSToSS;
1233: ctx->end = 0;
1234: ctx->destroy = VecScatterDestroy_SSToSS;
1235: ctx->copy = VecScatterCopy_SSToSS;
1236: ctx->view = VecScatterView_SSToSS;
1237: PetscInfo(xin,"Special case: processors only getting local values\n");
1238: goto functionend;
1239: }
1240: } else {
1241: MPI_Allreduce(&islocal,&cando,1,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)xin));
1242: }
1244: /* test for special case of all processors getting entire vector */
1245: /* contains check that PetscMPIInt can handle the sizes needed */
1246: totalv = PETSC_FALSE;
1247: if (ix_type == IS_STRIDE_ID && iy_type == IS_STRIDE_ID) {
1248: PetscInt i,nx,ny,to_first,to_step,from_first,from_step,N;
1249: PetscMPIInt *count = NULL,*displx;
1250: VecScatter_MPI_ToAll *sto = NULL;
1252: ISGetLocalSize(ix,&nx);
1253: ISStrideGetInfo(ix,&from_first,&from_step);
1254: ISGetLocalSize(iy,&ny);
1255: ISStrideGetInfo(iy,&to_first,&to_step);
1256: if (nx != ny) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Local scatter sizes don't match");
1257: VecGetSize(xin,&N);
1258: if (nx != N) totalv = PETSC_FALSE;
1259: else if (from_first == 0 && from_step == 1 && from_first == to_first && from_step == to_step) totalv = PETSC_TRUE;
1260: else totalv = PETSC_FALSE;
1261: MPI_Allreduce(&totalv,&cando,1,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)xin));
1263: #if defined(PETSC_USE_64BIT_INDICES)
1264: if (cando && (yin->map->N < PETSC_MPI_INT_MAX)) {
1265: #else
1266: if (cando) {
1267: #endif
1268: MPI_Comm_size(PetscObjectComm((PetscObject)ctx),&size);
1269: PetscMalloc3(1,VecScatter_MPI_ToAll,&sto,size,PetscMPIInt,&count,size,PetscMPIInt,&displx);
1270: range = xin->map->range;
1271: for (i=0; i<size; i++) {
1272: PetscMPIIntCast(range[i+1] - range[i],count+i);
1273: PetscMPIIntCast(range[i],displx+i);
1274: }
1275: sto->count = count;
1276: sto->displx = displx;
1277: sto->work1 = 0;
1278: sto->work2 = 0;
1279: sto->type = VEC_SCATTER_MPI_TOALL;
1280: ctx->todata = (void*)sto;
1281: ctx->fromdata = 0;
1282: ctx->begin = VecScatterBegin_MPI_ToAll;
1283: ctx->end = 0;
1284: ctx->destroy = VecScatterDestroy_MPI_ToAll;
1285: ctx->copy = VecScatterCopy_MPI_ToAll;
1286: ctx->view = VecScatterView_MPI_ToAll;
1287: PetscInfo(xin,"Special case: all processors get entire parallel vector\n");
1288: goto functionend;
1289: }
1290: } else {
1291: MPI_Allreduce(&totalv,&cando,1,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)xin));
1292: }
1294: /* test for special case of processor 0 getting entire vector */
1295: /* contains check that PetscMPIInt can handle the sizes needed */
1296: totalv = PETSC_FALSE;
1297: if (ix_type == IS_STRIDE_ID && iy_type == IS_STRIDE_ID) {
1298: PetscInt i,nx,ny,to_first,to_step,from_first,from_step,N;
1299: PetscMPIInt rank,*count = NULL,*displx;
1300: VecScatter_MPI_ToAll *sto = NULL;
1302: PetscObjectGetComm((PetscObject)xin,&comm);
1303: MPI_Comm_rank(comm,&rank);
1304: ISGetLocalSize(ix,&nx);
1305: ISStrideGetInfo(ix,&from_first,&from_step);
1306: ISGetLocalSize(iy,&ny);
1307: ISStrideGetInfo(iy,&to_first,&to_step);
1308: if (nx != ny) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Local scatter sizes don't match");
1309: if (!rank) {
1310: VecGetSize(xin,&N);
1311: if (nx != N) totalv = PETSC_FALSE;
1312: else if (from_first == 0 && from_step == 1 &&
1313: from_first == to_first && from_step == to_step) totalv = PETSC_TRUE;
1314: else totalv = PETSC_FALSE;
1315: } else {
1316: if (!nx) totalv = PETSC_TRUE;
1317: else totalv = PETSC_FALSE;
1318: }
1319: MPI_Allreduce(&totalv,&cando,1,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)xin));
1321: #if defined(PETSC_USE_64BIT_INDICES)
1322: if (cando && (yin->map->N < PETSC_MPI_INT_MAX)) {
1323: #else
1324: if (cando) {
1325: #endif
1326: MPI_Comm_size(PetscObjectComm((PetscObject)ctx),&size);
1327: PetscMalloc3(1,VecScatter_MPI_ToAll,&sto,size,PetscMPIInt,&count,size,PetscMPIInt,&displx);
1328: range = xin->map->range;
1329: for (i=0; i<size; i++) {
1330: PetscMPIIntCast(range[i+1] - range[i],count+i);
1331: PetscMPIIntCast(range[i],displx+i);
1332: }
1333: sto->count = count;
1334: sto->displx = displx;
1335: sto->work1 = 0;
1336: sto->work2 = 0;
1337: sto->type = VEC_SCATTER_MPI_TOONE;
1338: ctx->todata = (void*)sto;
1339: ctx->fromdata = 0;
1340: ctx->begin = VecScatterBegin_MPI_ToOne;
1341: ctx->end = 0;
1342: ctx->destroy = VecScatterDestroy_MPI_ToAll;
1343: ctx->copy = VecScatterCopy_MPI_ToAll;
1344: ctx->view = VecScatterView_MPI_ToAll;
1345: PetscInfo(xin,"Special case: processor zero gets entire parallel vector, rest get none\n");
1346: goto functionend;
1347: }
1348: } else {
1349: MPI_Allreduce(&totalv,&cando,1,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)xin));
1350: }
1352: PetscObjectTypeCompare((PetscObject)ix,ISBLOCK,&ixblock);
1353: PetscObjectTypeCompare((PetscObject)iy,ISBLOCK,&iyblock);
1354: PetscObjectTypeCompare((PetscObject)iy,ISSTRIDE,&iystride);
1355: if (ixblock) {
1356: /* special case block to block */
1357: if (iyblock) {
1358: PetscInt nx,ny,bsx,bsy;
1359: const PetscInt *idx,*idy;
1360: ISGetBlockSize(iy,&bsy);
1361: ISGetBlockSize(ix,&bsx);
1362: if (bsx == bsy && VecScatterOptimizedBS(bsx)) {
1363: ISBlockGetLocalSize(ix,&nx);
1364: ISBlockGetIndices(ix,&idx);
1365: ISBlockGetLocalSize(iy,&ny);
1366: ISBlockGetIndices(iy,&idy);
1367: if (nx != ny) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Local scatter sizes don't match");
1368: VecScatterCreate_PtoS(nx,idx,ny,idy,xin,yin,bsx,ctx);
1369: ISBlockRestoreIndices(ix,&idx);
1370: ISBlockRestoreIndices(iy,&idy);
1371: PetscInfo(xin,"Special case: blocked indices\n");
1372: goto functionend;
1373: }
1374: /* special case block to stride */
1375: } else if (iystride) {
1376: PetscInt ystart,ystride,ysize,bsx;
1377: ISStrideGetInfo(iy,&ystart,&ystride);
1378: ISGetLocalSize(iy,&ysize);
1379: ISGetBlockSize(ix,&bsx);
1380: /* see if stride index set is equivalent to block index set */
1381: if (VecScatterOptimizedBS(bsx) && ((ystart % bsx) == 0) && (ystride == 1) && ((ysize % bsx) == 0)) {
1382: PetscInt nx,il,*idy;
1383: const PetscInt *idx;
1384: ISBlockGetLocalSize(ix,&nx);
1385: ISBlockGetIndices(ix,&idx);
1386: if (ysize != bsx*nx) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Local scatter sizes don't match");
1387: PetscMalloc(nx*sizeof(PetscInt),&idy);
1388: if (nx) {
1389: idy[0] = ystart/bsx;
1390: for (il=1; il<nx; il++) idy[il] = idy[il-1] + 1;
1391: }
1392: VecScatterCreate_PtoS(nx,idx,nx,idy,xin,yin,bsx,ctx);
1393: PetscFree(idy);
1394: ISBlockRestoreIndices(ix,&idx);
1395: PetscInfo(xin,"Special case: blocked indices to stride\n");
1396: goto functionend;
1397: }
1398: }
1399: }
1400: /* left over general case */
1401: {
1402: PetscInt nx,ny;
1403: const PetscInt *idx,*idy;
1404: ISGetLocalSize(ix,&nx);
1405: ISGetIndices(ix,&idx);
1406: ISGetLocalSize(iy,&ny);
1407: ISGetIndices(iy,&idy);
1408: if (nx != ny) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Local scatter sizes don't match (%d %d)",nx,ny);
1409: VecScatterCreate_PtoS(nx,idx,ny,idy,xin,yin,1,ctx);
1410: ISRestoreIndices(ix,&idx);
1411: ISRestoreIndices(iy,&idy);
1412: PetscInfo(xin,"General case: MPI to Seq\n");
1413: goto functionend;
1414: }
1415: }
1416: /* ---------------------------------------------------------------------------*/
1417: if (xin_type == VEC_SEQ_ID && yin_type == VEC_MPI_ID) {
1418: /* ===========================================================================================================
1419: Check for special cases
1420: ==========================================================================================================*/
1421: /* special case local copy portion */
1422: islocal = PETSC_FALSE;
1423: if (ix_type == IS_STRIDE_ID && iy_type == IS_STRIDE_ID) {
1424: PetscInt nx,ny,to_first,to_step,from_step,start,end,from_first,min,max;
1425: VecScatter_Seq_Stride *from = NULL,*to = NULL;
1427: VecGetOwnershipRange(yin,&start,&end);
1428: ISGetLocalSize(ix,&nx);
1429: ISStrideGetInfo(ix,&from_first,&from_step);
1430: ISGetLocalSize(iy,&ny);
1431: ISStrideGetInfo(iy,&to_first,&to_step);
1432: if (nx != ny) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Local scatter sizes don't match");
1433: ISGetMinMax(iy,&min,&max);
1434: if (min >= start && max < end) islocal = PETSC_TRUE;
1435: else islocal = PETSC_FALSE;
1436: MPI_Allreduce(&islocal,&cando,1,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)yin));
1437: if (cando) {
1438: PetscMalloc2(1,VecScatter_Seq_Stride,&to,1,VecScatter_Seq_Stride,&from);
1439: to->n = nx;
1440: to->first = to_first-start;
1441: to->step = to_step;
1442: from->n = nx;
1443: from->first = from_first;
1444: from->step = from_step;
1445: to->type = VEC_SCATTER_SEQ_STRIDE;
1446: from->type = VEC_SCATTER_SEQ_STRIDE;
1447: ctx->todata = (void*)to;
1448: ctx->fromdata = (void*)from;
1449: ctx->begin = VecScatterBegin_SSToSS;
1450: ctx->end = 0;
1451: ctx->destroy = VecScatterDestroy_SSToSS;
1452: ctx->copy = VecScatterCopy_SSToSS;
1453: ctx->view = VecScatterView_SSToSS;
1454: PetscInfo(xin,"Special case: sequential stride to MPI stride\n");
1455: goto functionend;
1456: }
1457: } else {
1458: MPI_Allreduce(&islocal,&cando,1,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)yin));
1459: }
1460: /* special case block to stride */
1461: if (ix_type == IS_BLOCK_ID && iy_type == IS_STRIDE_ID) {
1462: PetscInt ystart,ystride,ysize,bsx;
1463: ISStrideGetInfo(iy,&ystart,&ystride);
1464: ISGetLocalSize(iy,&ysize);
1465: ISGetBlockSize(ix,&bsx);
1466: /* see if stride index set is equivalent to block index set */
1467: if (VecScatterOptimizedBS(bsx) && ((ystart % bsx) == 0) && (ystride == 1) && ((ysize % bsx) == 0)) {
1468: PetscInt nx,il,*idy;
1469: const PetscInt *idx;
1470: ISBlockGetLocalSize(ix,&nx);
1471: ISBlockGetIndices(ix,&idx);
1472: if (ysize != bsx*nx) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Local scatter sizes don't match");
1473: PetscMalloc(nx*sizeof(PetscInt),&idy);
1474: if (nx) {
1475: idy[0] = ystart/bsx;
1476: for (il=1; il<nx; il++) idy[il] = idy[il-1] + 1;
1477: }
1478: VecScatterCreate_StoP(nx,idx,nx,idy,xin,yin,bsx,ctx);
1479: PetscFree(idy);
1480: ISBlockRestoreIndices(ix,&idx);
1481: PetscInfo(xin,"Special case: Blocked indices to stride\n");
1482: goto functionend;
1483: }
1484: }
1486: /* general case */
1487: {
1488: PetscInt nx,ny;
1489: const PetscInt *idx,*idy;
1490: ISGetLocalSize(ix,&nx);
1491: ISGetIndices(ix,&idx);
1492: ISGetLocalSize(iy,&ny);
1493: ISGetIndices(iy,&idy);
1494: if (nx != ny) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Local scatter sizes don't match");
1495: VecScatterCreate_StoP(nx,idx,ny,idy,xin,yin,1,ctx);
1496: ISRestoreIndices(ix,&idx);
1497: ISRestoreIndices(iy,&idy);
1498: PetscInfo(xin,"General case: Seq to MPI\n");
1499: goto functionend;
1500: }
1501: }
1502: /* ---------------------------------------------------------------------------*/
1503: if (xin_type == VEC_MPI_ID && yin_type == VEC_MPI_ID) {
1504: /* no special cases for now */
1505: PetscInt nx,ny;
1506: const PetscInt *idx,*idy;
1507: ISGetLocalSize(ix,&nx);
1508: ISGetIndices(ix,&idx);
1509: ISGetLocalSize(iy,&ny);
1510: ISGetIndices(iy,&idy);
1511: if (nx != ny) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Local scatter sizes don't match");
1512: VecScatterCreate_PtoP(nx,idx,ny,idy,xin,yin,1,ctx);
1513: ISRestoreIndices(ix,&idx);
1514: ISRestoreIndices(iy,&idy);
1515: PetscInfo(xin,"General case: MPI to MPI\n");
1516: goto functionend;
1517: }
1519: functionend:
1520: *newctx = ctx;
1521: ISDestroy(&tix);
1522: ISDestroy(&tiy);
1524: PetscOptionsGetViewer(PetscObjectComm((PetscObject)ctx),((PetscObject)ctx)->prefix,"-vecscatter_view",&viewer,&format,&flag);
1525: if (flag) {
1526: PetscViewerPushFormat(viewer,format);
1527: VecScatterView(ctx,viewer);
1528: PetscViewerPopFormat(viewer);
1529: PetscViewerDestroy(&viewer);
1530: }
1531: return(0);
1532: }
1534: /* ------------------------------------------------------------------*/
1537: /*@
1538: VecScatterGetMerged - Returns true if the scatter is completed in the VecScatterBegin()
1539: and the VecScatterEnd() does nothing
1541: Not Collective
1543: Input Parameter:
1544: . ctx - scatter context created with VecScatterCreate()
1546: Output Parameter:
1547: . flg - PETSC_TRUE if the VecScatterBegin/End() are all done during the VecScatterBegin()
1549: Level: developer
1551: .seealso: VecScatterCreate(), VecScatterEnd(), VecScatterBegin()
1552: @*/
1553: PetscErrorCode VecScatterGetMerged(VecScatter ctx,PetscBool *flg)
1554: {
1557: *flg = ctx->beginandendtogether;
1558: return(0);
1559: }
1563: /*@
1564: VecScatterBegin - Begins a generalized scatter from one vector to
1565: another. Complete the scattering phase with VecScatterEnd().
1567: Neighbor-wise Collective on VecScatter and Vec
1569: Input Parameters:
1570: + inctx - scatter context generated by VecScatterCreate()
1571: . x - the vector from which we scatter
1572: . y - the vector to which we scatter
1573: . addv - either ADD_VALUES or INSERT_VALUES, with INSERT_VALUES mode any location
1574: not scattered to retains its old value; i.e. the vector is NOT first zeroed.
1575: - mode - the scattering mode, usually SCATTER_FORWARD. The available modes are:
1576: SCATTER_FORWARD or SCATTER_REVERSE
1579: Level: intermediate
1581: Options Database: See VecScatterCreate()
1583: Notes:
1584: The vectors x and y need not be the same vectors used in the call
1585: to VecScatterCreate(), but x must have the same parallel data layout
1586: as that passed in as the x to VecScatterCreate(), similarly for the y.
1587: Most likely they have been obtained from VecDuplicate().
1589: You cannot change the values in the input vector between the calls to VecScatterBegin()
1590: and VecScatterEnd().
1592: If you use SCATTER_REVERSE the two arguments x and y should be reversed, from
1593: the SCATTER_FORWARD.
1595: y[iy[i]] = x[ix[i]], for i=0,...,ni-1
1597: This scatter is far more general than the conventional
1598: scatter, since it can be a gather or a scatter or a combination,
1599: depending on the indices ix and iy. If x is a parallel vector and y
1600: is sequential, VecScatterBegin() can serve to gather values to a
1601: single processor. Similarly, if y is parallel and x sequential, the
1602: routine can scatter from one processor to many processors.
1604: Concepts: scatter^between vectors
1605: Concepts: gather^between vectors
1607: .seealso: VecScatterCreate(), VecScatterEnd()
1608: @*/
1609: PetscErrorCode VecScatterBegin(VecScatter inctx,Vec x,Vec y,InsertMode addv,ScatterMode mode)
1610: {
1612: #if defined(PETSC_USE_DEBUG)
1613: PetscInt to_n,from_n;
1614: #endif
1619: if (inctx->inuse) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE," Scatter ctx already in use");
1621: #if defined(PETSC_USE_DEBUG)
1622: /*
1623: Error checking to make sure these vectors match the vectors used
1624: to create the vector scatter context. -1 in the from_n and to_n indicate the
1625: vector lengths are unknown (for example with mapped scatters) and thus
1626: no error checking is performed.
1627: */
1628: if (inctx->from_n >= 0 && inctx->to_n >= 0) {
1629: VecGetLocalSize(x,&from_n);
1630: VecGetLocalSize(y,&to_n);
1631: if (mode & SCATTER_REVERSE) {
1632: if (to_n != inctx->from_n) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Vector wrong size %D for scatter %D (scatter reverse and vector to != ctx from size)",to_n,inctx->from_n);
1633: if (from_n != inctx->to_n) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Vector wrong size %D for scatter %D (scatter reverse and vector from != ctx to size)",from_n,inctx->to_n);
1634: } else {
1635: if (to_n != inctx->to_n) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Vector wrong size %D for scatter %D (scatter forward and vector to != ctx to size)",to_n,inctx->to_n);
1636: if (from_n != inctx->from_n) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Vector wrong size %D for scatter %D (scatter forward and vector from != ctx from size)",from_n,inctx->from_n);
1637: }
1638: }
1639: #endif
1641: inctx->inuse = PETSC_TRUE;
1642: PetscLogEventBarrierBegin(VEC_ScatterBarrier,0,0,0,0,PetscObjectComm((PetscObject)inctx));
1643: (*inctx->begin)(inctx,x,y,addv,mode);
1644: if (inctx->beginandendtogether && inctx->end) {
1645: inctx->inuse = PETSC_FALSE;
1646: (*inctx->end)(inctx,x,y,addv,mode);
1647: }
1648: PetscLogEventBarrierEnd(VEC_ScatterBarrier,0,0,0,0,PetscObjectComm((PetscObject)inctx));
1649: return(0);
1650: }
1652: /* --------------------------------------------------------------------*/
1655: /*@
1656: VecScatterEnd - Ends a generalized scatter from one vector to another. Call
1657: after first calling VecScatterBegin().
1659: Neighbor-wise Collective on VecScatter and Vec
1661: Input Parameters:
1662: + ctx - scatter context generated by VecScatterCreate()
1663: . x - the vector from which we scatter
1664: . y - the vector to which we scatter
1665: . addv - either ADD_VALUES or INSERT_VALUES.
1666: - mode - the scattering mode, usually SCATTER_FORWARD. The available modes are:
1667: SCATTER_FORWARD, SCATTER_REVERSE
1669: Level: intermediate
1671: Notes:
1672: If you use SCATTER_REVERSE the arguments x and y should be reversed, from the SCATTER_FORWARD.
1674: y[iy[i]] = x[ix[i]], for i=0,...,ni-1
1676: .seealso: VecScatterBegin(), VecScatterCreate()
1677: @*/
1678: PetscErrorCode VecScatterEnd(VecScatter ctx,Vec x,Vec y,InsertMode addv,ScatterMode mode)
1679: {
1686: ctx->inuse = PETSC_FALSE;
1687: if (!ctx->end) return(0);
1688: if (!ctx->beginandendtogether) {
1689: PetscLogEventBegin(VEC_ScatterEnd,ctx,x,y,0);
1690: (*(ctx)->end)(ctx,x,y,addv,mode);
1691: PetscLogEventEnd(VEC_ScatterEnd,ctx,x,y,0);
1692: }
1693: return(0);
1694: }
1698: /*@C
1699: VecScatterDestroy - Destroys a scatter context created by
1700: VecScatterCreate().
1702: Collective on VecScatter
1704: Input Parameter:
1705: . ctx - the scatter context
1707: Level: intermediate
1709: .seealso: VecScatterCreate(), VecScatterCopy()
1710: @*/
1711: PetscErrorCode VecScatterDestroy(VecScatter *ctx)
1712: {
1716: if (!*ctx) return(0);
1718: if ((*ctx)->inuse) SETERRQ(((PetscObject)(*ctx))->comm,PETSC_ERR_ARG_WRONGSTATE,"Scatter context is in use");
1719: if (--((PetscObject)(*ctx))->refct > 0) {*ctx = 0; return(0);}
1721: /* if memory was published with AMS then destroy it */
1722: PetscObjectAMSViewOff((PetscObject)(*ctx));
1723: (*(*ctx)->destroy)(*ctx);
1724: PetscHeaderDestroy(ctx);
1725: return(0);
1726: }
1730: /*@
1731: VecScatterCopy - Makes a copy of a scatter context.
1733: Collective on VecScatter
1735: Input Parameter:
1736: . sctx - the scatter context
1738: Output Parameter:
1739: . ctx - the context copy
1741: Level: advanced
1743: .seealso: VecScatterCreate(), VecScatterDestroy()
1744: @*/
1745: PetscErrorCode VecScatterCopy(VecScatter sctx,VecScatter *ctx)
1746: {
1752: if (!sctx->copy) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Cannot copy this type");
1753: PetscHeaderCreate(*ctx,_p_VecScatter,int,VEC_SCATTER_CLASSID,"VecScatter","VecScatter","Vec",PetscObjectComm((PetscObject)sctx),VecScatterDestroy,VecScatterView);
1754: (*ctx)->to_n = sctx->to_n;
1755: (*ctx)->from_n = sctx->from_n;
1756: (*sctx->copy)(sctx,*ctx);
1757: return(0);
1758: }
1761: /* ------------------------------------------------------------------*/
1764: /*@
1765: VecScatterView - Views a vector scatter context.
1767: Collective on VecScatter
1769: Input Parameters:
1770: + ctx - the scatter context
1771: - viewer - the viewer for displaying the context
1773: Level: intermediate
1775: @*/
1776: PetscErrorCode VecScatterView(VecScatter ctx,PetscViewer viewer)
1777: {
1782: if (!viewer) {
1783: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)ctx),&viewer);
1784: }
1786: if (ctx->view) {
1787: (*ctx->view)(ctx,viewer);
1788: }
1789: return(0);
1790: }
1794: /*@C
1795: VecScatterRemap - Remaps the "from" and "to" indices in a
1796: vector scatter context. FOR EXPERTS ONLY!
1798: Collective on VecScatter
1800: Input Parameters:
1801: + scat - vector scatter context
1802: . from - remapping for "from" indices (may be NULL)
1803: - to - remapping for "to" indices (may be NULL)
1805: Level: developer
1807: Notes: In the parallel case the todata is actually the indices
1808: from which the data is TAKEN! The from stuff is where the
1809: data is finally put. This is VERY VERY confusing!
1811: In the sequential case the todata is the indices where the
1812: data is put and the fromdata is where it is taken from.
1813: This is backwards from the paralllel case! CRY! CRY! CRY!
1815: @*/
1816: PetscErrorCode VecScatterRemap(VecScatter scat,PetscInt *rto,PetscInt *rfrom)
1817: {
1818: VecScatter_Seq_General *to,*from;
1819: VecScatter_MPI_General *mto;
1820: PetscInt i;
1827: from = (VecScatter_Seq_General*)scat->fromdata;
1828: mto = (VecScatter_MPI_General*)scat->todata;
1830: if (mto->type == VEC_SCATTER_MPI_TOALL) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Not for to all scatter");
1832: if (rto) {
1833: if (mto->type == VEC_SCATTER_MPI_GENERAL) {
1834: /* handle off processor parts */
1835: for (i=0; i<mto->starts[mto->n]; i++) mto->indices[i] = rto[mto->indices[i]];
1837: /* handle local part */
1838: to = &mto->local;
1839: for (i=0; i<to->n; i++) to->vslots[i] = rto[to->vslots[i]];
1840: } else if (from->type == VEC_SCATTER_SEQ_GENERAL) {
1841: for (i=0; i<from->n; i++) from->vslots[i] = rto[from->vslots[i]];
1842: } else if (from->type == VEC_SCATTER_SEQ_STRIDE) {
1843: VecScatter_Seq_Stride *sto = (VecScatter_Seq_Stride*)from;
1845: /* if the remapping is the identity and stride is identity then skip remap */
1846: if (sto->step == 1 && sto->first == 0) {
1847: for (i=0; i<sto->n; i++) {
1848: if (rto[i] != i) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Unable to remap such scatters");
1849: }
1850: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Unable to remap such scatters");
1851: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Unable to remap such scatters");
1852: }
1854: if (rfrom) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Unable to remap the FROM in scatters yet");
1856: /*
1857: Mark then vector lengths as unknown because we do not know the
1858: lengths of the remapped vectors
1859: */
1860: scat->from_n = -1;
1861: scat->to_n = -1;
1862: return(0);
1863: }