Actual source code: ex13.c

  1: static char help[] = "Test DMStagPopulateLocalToGlobalInjective.\n\n";

  3: #include <petscdm.h>
  4: #include <petscdmstag.h>

  6: static PetscErrorCode Test1(DM dm);
  7: static PetscErrorCode Test2_1d(DM dm);
  8: static PetscErrorCode Test2_2d(DM dm);
  9: static PetscErrorCode Test2_3d(DM dm);

 11: int main(int argc, char **argv)
 12: {
 13:   DM        dm;
 14:   PetscInt  dim;
 15:   PetscBool setSizes, useInjective;

 17:   /* Initialize PETSc and process command line arguments */
 18:   PetscFunctionBeginUser;
 19:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
 20:   dim = 2;
 21:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-dim", &dim, NULL));
 22:   setSizes = PETSC_FALSE;
 23:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-setsizes", &setSizes, NULL));
 24:   useInjective = PETSC_TRUE;
 25:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-useinjective", &useInjective, NULL));

 27:   /* Creation (normal) */
 28:   if (!setSizes) {
 29:     switch (dim) {
 30:     case 1:
 31:       PetscCall(DMStagCreate1d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, 3, 1, 1, DMSTAG_STENCIL_BOX, 1, NULL, &dm));
 32:       break;
 33:     case 2:
 34:       PetscCall(DMStagCreate2d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, 3, 2, PETSC_DECIDE, PETSC_DECIDE, 1, 1, 1, DMSTAG_STENCIL_BOX, 1, NULL, NULL, &dm));
 35:       break;
 36:     case 3:
 37:       PetscCall(DMStagCreate3d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, 3, 2, 4, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, 1, 1, 1, 1, DMSTAG_STENCIL_BOX, 1, NULL, NULL, NULL, &dm));
 38:       break;
 39:     default:
 40:       SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_SUP, "No support for dimension %" PetscInt_FMT, dim);
 41:     }
 42:   } else {
 43:     /* Creation (test providing decomp exactly)*/
 44:     PetscMPIInt size;
 45:     PetscInt    lx[4] = {2, 3, 4}, ranksx = 3, mx = 9;
 46:     PetscInt    ly[3] = {4, 5}, ranksy = 2, my = 9;
 47:     PetscInt    lz[2] = {6, 7}, ranksz = 2, mz = 13;

 49:     PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
 50:     switch (dim) {
 51:     case 1:
 52:       PetscCheck(size == ranksx, PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Must run on %" PetscInt_FMT " ranks with -dim 1 -setSizes", ranksx);
 53:       PetscCall(DMStagCreate1d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, mx, 1, 1, DMSTAG_STENCIL_BOX, 1, lx, &dm));
 54:       break;
 55:     case 2:
 56:       PetscCheck(size == ranksx * ranksy, PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Must run on %" PetscInt_FMT " ranks with -dim 2 -setSizes", ranksx * ranksy);
 57:       PetscCall(DMStagCreate2d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, mx, my, ranksx, ranksy, 1, 1, 1, DMSTAG_STENCIL_BOX, 1, lx, ly, &dm));
 58:       break;
 59:     case 3:
 60:       PetscCheck(size == ranksx * ranksy * ranksz, PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Must run on %" PetscInt_FMT " ranks with -dim 3 -setSizes", ranksx * ranksy * ranksz);
 61:       PetscCall(DMStagCreate3d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, mx, my, mz, ranksx, ranksy, ranksz, 1, 1, 1, 1, DMSTAG_STENCIL_BOX, 1, lx, ly, lz, &dm));
 62:       break;
 63:     default:
 64:       SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_SUP, "No support for dimension %" PetscInt_FMT, dim);
 65:     }
 66:   }

 68:   /* Setup */
 69:   PetscCall(DMSetFromOptions(dm));
 70:   PetscCall(DMSetUp(dm));

 72:   /* Populate Additional Injective Local-to-Global Map */
 73:   if (useInjective) PetscCall(DMStagPopulateLocalToGlobalInjective(dm));

 75:   /* Test: Make sure L2G inverts G2L */
 76:   PetscCall(Test1(dm));

 78:   /* Test: Make sure that G2L inverts L2G, on its domain */
 79:   PetscCall(DMGetDimension(dm, &dim));
 80:   switch (dim) {
 81:   case 1:
 82:     PetscCall(Test2_1d(dm));
 83:     break;
 84:   case 2:
 85:     PetscCall(Test2_2d(dm));
 86:     break;
 87:   case 3:
 88:     PetscCall(Test2_3d(dm));
 89:     break;
 90:   default:
 91:     SETERRQ(PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Not implemented for dimension %" PetscInt_FMT, dim);
 92:   }

 94:   /* Clean up and finalize PETSc */
 95:   PetscCall(DMDestroy(&dm));
 96:   PetscCall(PetscFinalize());
 97:   return 0;
 98: }

100: static PetscErrorCode Test1(DM dm)
101: {
102:   Vec         vecLocal, vecGlobal, vecGlobalCheck;
103:   PetscRandom rctx;
104:   PetscBool   equal;

106:   PetscFunctionBeginUser;
107:   PetscCall(DMCreateLocalVector(dm, &vecLocal));
108:   PetscCall(DMCreateGlobalVector(dm, &vecGlobal));
109:   PetscCall(PetscRandomCreate(PETSC_COMM_WORLD, &rctx));
110:   PetscCall(VecSetRandom(vecGlobal, rctx));
111:   PetscCall(VecSetRandom(vecLocal, rctx)); /* garbage */
112:   PetscCall(PetscRandomDestroy(&rctx));
113:   PetscCall(VecDuplicate(vecGlobal, &vecGlobalCheck));
114:   PetscCall(DMGlobalToLocal(dm, vecGlobal, INSERT_VALUES, vecLocal));
115:   PetscCall(DMLocalToGlobal(dm, vecLocal, INSERT_VALUES, vecGlobalCheck));
116:   PetscCall(VecEqual(vecGlobal, vecGlobalCheck, &equal));
117:   PetscCheck(equal, PetscObjectComm((PetscObject)dm), PETSC_ERR_PLIB, "Check failed - vectors should be bitwise identical");
118:   PetscCall(VecDestroy(&vecLocal));
119:   PetscCall(VecDestroy(&vecGlobal));
120:   PetscCall(VecDestroy(&vecGlobalCheck));
121:   PetscFunctionReturn(PETSC_SUCCESS);
122: }

124: /* Test function with positive values for positive arguments */
125: #define TEST_FUNCTION(i, j, k, idx, c) (8.33 * i + 7.343 * j + 1.234 * idx + 99.011 * c)

127: /* Helper function to check */
128: static PetscErrorCode CompareValues(PetscInt i, PetscInt j, PetscInt k, PetscInt c, PetscScalar val, PetscScalar valRef)
129: {
130:   PetscFunctionBeginUser;
131:   if (val != valRef && PetscAbsScalar(val - valRef) / PetscAbsScalar(valRef) > 10 * PETSC_MACHINE_EPSILON) {
132:     SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "(%" PetscInt_FMT ",%" PetscInt_FMT ",%" PetscInt_FMT ",%" PetscInt_FMT ") Value %.17g does not match the expected %.17g", i, j, k, c, (double)PetscRealPart(val), (double)PetscRealPart(valRef));
133:   }
134:   PetscFunctionReturn(PETSC_SUCCESS);
135: }

137: static PetscErrorCode Test2_1d(DM dm)
138: {
139:   Vec            vecLocal, vecLocalCheck, vecGlobal;
140:   PetscInt       i, startx, nx, nExtrax, dof0, dof1, c, idxLeft, idxElement;
141:   PetscScalar  **arr;
142:   const PetscInt j = -1, k = -1;

144:   PetscFunctionBeginUser;
145:   PetscCall(DMCreateLocalVector(dm, &vecLocal));
146:   PetscCall(VecSet(vecLocal, -1.0));
147:   PetscCall(DMStagGetCorners(dm, &startx, NULL, NULL, &nx, NULL, NULL, &nExtrax, NULL, NULL));
148:   PetscCall(DMStagGetDOF(dm, &dof0, &dof1, NULL, NULL));
149:   PetscCall(DMStagVecGetArray(dm, vecLocal, &arr));
150:   if (dof0 > 0) PetscCall(DMStagGetLocationSlot(dm, DMSTAG_LEFT, 0, &idxLeft));
151:   if (dof1 > 0) PetscCall(DMStagGetLocationSlot(dm, DMSTAG_ELEMENT, 0, &idxElement));
152:   for (i = startx; i < startx + nx + nExtrax; ++i) {
153:     for (c = 0; c < dof0; ++c) {
154:       const PetscScalar valRef = TEST_FUNCTION(i, 0, 0, idxLeft, c);
155:       arr[i][idxLeft + c]      = valRef;
156:     }
157:     if (i < startx + nx) {
158:       for (c = 0; c < dof1; ++c) {
159:         const PetscScalar valRef = TEST_FUNCTION(i, 0, 0, idxElement, c);
160:         arr[i][idxElement + c]   = valRef;
161:       }
162:     }
163:   }
164:   PetscCall(DMStagVecRestoreArray(dm, vecLocal, &arr));
165:   PetscCall(DMCreateGlobalVector(dm, &vecGlobal));
166:   PetscCall(DMLocalToGlobal(dm, vecLocal, INSERT_VALUES, vecGlobal));
167:   PetscCall(VecDuplicate(vecLocal, &vecLocalCheck));
168:   PetscCall(VecSet(vecLocalCheck, -1.0));
169:   PetscCall(DMGlobalToLocal(dm, vecGlobal, INSERT_VALUES, vecLocalCheck));
170:   PetscCall(DMStagVecGetArrayRead(dm, vecLocalCheck, &arr));
171:   for (i = startx; i < startx + nx + nExtrax; ++i) {
172:     for (c = 0; c < dof0; ++c) {
173:       const PetscScalar valRef = TEST_FUNCTION(i, 0, 0, idxLeft, c);
174:       const PetscScalar val    = arr[i][idxLeft + c];
175:       PetscCall(CompareValues(i, j, k, c, val, valRef));
176:     }
177:     if (i < startx + nx) {
178:       for (c = 0; c < dof1; ++c) {
179:         const PetscScalar valRef = TEST_FUNCTION(i, 0, 0, idxElement, c);
180:         const PetscScalar val    = arr[i][idxElement + c];
181:         PetscCall(CompareValues(i, j, k, c, val, valRef));
182:       }
183:     } else {
184:       for (c = 0; c < dof1; ++c) {
185:         const PetscScalar valRef = -1.0;
186:         const PetscScalar val    = arr[i][idxElement + c];
187:         PetscCall(CompareValues(i, j, k, c, val, valRef));
188:       }
189:     }
190:   }
191:   PetscCall(DMStagVecRestoreArrayRead(dm, vecLocalCheck, &arr));
192:   PetscCall(VecDestroy(&vecLocal));
193:   PetscCall(VecDestroy(&vecLocalCheck));
194:   PetscCall(VecDestroy(&vecGlobal));
195:   PetscFunctionReturn(PETSC_SUCCESS);
196: }

198: static PetscErrorCode Test2_2d(DM dm)
199: {
200:   Vec            vecLocal, vecLocalCheck, vecGlobal;
201:   PetscInt       i, j, startx, starty, nx, ny, nExtrax, nExtray, dof0, dof1, dof2, c, idxLeft, idxDown, idxDownLeft, idxElement;
202:   PetscScalar ***arr;
203:   const PetscInt k = -1;

205:   PetscFunctionBeginUser;
206:   PetscCall(DMCreateLocalVector(dm, &vecLocal));
207:   PetscCall(VecSet(vecLocal, -1.0));
208:   PetscCall(DMStagGetCorners(dm, &startx, &starty, NULL, &nx, &ny, NULL, &nExtrax, &nExtray, NULL));
209:   PetscCall(DMStagGetDOF(dm, &dof0, &dof1, &dof2, NULL));
210:   PetscCall(DMStagVecGetArray(dm, vecLocal, &arr));
211:   if (dof0 > 0) PetscCall(DMStagGetLocationSlot(dm, DMSTAG_DOWN_LEFT, 0, &idxDownLeft));
212:   if (dof1 > 0) {
213:     PetscCall(DMStagGetLocationSlot(dm, DMSTAG_LEFT, 0, &idxLeft));
214:     PetscCall(DMStagGetLocationSlot(dm, DMSTAG_DOWN, 0, &idxDown));
215:   }
216:   if (dof2 > 0) PetscCall(DMStagGetLocationSlot(dm, DMSTAG_ELEMENT, 0, &idxElement));
217:   for (j = starty; j < starty + ny + nExtray; ++j) {
218:     for (i = startx; i < startx + nx + nExtrax; ++i) {
219:       for (c = 0; c < dof0; ++c) {
220:         const PetscScalar valRef   = TEST_FUNCTION(i, j, 0, idxDownLeft, c);
221:         arr[j][i][idxDownLeft + c] = valRef;
222:       }
223:       if (j < starty + ny) {
224:         for (c = 0; c < dof1; ++c) {
225:           const PetscScalar valRef = TEST_FUNCTION(i, j, 0, idxLeft, c);
226:           arr[j][i][idxLeft + c]   = valRef;
227:         }
228:       }
229:       if (i < startx + nx) {
230:         for (c = 0; c < dof1; ++c) {
231:           const PetscScalar valRef = TEST_FUNCTION(i, j, 0, idxDown, c);
232:           arr[j][i][idxDown + c]   = valRef;
233:         }
234:       }
235:       if (i < startx + nx && j < starty + ny) {
236:         for (c = 0; c < dof2; ++c) {
237:           const PetscScalar valRef  = TEST_FUNCTION(i, j, 0, idxElement, c);
238:           arr[j][i][idxElement + c] = valRef;
239:         }
240:       }
241:     }
242:   }
243:   PetscCall(DMStagVecRestoreArray(dm, vecLocal, &arr));
244:   PetscCall(DMCreateGlobalVector(dm, &vecGlobal));
245:   PetscCall(DMLocalToGlobal(dm, vecLocal, INSERT_VALUES, vecGlobal));
246:   PetscCall(VecDuplicate(vecLocal, &vecLocalCheck));
247:   PetscCall(VecSet(vecLocalCheck, -1.0));
248:   PetscCall(DMGlobalToLocal(dm, vecGlobal, INSERT_VALUES, vecLocalCheck));
249:   PetscCall(DMStagVecGetArrayRead(dm, vecLocalCheck, &arr));
250:   for (j = starty; j < starty + ny + nExtray; ++j) {
251:     for (i = startx; i < startx + nx + nExtrax; ++i) {
252:       for (c = 0; c < dof0; ++c) {
253:         const PetscScalar valRef = TEST_FUNCTION(i, j, 0, idxDownLeft, c);
254:         const PetscScalar val    = arr[j][i][idxDownLeft + c];
255:         PetscCall(CompareValues(i, j, k, c, val, valRef));
256:       }
257:       if (j < starty + ny) {
258:         for (c = 0; c < dof1; ++c) {
259:           const PetscScalar valRef = TEST_FUNCTION(i, j, 0, idxLeft, c);
260:           const PetscScalar val    = arr[j][i][idxLeft + c];
261:           PetscCall(CompareValues(i, j, k, c, val, valRef));
262:         }
263:       } else {
264:         for (c = 0; c < dof1; ++c) {
265:           const PetscScalar valRef = -1.0;
266:           const PetscScalar val    = arr[j][i][idxLeft + c];
267:           PetscCall(CompareValues(i, j, k, c, val, valRef));
268:         }
269:       }
270:       if (i < startx + nx) {
271:         for (c = 0; c < dof1; ++c) {
272:           const PetscScalar valRef = TEST_FUNCTION(i, j, 0, idxDown, c);
273:           const PetscScalar val    = arr[j][i][idxDown + c];
274:           PetscCall(CompareValues(i, j, k, c, val, valRef));
275:         }
276:       } else {
277:         for (c = 0; c < dof1; ++c) {
278:           const PetscScalar valRef = -1.0;
279:           const PetscScalar val    = arr[j][i][idxDown + c];
280:           PetscCall(CompareValues(i, j, k, c, val, valRef));
281:         }
282:       }
283:       if (i < startx + nx && j < starty + ny) {
284:         for (c = 0; c < dof2; ++c) {
285:           const PetscScalar valRef = TEST_FUNCTION(i, j, 0, idxElement, c);
286:           const PetscScalar val    = arr[j][i][idxElement + c];
287:           PetscCall(CompareValues(i, j, k, c, val, valRef));
288:         }
289:       } else {
290:         for (c = 0; c < dof2; ++c) {
291:           const PetscScalar valRef = -1.0;
292:           const PetscScalar val    = arr[j][i][idxElement + c];
293:           PetscCall(CompareValues(i, j, k, c, val, valRef));
294:         }
295:       }
296:     }
297:   }
298:   PetscCall(DMStagVecRestoreArrayRead(dm, vecLocalCheck, &arr));
299:   PetscCall(VecDestroy(&vecLocal));
300:   PetscCall(VecDestroy(&vecLocalCheck));
301:   PetscCall(VecDestroy(&vecGlobal));
302:   PetscFunctionReturn(PETSC_SUCCESS);
303: }

305: static PetscErrorCode Test2_3d(DM dm)
306: {
307:   Vec             vecLocal, vecLocalCheck, vecGlobal;
308:   PetscInt        i, j, k, startx, starty, startz, nx, ny, nz, nExtrax, nExtray, nExtraz, dof0, dof1, dof2, dof3, c, idxLeft, idxDown, idxDownLeft, idxBackDownLeft, idxBackDown, idxBack, idxBackLeft, idxElement;
309:   PetscScalar ****arr;

311:   PetscFunctionBeginUser;
312:   PetscCall(DMCreateLocalVector(dm, &vecLocal));
313:   PetscCall(VecSet(vecLocal, -1.0));
314:   PetscCall(DMStagGetCorners(dm, &startx, &starty, &startz, &nx, &ny, &nz, &nExtrax, &nExtray, &nExtraz));
315:   PetscCall(DMStagGetDOF(dm, &dof0, &dof1, &dof2, &dof3));
316:   PetscCall(DMStagVecGetArray(dm, vecLocal, &arr));
317:   if (dof0 > 0) PetscCall(DMStagGetLocationSlot(dm, DMSTAG_BACK_DOWN_LEFT, 0, &idxBackDownLeft));
318:   if (dof1 > 0) {
319:     PetscCall(DMStagGetLocationSlot(dm, DMSTAG_BACK_LEFT, 0, &idxBackLeft));
320:     PetscCall(DMStagGetLocationSlot(dm, DMSTAG_BACK_DOWN, 0, &idxBackDown));
321:     PetscCall(DMStagGetLocationSlot(dm, DMSTAG_DOWN_LEFT, 0, &idxDownLeft));
322:   }
323:   if (dof2 > 0) {
324:     PetscCall(DMStagGetLocationSlot(dm, DMSTAG_LEFT, 0, &idxLeft));
325:     PetscCall(DMStagGetLocationSlot(dm, DMSTAG_DOWN, 0, &idxDown));
326:     PetscCall(DMStagGetLocationSlot(dm, DMSTAG_BACK, 0, &idxBack));
327:   }
328:   if (dof3 > 0) PetscCall(DMStagGetLocationSlot(dm, DMSTAG_ELEMENT, 0, &idxElement));
329:   for (k = startz; k < startz + nz + nExtraz; ++k) {
330:     for (j = starty; j < starty + ny + nExtray; ++j) {
331:       for (i = startx; i < startx + nx + nExtrax; ++i) {
332:         for (c = 0; c < dof0; ++c) {
333:           const PetscScalar valRef          = TEST_FUNCTION(i, j, k, idxBackDownLeft, c);
334:           arr[k][j][i][idxBackDownLeft + c] = valRef;
335:         }
336:         if (k < startz + nz) {
337:           for (c = 0; c < dof1; ++c) {
338:             const PetscScalar valRef      = TEST_FUNCTION(i, j, k, idxDownLeft, c);
339:             arr[k][j][i][idxDownLeft + c] = valRef;
340:           }
341:         }
342:         if (j < starty + ny) {
343:           for (c = 0; c < dof1; ++c) {
344:             const PetscScalar valRef      = TEST_FUNCTION(i, j, k, idxBackLeft, c);
345:             arr[k][j][i][idxBackLeft + c] = valRef;
346:           }
347:         }
348:         if (i < startx + nx) {
349:           for (c = 0; c < dof1; ++c) {
350:             const PetscScalar valRef      = TEST_FUNCTION(i, j, k, idxBackDown, c);
351:             arr[k][j][i][idxBackDown + c] = valRef;
352:           }
353:         }
354:         if (j < starty + ny && k < startz + nz) {
355:           for (c = 0; c < dof2; ++c) {
356:             const PetscScalar valRef  = TEST_FUNCTION(i, j, k, idxLeft, c);
357:             arr[k][j][i][idxLeft + c] = valRef;
358:           }
359:         }
360:         if (i < startx + nx && k < startz + nz) {
361:           for (c = 0; c < dof2; ++c) {
362:             const PetscScalar valRef  = TEST_FUNCTION(i, j, k, idxDown, c);
363:             arr[k][j][i][idxDown + c] = valRef;
364:           }
365:         }
366:         if (i < startx + nx && j < starty + ny) {
367:           for (c = 0; c < dof2; ++c) {
368:             const PetscScalar valRef  = TEST_FUNCTION(i, j, k, idxBack, c);
369:             arr[k][j][i][idxBack + c] = valRef;
370:           }
371:         }
372:         if (i < startx + nx && j < starty + ny && k < startz + nz) {
373:           for (c = 0; c < dof3; ++c) {
374:             const PetscScalar valRef     = TEST_FUNCTION(i, j, k, idxElement, c);
375:             arr[k][j][i][idxElement + c] = valRef;
376:           }
377:         }
378:       }
379:     }
380:   }
381:   PetscCall(DMStagVecRestoreArray(dm, vecLocal, &arr));
382:   PetscCall(DMCreateGlobalVector(dm, &vecGlobal));
383:   PetscCall(DMLocalToGlobal(dm, vecLocal, INSERT_VALUES, vecGlobal));
384:   PetscCall(VecDuplicate(vecLocal, &vecLocalCheck));
385:   PetscCall(VecSet(vecLocalCheck, -1.0));
386:   PetscCall(DMGlobalToLocal(dm, vecGlobal, INSERT_VALUES, vecLocalCheck));
387:   PetscCall(DMStagVecGetArrayRead(dm, vecLocalCheck, &arr));
388:   for (k = startz; k < startz + nz + nExtraz; ++k) {
389:     for (j = starty; j < starty + ny + nExtray; ++j) {
390:       for (i = startx; i < startx + nx + nExtrax; ++i) {
391:         for (c = 0; c < dof0; ++c) {
392:           const PetscScalar valRef = TEST_FUNCTION(i, j, k, idxBackDownLeft, c);
393:           const PetscScalar val    = arr[k][j][i][idxBackDownLeft + c];
394:           PetscCall(CompareValues(i, j, k, c, val, valRef));
395:         }
396:         if (k < startz + nz) {
397:           for (c = 0; c < dof1; ++c) {
398:             const PetscScalar valRef = TEST_FUNCTION(i, j, k, idxDownLeft, c);
399:             const PetscScalar val    = arr[k][j][i][idxDownLeft + c];
400:             PetscCall(CompareValues(i, j, k, c, val, valRef));
401:           }
402:         } else {
403:           for (c = 0; c < dof1; ++c) {
404:             const PetscScalar valRef = -1.0;
405:             const PetscScalar val    = arr[k][j][i][idxDownLeft + c];
406:             PetscCall(CompareValues(i, j, k, c, val, valRef));
407:           }
408:         }
409:         if (j < starty + ny) {
410:           for (c = 0; c < dof1; ++c) {
411:             const PetscScalar valRef = TEST_FUNCTION(i, j, k, idxBackLeft, c);
412:             const PetscScalar val    = arr[k][j][i][idxBackLeft + c];
413:             PetscCall(CompareValues(i, j, k, c, val, valRef));
414:           }
415:         } else {
416:           for (c = 0; c < dof1; ++c) {
417:             const PetscScalar valRef = -1.0;
418:             const PetscScalar val    = arr[k][j][i][idxBackLeft + c];
419:             PetscCall(CompareValues(i, j, k, c, val, valRef));
420:           }
421:         }
422:         if (i < startx + nx) {
423:           for (c = 0; c < dof1; ++c) {
424:             const PetscScalar valRef = TEST_FUNCTION(i, j, k, idxBackDown, c);
425:             const PetscScalar val    = arr[k][j][i][idxBackDown + c];
426:             PetscCall(CompareValues(i, j, k, c, val, valRef));
427:           }
428:         } else {
429:           for (c = 0; c < dof1; ++c) {
430:             const PetscScalar valRef = -1.0;
431:             const PetscScalar val    = arr[k][j][i][idxBackDown + c];
432:             PetscCall(CompareValues(i, j, k, c, val, valRef));
433:           }
434:         }
435:         if (j < starty + ny && k < startz + nz) {
436:           for (c = 0; c < dof2; ++c) {
437:             const PetscScalar valRef = TEST_FUNCTION(i, j, k, idxLeft, c);
438:             const PetscScalar val    = arr[k][j][i][idxLeft + c];
439:             PetscCall(CompareValues(i, j, k, c, val, valRef));
440:           }
441:         } else {
442:           for (c = 0; c < dof2; ++c) {
443:             const PetscScalar valRef = -1.0;
444:             const PetscScalar val    = arr[k][j][i][idxLeft + c];
445:             PetscCall(CompareValues(i, j, k, c, val, valRef));
446:           }
447:         }
448:         if (i < startx + nx && k < startz + nz) {
449:           for (c = 0; c < dof2; ++c) {
450:             const PetscScalar valRef = TEST_FUNCTION(i, j, k, idxDown, c);
451:             const PetscScalar val    = arr[k][j][i][idxDown + c];
452:             PetscCall(CompareValues(i, j, k, c, val, valRef));
453:           }
454:         } else {
455:           for (c = 0; c < dof2; ++c) {
456:             const PetscScalar valRef = -1.0;
457:             const PetscScalar val    = arr[k][j][i][idxDown + c];
458:             PetscCall(CompareValues(i, j, k, c, val, valRef));
459:           }
460:         }
461:         if (i < startx + nx && j < starty + ny) {
462:           for (c = 0; c < dof2; ++c) {
463:             const PetscScalar valRef = TEST_FUNCTION(i, j, k, idxBack, c);
464:             const PetscScalar val    = arr[k][j][i][idxBack + c];
465:             PetscCall(CompareValues(i, j, k, c, val, valRef));
466:           }
467:         } else {
468:           for (c = 0; c < dof2; ++c) {
469:             const PetscScalar valRef = -1.0;
470:             const PetscScalar val    = arr[k][j][i][idxBack + c];
471:             PetscCall(CompareValues(i, j, k, c, val, valRef));
472:           }
473:         }
474:         if (i < startx + nx && j < starty + ny && k < startz + nz) {
475:           for (c = 0; c < dof3; ++c) {
476:             const PetscScalar valRef = TEST_FUNCTION(i, j, k, idxElement, c);
477:             const PetscScalar val    = arr[k][j][i][idxElement + c];
478:             PetscCall(CompareValues(i, j, k, c, val, valRef));
479:           }
480:         } else {
481:           for (c = 0; c < dof3; ++c) {
482:             const PetscScalar valRef = -1.0;
483:             const PetscScalar val    = arr[k][j][i][idxElement + c];
484:             PetscCall(CompareValues(i, j, k, c, val, valRef));
485:           }
486:         }
487:       }
488:     }
489:   }
490:   PetscCall(DMStagVecRestoreArrayRead(dm, vecLocalCheck, &arr));
491:   PetscCall(VecDestroy(&vecLocal));
492:   PetscCall(VecDestroy(&vecLocalCheck));
493:   PetscCall(VecDestroy(&vecGlobal));
494:   PetscFunctionReturn(PETSC_SUCCESS);
495: }
496: #undef TEST_FUNCTION

498: /*TEST

500:    testset:
501:       suffix: periodic_1d_seq
502:       nsize: 1
503:       args: -dm_view -dim 1 -stag_grid_x 4 -stag_boundary_type_x periodic -stag_stencil_width {{0 1 2}separate output}

505:    testset:
506:       suffix: ghosted_1d_seq
507:       nsize: 1
508:       args: -dm_view -dim 1 -stag_grid_x 4 -stag_boundary_type_x ghosted -stag_stencil_width {{0 1 2}separate output}

510:    testset:
511:       suffix: none_1d_seq
512:       nsize: 1
513:       args: -dm_view -dim 1 -stag_grid_x 4 -stag_boundary_type_x ghosted -stag_stencil_width {{0 1 2}separate output}

515:    testset:
516:       suffix: periodic_1d_par
517:       nsize: 3
518:       args: -dm_view -dim 1 -setsizes -stag_boundary_type_x periodic -stag_stencil_width {{0 1 2}separate output}

520:    testset:
521:       suffix: ghosted_1d_par
522:       nsize: 3
523:       args: -dm_view -dim 1 -setsizes -stag_boundary_type_x ghosted -stag_stencil_width {{0 1 2}separate output}

525:    testset:
526:       suffix: none_1d_par
527:       nsize: 3
528:       args: -dm_view -dim 1 -setsizes -stag_boundary_type_x ghosted -stag_stencil_width {{0 1 2}separate output}

530:    testset:
531:       suffix: periodic_periodic_2d_seq
532:       nsize: 1
533:       args: -dm_view -dim 2 -stag_grid_x 4 -stag_grid_y 5 -stag_boundary_type_x periodic -stag_boundary_type_y periodic -stag_stencil_width {{0 1 2}separate output}

535:    testset:
536:       suffix: periodic_ghosted_2d_seq
537:       nsize: 1
538:       args: -dm_view -dim 2 -stag_grid_x 4 -stag_grid_y 5 -stag_boundary_type_x periodic -stag_boundary_type_y ghosted -stag_stencil_width {{0 1 2}separate output}

540:    testset:
541:       suffix: none_none_2d_seq
542:       nsize: 1
543:       args: -dm_view -dim 2 -stag_grid_x 4 -stag_grid_y 5 -stag_boundary_type_x none -stag_boundary_type_y none -stag_stencil_width {{0 1 2}separate output}

545:    testset:
546:       suffix: none_ghosted_2d_seq
547:       nsize: 1
548:       args: -dm_view -dim 2 -stag_grid_x 4 -stag_grid_y 5 -stag_boundary_type_x none -stag_boundary_type_y ghosted -stag_stencil_width {{0 1 2}separate output}

550:    testset:
551:       suffix: none_periodic_2d_seq
552:       nsize: 1
553:       args: -dm_view -dim 2 -stag_grid_x 4 -stag_grid_y 5 -stag_boundary_type_x none -stag_boundary_type_y periodic -stag_stencil_width {{0 1 2}separate output}

555:    testset:
556:       suffix: periodic_periodic_2d_par
557:       nsize: 6
558:       args: -dm_view -dim 2 -setsizes -stag_boundary_type_x periodic -stag_boundary_type_y periodic -stag_stencil_width {{0 1 2}separate output}

560:    testset:
561:       suffix: periodic_ghosted_2d_par
562:       nsize: 6
563:       args: -dm_view -dim 2 -setsizes -stag_boundary_type_x periodic -stag_boundary_type_y ghosted -stag_stencil_width {{0 1 2}separate output}

565:    testset:
566:       suffix: none_none_2d_par
567:       nsize: 6
568:       args: -dm_view -dim 2 -setsizes -stag_boundary_type_x none -stag_boundary_type_y none -stag_stencil_width {{0 1 2}separate output}

570:    testset:
571:       suffix: none_ghosted_2d_par
572:       nsize: 6
573:       args: -dm_view -dim 2 -setsizes -stag_boundary_type_x none -stag_boundary_type_y ghosted -stag_stencil_width {{0 1 2}separate output}

575:    testset:
576:       suffix: none_periodic_2d_par
577:       nsize: 6
578:       args: -dm_view -dim 2 -setsizes -stag_boundary_type_x none -stag_boundary_type_y periodic -stag_stencil_width {{0 1 2}separate output}

580:    testset:
581:       suffix: periodic_periodic_periodic_3d_seq
582:       nsize: 1
583:       args: -dm_view -dim 3 -stag_grid_x 4 -stag_grid_y 5 -stag_grid_z 3 -stag_boundary_type_x periodic -stag_boundary_type_y periodic -stag_boundary_type_z periodic -stag_stencil_width {{0 1 2}separate output}

585:    testset:
586:       suffix: periodic_ghosted_periodic_3d_seq
587:       nsize: 1
588:       args: -dm_view -dim 3 -stag_grid_x 4 -stag_grid_y 5 -stag_grid_z 3 -stag_boundary_type_x periodic -stag_boundary_type_y ghosted -stag_boundary_type_z periodic -stag_stencil_width {{0 1 2}separate output}

590:    testset:
591:       suffix: none_periodic_ghosted_3d_seq
592:       nsize: 1
593:       args: -dm_view -dim 3 -stag_grid_x 4 -stag_grid_y 5 -stag_grid_z 3 -stag_boundary_type_x none -stag_boundary_type_y periodic -stag_boundary_type_z ghosted -stag_stencil_width {{0 1 2}separate output}

595:    testset:
596:       suffix: none_none_none_3d_seq
597:       nsize: 1
598:       args: -dm_view -dim 3 -stag_grid_x 4 -stag_grid_y 5 -stag_grid_z 3 -stag_boundary_type_x none -stag_boundary_type_y none -stag_boundary_type_z none -stag_stencil_width {{0 1 2}separate output}

600:    testset:
601:       suffix: periodic_periodic_periodic_3d_par
602:       nsize: 12
603:       args: -dm_view -dim 3 -setsizes -stag_boundary_type_x periodic -stag_boundary_type_y periodic -stag_boundary_type_z periodic -stag_stencil_width {{0 1 2}separate output}

605:    testset:
606:       suffix: periodic_ghosted_ghosted_3d_par
607:       nsize: 12
608:       args: -dm_view -dim 3 -setsizes -stag_boundary_type_x periodic -stag_boundary_type_y ghosted -stag_boundary_type_z ghosted -stag_stencil_width {{0 1 2}separate output}

610:    testset:
611:       suffix: ghosted_periodic_periodic_3d_par
612:       nsize: 12
613:       args: -dm_view -dim 3 -setsizes -stag_boundary_type_x ghosted -stag_boundary_type_y periodic -stag_boundary_type_z periodic -stag_stencil_width {{0 1 2}separate output}

615:    testset:
616:       suffix: none_none_none_3d_par
617:       nsize: 12
618:       args: -dm_view -dim 3 -setsizes -stag_boundary_type_x none -stag_boundary_type_y none -stag_boundary_type_z none -stag_stencil_width {{0 1 2}separate output}

620:    test:
621:       suffix: periodic_none_none_3d_skinny_seq
622:       nsize: 1
623:       args: -dm_view -dim 3 -stag_boundary_type_x periodic -stag_boundary_type_y none -stag_boundary_type_z none -stag_grid_x 3 -stag_grid_y 6 -stag_grid_z 5 -stag_stencil_width 1 -useinjective 0

625:    test:
626:       suffix: periodic_none_none_3d_skinny_par
627:       nsize: 4
628:       args: -dm_view -dim 3 -stag_boundary_type_x periodic -stag_boundary_type_y none -stag_boundary_type_z none -stag_grid_x 3 -stag_grid_y 6 -stag_grid_z 5 -stag_ranks_x 1 -stag_ranks_y 2 -stag_ranks_z 2 -stag_stencil_width 1 -useinjective 0

630: TEST*/