Actual source code: ex5.c

  1: static char help[] = "Test PetscSFFCompose and PetscSFCreateStridedSF when the ilocal arrays are not identity nor dense\n\n";

  3: #include <petsc.h>
  4: #include <petscsf.h>

  6: static PetscErrorCode TestVector(PetscSF sf, const char *sfname)
  7: {
  8:   PetscInt mr, ml;
  9:   MPI_Comm comm;

 11:   PetscFunctionBeginUser;
 12:   comm = PetscObjectComm((PetscObject)sf);
 13:   PetscCall(PetscSFGetGraph(sf, &mr, NULL, NULL, NULL));
 14:   PetscCall(PetscSFGetLeafRange(sf, NULL, &ml));
 15:   for (PetscInt bs = 1; bs < 3; bs++) {
 16:     for (PetscInt r = 0; r < 2; r++) {
 17:       for (PetscInt l = 0; l < 2; l++) {
 18:         PetscSF   vsf;
 19:         PetscInt *rdata, *ldata, *rdatav, *ldatav;
 20:         PetscInt  ldr = PETSC_DECIDE;
 21:         PetscInt  ldl = PETSC_DECIDE;
 22:         PetscBool flg;

 24:         if (r == 1) ldr = mr;
 25:         if (r == 2) ldr = mr + 7;
 26:         if (l == 1) ldl = ml + 1;
 27:         if (l == 2) ldl = ml + 5;

 29:         PetscCall(PetscSFCreateStridedSF(sf, bs, ldr, ldl, &vsf));
 30:         if (ldr == PETSC_DECIDE) ldr = mr;
 31:         if (ldl == PETSC_DECIDE) ldl = ml + 1;

 33:         PetscCall(PetscCalloc4(bs * ldr, &rdata, bs * ldl, &ldata, bs * ldr, &rdatav, bs * ldl, &ldatav));
 34:         for (PetscInt i = 0; i < bs * ldr; i++) rdata[i] = i + 1;

 36:         for (PetscInt i = 0; i < bs; i++) {
 37:           PetscCall(PetscSFBcastBegin(sf, MPIU_INT, PetscSafePointerPlusOffset(rdata, i * ldr), PetscSafePointerPlusOffset(ldata, i * ldl), MPI_REPLACE));
 38:           PetscCall(PetscSFBcastEnd(sf, MPIU_INT, PetscSafePointerPlusOffset(rdata, i * ldr), PetscSafePointerPlusOffset(ldata, i * ldl), MPI_REPLACE));
 39:         }
 40:         PetscCall(PetscSFBcastBegin(vsf, MPIU_INT, rdata, ldatav, MPI_REPLACE));
 41:         PetscCall(PetscSFBcastEnd(vsf, MPIU_INT, rdata, ldatav, MPI_REPLACE));
 42:         PetscCall(PetscArraycmp(ldata, ldatav, bs * ldl, &flg));

 44:         PetscCallMPI(MPIU_Allreduce(MPI_IN_PLACE, &flg, 1, MPIU_BOOL, MPI_LAND, comm));
 45:         if (!flg) {
 46:           PetscCall(PetscPrintf(comm, "Error with Bcast on %s: block size %" PetscInt_FMT ", ldr %" PetscInt_FMT ", ldl %" PetscInt_FMT "\n", sfname, bs, ldr, ldl));
 47:           PetscCall(PetscPrintf(comm, "Single SF\n"));
 48:           PetscCall(PetscIntView(bs * ldl, ldata, PETSC_VIEWER_STDOUT_(comm)));
 49:           PetscCall(PetscPrintf(comm, "Vector SF\n"));
 50:           PetscCall(PetscIntView(bs * ldl, ldatav, PETSC_VIEWER_STDOUT_(comm)));
 51:         }
 52:         PetscCall(PetscArrayzero(rdata, bs * ldr));
 53:         PetscCall(PetscArrayzero(rdatav, bs * ldr));

 55:         for (PetscInt i = 0; i < bs; i++) {
 56:           PetscCall(PetscSFReduceBegin(sf, MPIU_INT, PetscSafePointerPlusOffset(ldata, i * ldl), PetscSafePointerPlusOffset(rdata, i * ldr), MPI_SUM));
 57:           PetscCall(PetscSFReduceEnd(sf, MPIU_INT, PetscSafePointerPlusOffset(ldata, i * ldl), PetscSafePointerPlusOffset(rdata, i * ldr), MPI_SUM));
 58:         }
 59:         PetscCall(PetscSFReduceBegin(vsf, MPIU_INT, ldata, rdatav, MPI_SUM));
 60:         PetscCall(PetscSFReduceEnd(vsf, MPIU_INT, ldata, rdatav, MPI_SUM));
 61:         PetscCall(PetscArraycmp(rdata, rdatav, bs * ldr, &flg));
 62:         PetscCallMPI(MPIU_Allreduce(MPI_IN_PLACE, &flg, 1, MPIU_BOOL, MPI_LAND, comm));
 63:         if (!flg) {
 64:           PetscCall(PetscPrintf(comm, "Error with Reduce on %s: block size %" PetscInt_FMT ", ldr %" PetscInt_FMT ", ldl %" PetscInt_FMT "\n", sfname, bs, ldr, ldl));
 65:           PetscCall(PetscPrintf(comm, "Single SF\n"));
 66:           PetscCall(PetscIntView(bs * ldr, rdata, PETSC_VIEWER_STDOUT_(comm)));
 67:           PetscCall(PetscPrintf(comm, "Vector SF\n"));
 68:           PetscCall(PetscIntView(bs * ldr, rdatav, PETSC_VIEWER_STDOUT_(comm)));
 69:         }
 70:         PetscCall(PetscFree4(rdata, ldata, rdatav, ldatav));
 71:         PetscCall(PetscSFDestroy(&vsf));
 72:       }
 73:     }
 74:   }
 75:   PetscFunctionReturn(PETSC_SUCCESS);
 76: }

 78: int main(int argc, char **argv)
 79: {
 80:   PetscSF      sfA, sfB, sfBA, sfAAm, sfBBm, sfAm, sfBm;
 81:   PetscInt     nrootsA, nleavesA, nrootsB, nleavesB;
 82:   PetscInt    *ilocalA, *ilocalB;
 83:   PetscSFNode *iremoteA, *iremoteB;
 84:   PetscMPIInt  rank, size;
 85:   PetscInt     i, m, n, k, nl = 2, mA, mB, nldataA, nldataB;
 86:   PetscInt    *rdA, *rdB, *ldA, *ldB;
 87:   PetscBool    inverse = PETSC_FALSE, test_vector = PETSC_TRUE;

 89:   PetscFunctionBeginUser;
 90:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
 91:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-nl", &nl, NULL));
 92:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-explicit_inverse", &inverse, NULL));
 93:   PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
 94:   PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));

 96:   PetscCall(PetscSFCreate(PETSC_COMM_WORLD, &sfA));
 97:   PetscCall(PetscSFCreate(PETSC_COMM_WORLD, &sfB));
 98:   PetscCall(PetscSFSetFromOptions(sfA));
 99:   PetscCall(PetscSFSetFromOptions(sfB));

101:   // disable vector tests with linux-misc-32bit and sftype window
102: #if (PETSC_SIZEOF_SIZE_T == 4)
103:   {
104:     PetscBool iswin;

106:     PetscCall(PetscObjectTypeCompare((PetscObject)sfA, PETSCSFWINDOW, &iswin));
107:     if (iswin) test_vector = PETSC_FALSE;
108:   }
109: #endif
110:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-test_vector", &test_vector, NULL));

112:   n = 4 * nl * size;
113:   m = 2 * nl;
114:   k = nl;

116:   nldataA = rank == 0 ? n : 0;
117:   nldataB = 3 * nl;

119:   nrootsA  = m;
120:   nleavesA = rank == 0 ? size * m : 0;
121:   nrootsB  = rank == 0 ? n : 0;
122:   nleavesB = k;

124:   PetscCall(PetscMalloc1(nleavesA, &ilocalA));
125:   PetscCall(PetscMalloc1(nleavesA, &iremoteA));
126:   PetscCall(PetscMalloc1(nleavesB, &ilocalB));
127:   PetscCall(PetscMalloc1(nleavesB, &iremoteB));

129:   /* sf A bcast is equivalent to a sparse gather on process 0
130:      process 0 receives data in the middle [nl,3*nl] of the leaf data array for A */
131:   for (i = 0; i < nleavesA; i++) {
132:     iremoteA[i].rank  = (PetscMPIInt)(i / m);
133:     iremoteA[i].index = i % m;
134:     ilocalA[i]        = nl + i / m * 4 * nl + i % m;
135:   }

137:   /* sf B bcast is equivalent to a sparse scatter from process 0
138:      process 0 sends data from [nl,2*nl] of the leaf data array for A
139:      each process receives, in reverse order, in the middle [nl,2*nl] of the leaf data array for B */
140:   for (i = 0; i < nleavesB; i++) {
141:     iremoteB[i].rank  = 0;
142:     iremoteB[i].index = rank * 4 * nl + nl + i % m;
143:     ilocalB[i]        = 2 * nl - i - 1;
144:   }
145:   PetscCall(PetscSFSetGraph(sfA, nrootsA, nleavesA, ilocalA, PETSC_OWN_POINTER, iremoteA, PETSC_OWN_POINTER));
146:   PetscCall(PetscSFSetGraph(sfB, nrootsB, nleavesB, ilocalB, PETSC_OWN_POINTER, iremoteB, PETSC_OWN_POINTER));
147:   PetscCall(PetscSFSetUp(sfA));
148:   PetscCall(PetscSFSetUp(sfB));
149:   PetscCall(PetscObjectSetName((PetscObject)sfA, "sfA"));
150:   PetscCall(PetscObjectSetName((PetscObject)sfB, "sfB"));
151:   PetscCall(PetscSFViewFromOptions(sfA, NULL, "-view"));
152:   PetscCall(PetscSFViewFromOptions(sfB, NULL, "-view"));

154:   if (test_vector) PetscCall(TestVector(sfA, "sfA"));
155:   if (test_vector) PetscCall(TestVector(sfB, "sfB"));

157:   PetscCall(PetscSFGetLeafRange(sfA, NULL, &mA));
158:   PetscCall(PetscSFGetLeafRange(sfB, NULL, &mB));
159:   PetscCall(PetscMalloc2(nrootsA, &rdA, nldataA, &ldA));
160:   PetscCall(PetscMalloc2(nrootsB, &rdB, nldataB, &ldB));
161:   for (i = 0; i < nrootsA; i++) rdA[i] = m * rank + i;
162:   for (i = 0; i < nldataA; i++) ldA[i] = -1;
163:   for (i = 0; i < nldataB; i++) ldB[i] = -1;

165:   PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "BcastB(BcastA)\n"));
166:   PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "A: root data\n"));
167:   PetscCall(PetscIntView(nrootsA, rdA, PETSC_VIEWER_STDOUT_WORLD));
168:   PetscCall(PetscSFBcastBegin(sfA, MPIU_INT, rdA, ldA, MPI_REPLACE));
169:   PetscCall(PetscSFBcastEnd(sfA, MPIU_INT, rdA, ldA, MPI_REPLACE));
170:   PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "A: leaf data (all)\n"));
171:   PetscCall(PetscIntView(nldataA, ldA, PETSC_VIEWER_STDOUT_WORLD));
172:   PetscCall(PetscSFBcastBegin(sfB, MPIU_INT, ldA, ldB, MPI_REPLACE));
173:   PetscCall(PetscSFBcastEnd(sfB, MPIU_INT, ldA, ldB, MPI_REPLACE));
174:   PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "B: leaf data (all)\n"));
175:   PetscCall(PetscIntView(nldataB, ldB, PETSC_VIEWER_STDOUT_WORLD));

177:   PetscCall(PetscSFCompose(sfA, sfB, &sfBA));
178:   PetscCall(PetscSFSetFromOptions(sfBA));
179:   PetscCall(PetscSFSetUp(sfBA));
180:   PetscCall(PetscObjectSetName((PetscObject)sfBA, "sfBA"));
181:   PetscCall(PetscSFViewFromOptions(sfBA, NULL, "-view"));
182:   if (test_vector) PetscCall(TestVector(sfBA, "sfBA"));

184:   for (i = 0; i < nldataB; i++) ldB[i] = -1;
185:   PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "BcastBA\n"));
186:   PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "BA: root data\n"));
187:   PetscCall(PetscIntView(nrootsA, rdA, PETSC_VIEWER_STDOUT_WORLD));
188:   PetscCall(PetscSFBcastBegin(sfBA, MPIU_INT, rdA, ldB, MPI_REPLACE));
189:   PetscCall(PetscSFBcastEnd(sfBA, MPIU_INT, rdA, ldB, MPI_REPLACE));
190:   PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "BA: leaf data (all)\n"));
191:   PetscCall(PetscIntView(nldataB, ldB, PETSC_VIEWER_STDOUT_WORLD));

193:   PetscCall(PetscSFCreateInverseSF(sfA, &sfAm));
194:   PetscCall(PetscSFSetFromOptions(sfAm));
195:   PetscCall(PetscObjectSetName((PetscObject)sfAm, "sfAm"));
196:   PetscCall(PetscSFViewFromOptions(sfAm, NULL, "-view"));
197:   if (test_vector) PetscCall(TestVector(sfAm, "sfAm"));

199:   if (!inverse) {
200:     PetscCall(PetscSFComposeInverse(sfA, sfA, &sfAAm));
201:   } else {
202:     PetscCall(PetscSFCompose(sfA, sfAm, &sfAAm));
203:   }
204:   PetscCall(PetscSFSetFromOptions(sfAAm));
205:   PetscCall(PetscSFSetUp(sfAAm));
206:   PetscCall(PetscObjectSetName((PetscObject)sfAAm, "sfAAm"));
207:   PetscCall(PetscSFViewFromOptions(sfAAm, NULL, "-view"));
208:   if (test_vector) PetscCall(TestVector(sfAAm, "sfAAm"));

210:   PetscCall(PetscSFCreateInverseSF(sfB, &sfBm));
211:   PetscCall(PetscSFSetFromOptions(sfBm));
212:   PetscCall(PetscObjectSetName((PetscObject)sfBm, "sfBm"));
213:   PetscCall(PetscSFViewFromOptions(sfBm, NULL, "-view"));
214:   if (test_vector) PetscCall(TestVector(sfBm, "sfBm"));

216:   if (!inverse) {
217:     PetscCall(PetscSFComposeInverse(sfB, sfB, &sfBBm));
218:   } else {
219:     PetscCall(PetscSFCompose(sfB, sfBm, &sfBBm));
220:   }
221:   PetscCall(PetscSFSetFromOptions(sfBBm));
222:   PetscCall(PetscSFSetUp(sfBBm));
223:   PetscCall(PetscObjectSetName((PetscObject)sfBBm, "sfBBm"));
224:   PetscCall(PetscSFViewFromOptions(sfBBm, NULL, "-view"));
225:   if (test_vector) PetscCall(TestVector(sfBBm, "sfBBm"));

227:   PetscCall(PetscFree2(rdA, ldA));
228:   PetscCall(PetscFree2(rdB, ldB));

230:   PetscCall(PetscSFDestroy(&sfA));
231:   PetscCall(PetscSFDestroy(&sfB));
232:   PetscCall(PetscSFDestroy(&sfBA));
233:   PetscCall(PetscSFDestroy(&sfAm));
234:   PetscCall(PetscSFDestroy(&sfBm));
235:   PetscCall(PetscSFDestroy(&sfAAm));
236:   PetscCall(PetscSFDestroy(&sfBBm));

238:   PetscCall(PetscFinalize());
239:   return 0;
240: }

242: /*TEST

244:    test:
245:      suffix: 1
246:      args: -view -explicit_inverse {{0 1}}

248:    test:
249:      nsize: 7
250:      filter: grep -v "type" | grep -v "sort"
251:      suffix: 2
252:      args: -view -nl 5 -explicit_inverse {{0 1}}

254:    # we cannot test for -sf_window_flavor dynamic because SFCompose with sparse leaves may change the root data pointer only locally, and this is not supported by the dynamic case
255:    test:
256:      nsize: 7
257:      suffix: 2_window
258:      filter: grep -v "type" | grep -v "sort"
259:      output_file: output/ex5_2.out
260:      args: -view -nl 5 -explicit_inverse {{0 1}} -sf_type window -sf_window_sync {{fence lock active}} -sf_window_flavor {{create allocate}}
261:      requires: defined(PETSC_HAVE_MPI_ONE_SIDED) defined(PETSC_HAVE_MPI_FEATURE_DYNAMIC_WINDOW)

263:    # The nightly test suite with MPICH uses ch3:sock, which is broken when winsize == 0 in some of the processes
264:    test:
265:      nsize: 7
266:      suffix: 2_window_shared
267:      filter: grep -v "type" | grep -v "sort"
268:      output_file: output/ex5_2.out
269:      args: -view -nl 5 -explicit_inverse {{0 1}} -sf_type window -sf_window_sync {{fence lock active}} -sf_window_flavor shared
270:      requires: defined(PETSC_HAVE_MPI_PROCESS_SHARED_MEMORY) !defined(PETSC_HAVE_MPICH) defined(PETSC_HAVE_MPI_ONE_SIDED)

272: TEST*/