Actual source code: ex2.c
petsc-3.11.4 2019-09-28
1: static const char help[] = "Test overlapped communication on a single star forest (PetscSF)\n\n";
3: #include <petscvec.h>
4: #include <petscsf.h>
5: #include <petscviewer.h>
7: int main(int argc, char **argv)
8: {
9: PetscInt ierr;
10: PetscSF sf;
11: Vec A,B,Aout,Bout;
12: PetscScalar *bufA,*bufB,*bufAout,*bufBout;
13: PetscInt a,b,aout,bout,*bufa,*bufb,*bufaout,*bufbout;
14: PetscMPIInt rank,size;
15: PetscInt i,*ilocal,nroots,nleaves;
16: PetscSFNode *iremote;
18: PetscInitialize(&argc,&argv,NULL,help);if (ierr) return ierr;
19: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
20: MPI_Comm_size(PETSC_COMM_WORLD,&size);
22: if (size != 2) SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_USER, "Only coded for two MPI processes\n");
24: PetscSFCreate(PETSC_COMM_WORLD,&sf);
25: PetscSFSetFromOptions(sf);
27: nleaves = 2;
28: nroots = 1;
29: PetscMalloc1(nleaves,&ilocal);
31: for (i = 0; i<nleaves; i++) {
32: ilocal[i] = i;
33: }
35: PetscMalloc1(nleaves,&iremote);
36: if (rank == 0) {
37: iremote[0].rank = 0;
38: iremote[0].index = 0;
39: iremote[1].rank = 1;
40: iremote[1].index = 0;
41: } else {
42: iremote[0].rank = 1;
43: iremote[0].index = 0;
44: iremote[1].rank = 0;
45: iremote[1].index = 0;
46: }
47: PetscSFSetGraph(sf,nroots,nleaves,ilocal,PETSC_OWN_POINTER,iremote,PETSC_OWN_POINTER);
48: PetscSFSetUp(sf);
49: PetscSFView(sf,PETSC_VIEWER_STDOUT_WORLD);
50: VecCreate(PETSC_COMM_WORLD,&A);
51: VecSetSizes(A,2,PETSC_DETERMINE);
52: VecSetFromOptions(A);
53: VecSetUp(A);
55: VecDuplicate(A,&B);
56: VecDuplicate(A,&Aout);
57: VecDuplicate(A,&Bout);
58: VecGetArray(A,&bufA);
59: VecGetArray(B,&bufB);
60: for (i=0; i<2; i++) {
61: bufA[i] = (PetscScalar)rank;
62: bufB[i] = (PetscScalar)(rank) + 10.0;
63: }
64: VecRestoreArray(A,&bufA);
65: VecRestoreArray(B,&bufB);
67: VecGetArrayRead(A,(const PetscScalar**)&bufA);
68: VecGetArrayRead(B,(const PetscScalar**)&bufB);
69: VecGetArray(Aout,&bufAout);
70: VecGetArray(Bout,&bufBout);
72: /* Testing overlapped PetscSFBcast with different rootdata and leafdata */
73: PetscSFBcastBegin(sf,MPIU_SCALAR,bufA,bufAout);
74: PetscSFBcastBegin(sf,MPIU_SCALAR,bufB,bufBout);
75: PetscSFBcastEnd (sf,MPIU_SCALAR,bufA,bufAout);
76: PetscSFBcastEnd (sf,MPIU_SCALAR,bufB,bufBout);
77: VecRestoreArrayRead(A,(const PetscScalar**)&bufA);
78: VecRestoreArrayRead(B,(const PetscScalar**)&bufB);
79: VecRestoreArray(Aout,&bufAout);
80: VecRestoreArray(Bout,&bufBout);
82: VecView(Aout,PETSC_VIEWER_STDOUT_WORLD);
83: VecView(Bout,PETSC_VIEWER_STDOUT_WORLD);
84: VecDestroy(&A);
85: VecDestroy(&B);
86: VecDestroy(&Aout);
87: VecDestroy(&Bout);
88: PetscSFDestroy(&sf);
90: /* Another very simple graph: rank 0 has one root, zero leaves; rank 1 has zero roots, one leave.
91: Zero roots or leaves will result in NULL rootdata or leafdata. Therefore, one can not use that
92: as key to identify pending communications.
93: */
94: if (!rank) {
95: nroots = 1;
96: nleaves = 0;
97: } else {
98: nroots = 0;
99: nleaves = 1;
100: }
102: PetscMalloc1(nleaves,&ilocal);
103: PetscMalloc1(nleaves,&iremote);
104: if (rank) {
105: ilocal[0] = 0;
106: iremote[0].rank = 0;
107: iremote[0].index = 0;
108: }
110: PetscSFCreate(PETSC_COMM_WORLD,&sf);
111: PetscSFSetFromOptions(sf);
112: PetscSFSetGraph(sf,nroots,nleaves,ilocal,PETSC_OWN_POINTER,iremote,PETSC_OWN_POINTER);
113: PetscSFSetUp(sf);
115: if (!rank) {
116: a = 1;
117: b = 10;
118: bufa = &a;
119: bufb = &b;
120: bufaout = NULL;
121: bufbout = NULL;
122: } else {
123: bufa = NULL;
124: bufb = NULL;
125: bufaout = &aout;
126: bufbout = &bout;
127: }
129: /* Do Bcast to challenge PetscSFBcast if it uses rootdata to identify pending communications.
130: Rank 0 Rank 1
131: rootdata: bufa=XXX (1) bufa=NULL
132: bufb=YYY (10) bufb=NULL
134: leafdata: bufaout=NULL bufaout=WWW
135: bufbout=NULL bufbout=ZZZ
136: */
137: PetscSFBcastBegin(sf,MPIU_INT,bufa,bufaout);
138: PetscSFBcastBegin(sf,MPIU_INT,bufb,bufbout);
139: PetscSFBcastEnd (sf,MPIU_INT,bufa,bufaout);
140: PetscSFBcastEnd (sf,MPIU_INT,bufb,bufbout);
141: if (rank) {PetscSynchronizedPrintf(PETSC_COMM_WORLD,"On rank 1, aout=%D, bout=%D\n",aout,bout);} /* On rank 1, aout=1, bout=10 */
142: PetscSynchronizedFlush(PETSC_COMM_WORLD,NULL);
144: /* Do Reduce to challenge PetscSFReduce if it uses leafdata to identify pending communications.
145: Rank 0 Rank 1
146: rootdata: bufa=XXX (1) bufa=NULL
147: bufb=YYY (10) bufb=NULL
149: leafdata: bufaout=NULL bufaout=WWW (1)
150: bufbout=NULL bufbout=ZZZ (10)
151: */
152: PetscSFReduceBegin(sf,MPIU_INT,bufaout,bufa,MPI_SUM);
153: PetscSFReduceBegin(sf,MPIU_INT,bufbout,bufb,MPI_SUM);
154: PetscSFReduceEnd (sf,MPIU_INT,bufaout,bufa,MPI_SUM);
155: PetscSFReduceEnd (sf,MPIU_INT,bufbout,bufb,MPI_SUM);
156: if (!rank) {PetscSynchronizedPrintf(PETSC_COMM_WORLD,"On rank 0, a =%D, b =%D\n",a,b);} /* On rank 0, a=2, b=20 */
157: PetscSynchronizedFlush(PETSC_COMM_WORLD,NULL);
159: /* A difficult case PETSc can not handle correctly. So comment out. */
160: #if 0
161: /* Note on rank 0, the following overlapped Bcast have the same rootdata (bufa) and leafdata (NULL).
162: It challenges PetscSF if it uses (rootdata, leafdata) as key to identify pending communications.
163: Rank 0 Rank 1
164: rootdata: bufa=XXX (2) bufa=NULL
166: leafdata: bufaout=NULL bufaout=WWW (1)
167: bufbout=NULL bufbout=ZZZ (10)
168: */
169: PetscSFBcastBegin(sf,MPIU_INT,bufa,bufaout);
170: PetscSFBcastBegin(sf,MPIU_INT,bufa,bufbout);
171: PetscSFBcastEnd (sf,MPIU_INT,bufa,bufaout); /* Buggy PetscSF could match a wrong PetscSFBcastBegin */
172: PetscSFBcastEnd (sf,MPIU_INT,bufa,bufbout);
173: if (rank) {PetscSynchronizedPrintf(PETSC_COMM_WORLD,"On rank 1, aout=%D, bout=%D\n",aout,bout);} /* On rank 1, aout=2, bout=2 */
175: /* The above Bcasts may successfully populate leafdata with correct values. But communinication contexts (i.e., the links,
176: each with a unique MPI tag) may have retired in different order in the ->avail list. Processes doing the following
177: PetscSFReduce will get tag-unmatched links from the ->avail list, resulting in dead lock.
178: */
179: PetscSFReduceBegin(sf,MPIU_INT,bufaout,bufa,MPI_SUM);
180: PetscSFReduceEnd (sf,MPIU_INT,bufaout,bufa,MPI_SUM);
181: if (!rank) {PetscSynchronizedPrintf(PETSC_COMM_WORLD,"On rank 0, a =%D, b =%D\n",a,b);} /* On rank 0, a=4, b=20 */
182: PetscSynchronizedFlush(PETSC_COMM_WORLD,NULL);
183: #endif
185: PetscSFDestroy(&sf);
186: PetscFinalize();
187: return ierr;
188: }
191: /*TEST
193: test:
194: suffix: basic
195: nsize: 2
196: args: -sf_type basic
198: test:
199: suffix: window
200: nsize: 2
201: args: -sf_type window
202: requires: define(PETSC_HAVE_MPI_WIN_CREATE)
204: TEST*/