Actual source code: ex8.c
1: static char help[] = "Demonstrates BuildTwoSided functions.\n";
3: #include <petscsys.h>
5: typedef struct {
6: PetscInt rank;
7: PetscScalar value;
8: char ok[3];
9: } Unit;
11: static PetscErrorCode MakeDatatype(MPI_Datatype *dtype)
12: {
13: MPI_Datatype dtypes[3],tmptype;
14: PetscMPIInt lengths[3];
15: MPI_Aint displs[3];
16: Unit dummy;
18: dtypes[0] = MPIU_INT;
19: dtypes[1] = MPIU_SCALAR;
20: dtypes[2] = MPI_CHAR;
21: lengths[0] = 1;
22: lengths[1] = 1;
23: lengths[2] = 3;
24: /* Curse the evil beings that made std::complex a non-POD type. */
25: displs[0] = (char*)&dummy.rank - (char*)&dummy; /* offsetof(Unit,rank); */
26: displs[1] = (char*)&dummy.value - (char*)&dummy; /* offsetof(Unit,value); */
27: displs[2] = (char*)&dummy.ok - (char*)&dummy; /* offsetof(Unit,ok); */
28: MPI_Type_create_struct(3,lengths,displs,dtypes,&tmptype);
29: MPI_Type_commit(&tmptype);
30: MPI_Type_create_resized(tmptype,0,sizeof(Unit),dtype);
31: MPI_Type_commit(dtype);
32: MPI_Type_free(&tmptype);
33: {
34: MPI_Aint lb,extent;
35: MPI_Type_get_extent(*dtype,&lb,&extent);
37: }
38: return 0;
39: }
41: struct FCtx {
42: PetscMPIInt rank;
43: PetscMPIInt nto;
44: PetscMPIInt *toranks;
45: Unit *todata;
46: PetscSegBuffer seg;
47: };
49: static PetscErrorCode FSend(MPI_Comm comm,const PetscMPIInt tag[],PetscMPIInt tonum,PetscMPIInt rank,void *todata,MPI_Request req[],void *ctx)
50: {
51: struct FCtx *fctx = (struct FCtx*)ctx;
55: MPI_Isend(&fctx->todata[tonum].rank,1,MPIU_INT,rank,tag[0],comm,&req[0]);
56: MPI_Isend(&fctx->todata[tonum].value,1,MPIU_SCALAR,rank,tag[1],comm,&req[1]);
57: return 0;
58: }
60: static PetscErrorCode FRecv(MPI_Comm comm,const PetscMPIInt tag[],PetscMPIInt rank,void *fromdata,MPI_Request req[],void *ctx)
61: {
62: struct FCtx *fctx = (struct FCtx*)ctx;
63: Unit *buf;
66: PetscSegBufferGet(fctx->seg,1,&buf);
67: MPI_Irecv(&buf->rank,1,MPIU_INT,rank,tag[0],comm,&req[0]);
68: MPI_Irecv(&buf->value,1,MPIU_SCALAR,rank,tag[1],comm,&req[1]);
69: buf->ok[0] = 'o';
70: buf->ok[1] = 'k';
71: buf->ok[2] = 0;
72: return 0;
73: }
75: int main(int argc,char **argv)
76: {
77: PetscMPIInt rank,size,*toranks,*fromranks,nto,nfrom;
78: PetscInt i,n;
79: PetscBool verbose,build_twosided_f;
80: Unit *todata,*fromdata;
81: MPI_Datatype dtype;
83: PetscInitialize(&argc,&argv,(char*)0,help);
84: MPI_Comm_size(PETSC_COMM_WORLD,&size);
85: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
87: verbose = PETSC_FALSE;
88: PetscOptionsGetBool(NULL,NULL,"-verbose",&verbose,NULL);
89: build_twosided_f = PETSC_FALSE;
90: PetscOptionsGetBool(NULL,NULL,"-build_twosided_f",&build_twosided_f,NULL);
92: for (i=1,nto=0; i<size; i*=2) nto++;
93: PetscMalloc2(nto,&todata,nto,&toranks);
94: for (n=0,i=1; i<size; n++,i*=2) {
95: toranks[n] = (rank+i) % size;
96: todata[n].rank = (rank+i) % size;
97: todata[n].value = (PetscScalar)rank;
98: todata[n].ok[0] = 'o';
99: todata[n].ok[1] = 'k';
100: todata[n].ok[2] = 0;
101: }
102: if (verbose) {
103: for (i=0; i<nto; i++) {
104: PetscSynchronizedPrintf(PETSC_COMM_WORLD,"[%d] TO %d: {%" PetscInt_FMT ", %g, \"%s\"}\n",rank,toranks[i],todata[i].rank,(double)PetscRealPart(todata[i].value),todata[i].ok);
105: }
106: PetscSynchronizedFlush(PETSC_COMM_WORLD,PETSC_STDOUT);
107: }
109: MakeDatatype(&dtype);
111: if (build_twosided_f) {
112: struct FCtx fctx;
113: PetscMPIInt *todummy,*fromdummy;
114: fctx.rank = rank;
115: fctx.nto = nto;
116: fctx.toranks = toranks;
117: fctx.todata = todata;
118: PetscSegBufferCreate(sizeof(Unit),1,&fctx.seg);
119: PetscMalloc1(nto,&todummy);
120: for (i=0; i<nto; i++) todummy[i] = rank;
121: PetscCommBuildTwoSidedF(PETSC_COMM_WORLD,1,MPI_INT,nto,toranks,todummy,&nfrom,&fromranks,&fromdummy,2,FSend,FRecv,&fctx);
122: PetscFree(todummy);
123: PetscFree(fromdummy);
124: PetscSegBufferExtractAlloc(fctx.seg,&fromdata);
125: PetscSegBufferDestroy(&fctx.seg);
126: } else {
127: PetscCommBuildTwoSided(PETSC_COMM_WORLD,1,dtype,nto,toranks,todata,&nfrom,&fromranks,&fromdata);
128: }
129: MPI_Type_free(&dtype);
131: if (verbose) {
132: PetscInt *iranks,*iperm;
133: PetscMalloc2(nfrom,&iranks,nfrom,&iperm);
134: for (i=0; i<nfrom; i++) {
135: iranks[i] = fromranks[i];
136: iperm[i] = i;
137: }
138: /* Receive ordering is non-deterministic in general, so sort to make verbose output deterministic. */
139: PetscSortIntWithPermutation(nfrom,iranks,iperm);
140: for (i=0; i<nfrom; i++) {
141: PetscInt ip = iperm[i];
142: PetscSynchronizedPrintf(PETSC_COMM_WORLD,"[%d] FROM %d: {%" PetscInt_FMT ", %g, \"%s\"}\n",rank,fromranks[ip],fromdata[ip].rank,(double)PetscRealPart(fromdata[ip].value),fromdata[ip].ok);
143: }
144: PetscSynchronizedFlush(PETSC_COMM_WORLD,PETSC_STDOUT);
145: PetscFree2(iranks,iperm);
146: }
149: for (i=1; i<size; i*=2) {
150: PetscMPIInt expected_rank = (rank-i+size)%size;
151: PetscBool flg;
152: for (n=0; n<nfrom; n++) {
153: if (expected_rank == fromranks[n]) goto found;
154: }
155: SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"[%d] Could not find expected from rank %d",rank,expected_rank);
156: found:
158: PetscStrcmp(fromdata[n].ok,"ok",&flg);
160: }
161: PetscFree2(todata,toranks);
162: PetscFree(fromdata);
163: PetscFree(fromranks);
164: PetscFinalize();
165: return 0;
166: }
168: /*TEST
170: test:
171: nsize: 4
172: args: -verbose -build_twosided allreduce
174: test:
175: suffix: f
176: nsize: 4
177: args: -verbose -build_twosided_f -build_twosided allreduce
178: output_file: output/ex8_1.out
180: test:
181: suffix: f_ibarrier
182: nsize: 4
183: args: -verbose -build_twosided_f -build_twosided ibarrier
184: output_file: output/ex8_1.out
185: requires: defined(PETSC_HAVE_MPI_NONBLOCKING_COLLECTIVES)
187: test:
188: suffix: ibarrier
189: nsize: 4
190: args: -verbose -build_twosided ibarrier
191: output_file: output/ex8_1.out
192: requires: defined(PETSC_HAVE_MPI_NONBLOCKING_COLLECTIVES)
194: test:
195: suffix: redscatter
196: requires: mpi_reduce_scatter_block
197: nsize: 4
198: args: -verbose -build_twosided redscatter
199: output_file: output/ex8_1.out
201: TEST*/