Actual source code: vpscat.h
2: /*
3: Defines the methods VecScatterBegin/End_1,2,......
4: This is included by vpscat.c with different values for BS
6: This is a terrible way of doing "templates" in C.
7: */
8: #define PETSCMAP1_a(a,b) a ## _ ## b
9: #define PETSCMAP1_b(a,b) PETSCMAP1_a(a,b)
10: #define PETSCMAP1(a) PETSCMAP1_b(a,BS)
14: PetscErrorCode PETSCMAP1(VecScatterBegin)(VecScatter ctx,Vec xin,Vec yin,InsertMode addv,ScatterMode mode)
15: {
16: VecScatter_MPI_General *to,*from;
17: PetscScalar *xv,*yv,*svalues;
18: MPI_Request *rwaits,*swaits;
19: PetscErrorCode ierr;
20: PetscInt i,*indices,*sstarts,nrecvs,nsends,bs;
23: CHKMEMQ;
24: VecGetArray(xin,&xv);
25: if (xin != yin) {VecGetArray(yin,&yv);} else {yv = xv;}
26: if (mode & SCATTER_REVERSE) {
27: to = (VecScatter_MPI_General*)ctx->fromdata;
28: from = (VecScatter_MPI_General*)ctx->todata;
29: rwaits = from->rev_requests;
30: swaits = to->rev_requests;
31: } else {
32: to = (VecScatter_MPI_General*)ctx->todata;
33: from = (VecScatter_MPI_General*)ctx->fromdata;
34: rwaits = from->requests;
35: swaits = to->requests;
36: }
37: bs = to->bs;
38: svalues = to->values;
39: nrecvs = from->n;
40: nsends = to->n;
41: indices = to->indices;
42: sstarts = to->starts;
44: if (!(mode & SCATTER_LOCAL)) {
46: if (!from->use_readyreceiver && !to->sendfirst && !to->use_alltoallv & !to->use_window) {
47: /* post receives since they were not previously posted */
48: if (nrecvs) {MPI_Startall_irecv(from->starts[nrecvs]*bs,nrecvs,rwaits);}
49: }
51: #if defined(PETSC_HAVE_MPI_ALLTOALLW)
52: if (to->use_alltoallw && addv == INSERT_VALUES) {
53: MPI_Alltoallw(xv,to->wcounts,to->wdispls,to->types,yv,from->wcounts,from->wdispls,from->types,ctx->comm);
54: } else
55: #endif
56: if (ctx->packtogether || to->use_alltoallv) {
57: /* this version packs all the messages together and sends, when -vecscatter_packtogether used */
58: PETSCMAP1(Pack)(sstarts[nsends],indices,xv,svalues);
59: if (to->use_alltoallv) {
60: MPI_Alltoallv(to->values,to->counts,to->displs,MPIU_SCALAR,from->values,from->counts,from->displs,MPIU_SCALAR,ctx->comm);
61: #if defined(PETSC_HAVE_MPI_WIN_CREATE)
62: } else if (to->use_window) {
63: PetscInt cnt;
65: MPI_Win_fence(0,from->window);
66: for (i=0; i<nsends; i++) {
67: cnt = bs*(to->starts[i+1]-to->starts[i]);
68: MPI_Put(to->values+to->starts[i],cnt,MPIU_SCALAR,to->procs[i],0/*disp*/,cnt,MPIU_SCALAR,from->window);
69: }
70: MPI_Win_fence(0,from->window);
71: #endif
72: } else if (nsends) {
73: MPI_Startall_isend(to->starts[to->n],nsends,swaits);
74: }
75: } else {
76: /* this version packs and sends one at a time */
77: for (i=0; i<nsends; i++) {
78: PETSCMAP1(Pack)(sstarts[i+1]-sstarts[i],indices + sstarts[i],xv,svalues + bs*sstarts[i]);
79: MPI_Start_isend(sstarts[i+1]-sstarts[i],swaits+i);
80: }
81: }
83: if (!from->use_readyreceiver && to->sendfirst && !to->use_alltoallv && !to->use_window) {
84: /* post receives since they were not previously posted */
85: if (nrecvs) {MPI_Startall_irecv(from->starts[nrecvs]*bs,nrecvs,rwaits);}
86: }
87: }
89: /* take care of local scatters */
90: if (to->local.n) {
91: if (to->local.is_copy && addv == INSERT_VALUES) {
92: PetscMemcpy(yv + from->local.copy_start,xv + to->local.copy_start,to->local.copy_length);
93: } else {
94: PETSCMAP1(Scatter)(to->local.n,to->local.vslots,xv,from->local.vslots,yv,addv);
95: }
96: }
97: VecRestoreArray(xin,&xv);
98: if (xin != yin) {VecRestoreArray(yin,&yv);}
99: CHKMEMQ;
100: return(0);
101: }
103: /* --------------------------------------------------------------------------------------*/
107: PetscErrorCode PETSCMAP1(VecScatterEnd)(VecScatter ctx,Vec xin,Vec yin,InsertMode addv,ScatterMode mode)
108: {
109: VecScatter_MPI_General *to,*from;
110: PetscScalar *rvalues,*yv;
111: PetscErrorCode ierr;
112: PetscInt nrecvs,nsends,*indices,count,*rstarts,bs;
113: PetscMPIInt imdex;
114: MPI_Request *rwaits,*swaits;
115: MPI_Status xrstatus,*rstatus,*sstatus;
118: CHKMEMQ;
119: if (mode & SCATTER_LOCAL) return(0);
120: VecGetArray(yin,&yv);
122: to = (VecScatter_MPI_General*)ctx->todata;
123: from = (VecScatter_MPI_General*)ctx->fromdata;
124: rwaits = from->requests;
125: swaits = to->requests;
126: sstatus = to->sstatus; /* sstatus and rstatus are always stored in to */
127: rstatus = to->rstatus;
128: if (mode & SCATTER_REVERSE) {
129: to = (VecScatter_MPI_General*)ctx->fromdata;
130: from = (VecScatter_MPI_General*)ctx->todata;
131: rwaits = from->rev_requests;
132: swaits = to->rev_requests;
133: }
134: bs = from->bs;
135: rvalues = from->values;
136: nrecvs = from->n;
137: nsends = to->n;
138: indices = from->indices;
139: rstarts = from->starts;
141: if (ctx->packtogether || (to->use_alltoallw && (addv != INSERT_VALUES)) || (to->use_alltoallv && !to->use_alltoallw)) {
142: if (nrecvs && !to->use_alltoallv && !to->use_window) {MPI_Waitall(nrecvs,rwaits,rstatus);}
143: PETSCMAP1(UnPack)(from->starts[from->n],from->values,indices,yv,addv);
144: } else if (!to->use_alltoallw) {
145: /* unpack one at a time */
146: count = nrecvs;
147: while (count) {
148: MPI_Waitany(nrecvs,rwaits,&imdex,&xrstatus);
149: /* unpack receives into our local space */
150: PETSCMAP1(UnPack)(rstarts[imdex+1] - rstarts[imdex],rvalues + bs*rstarts[imdex],indices + rstarts[imdex],yv,addv);
151: count--;
152: }
153: }
154: if (from->use_readyreceiver) {
155: if (nrecvs) {MPI_Startall_irecv(from->starts[nrecvs]*bs,nrecvs,rwaits);}
156: MPI_Barrier(ctx->comm);
157: }
159: /* wait on sends */
160: if (nsends && !to->use_alltoallv && !to->use_window) {MPI_Waitall(nsends,swaits,sstatus);}
161: VecRestoreArray(yin,&yv);
162: CHKMEMQ;
163: return(0);
164: }
166: #undef PETSCMAP1_a
167: #undef PETSCMAP1_b
168: #undef PETSCMAP1
169: #undef BS