xref: /petsc/src/dm/impls/stag/tests/ex13.c (revision 834855d6effb0d027771461c8e947ee1ce5a1e17) 
1 static char help[] = "Test DMStagPopulateLocalToGlobalInjective.\n\n";
2 
3 #include <petscdm.h>
4 #include <petscdmstag.h>
5 
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);
10 
main(int argc,char ** argv)11 int main(int argc, char **argv)
12 {
13   DM        dm;
14   PetscInt  dim;
15   PetscBool setSizes, useInjective;
16 
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));
26 
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;
48 
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   }
67 
68   /* Setup */
69   PetscCall(DMSetFromOptions(dm));
70   PetscCall(DMSetUp(dm));
71 
72   /* Populate Additional Injective Local-to-Global Map */
73   if (useInjective) PetscCall(DMStagPopulateLocalToGlobalInjective(dm));
74 
75   /* Test: Make sure L2G inverts G2L */
76   PetscCall(Test1(dm));
77 
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   }
93 
94   /* Clean up and finalize PETSc */
95   PetscCall(DMDestroy(&dm));
96   PetscCall(PetscFinalize());
97   return 0;
98 }
99 
Test1(DM dm)100 static PetscErrorCode Test1(DM dm)
101 {
102   Vec         vecLocal, vecGlobal, vecGlobalCheck;
103   PetscRandom rctx;
104   PetscBool   equal;
105 
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 }
123 
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)
126 
127 /* Helper function to check */
CompareValues(PetscInt i,PetscInt j,PetscInt k,PetscInt c,PetscScalar val,PetscScalar valRef)128 static PetscErrorCode CompareValues(PetscInt i, PetscInt j, PetscInt k, PetscInt c, PetscScalar val, PetscScalar valRef)
129 {
130   PetscFunctionBeginUser;
131   PetscCheck(val == valRef || PetscAbsScalar(val - valRef) / PetscAbsScalar(valRef) <= 10 * PETSC_MACHINE_EPSILON, 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));
132   PetscFunctionReturn(PETSC_SUCCESS);
133 }
134 
Test2_1d(DM dm)135 static PetscErrorCode Test2_1d(DM dm)
136 {
137   Vec            vecLocal, vecLocalCheck, vecGlobal;
138   PetscInt       i, startx, nx, nExtrax, dof0, dof1, c, idxLeft, idxElement;
139   PetscScalar  **arr;
140   const PetscInt j = -1, k = -1;
141 
142   PetscFunctionBeginUser;
143   PetscCall(DMCreateLocalVector(dm, &vecLocal));
144   PetscCall(VecSet(vecLocal, -1.0));
145   PetscCall(DMStagGetCorners(dm, &startx, NULL, NULL, &nx, NULL, NULL, &nExtrax, NULL, NULL));
146   PetscCall(DMStagGetDOF(dm, &dof0, &dof1, NULL, NULL));
147   PetscCall(DMStagVecGetArray(dm, vecLocal, &arr));
148   if (dof0 > 0) PetscCall(DMStagGetLocationSlot(dm, DMSTAG_LEFT, 0, &idxLeft));
149   if (dof1 > 0) PetscCall(DMStagGetLocationSlot(dm, DMSTAG_ELEMENT, 0, &idxElement));
150   for (i = startx; i < startx + nx + nExtrax; ++i) {
151     for (c = 0; c < dof0; ++c) {
152       const PetscScalar valRef = TEST_FUNCTION(i, 0, 0, idxLeft, c);
153       arr[i][idxLeft + c]      = valRef;
154     }
155     if (i < startx + nx) {
156       for (c = 0; c < dof1; ++c) {
157         const PetscScalar valRef = TEST_FUNCTION(i, 0, 0, idxElement, c);
158         arr[i][idxElement + c]   = valRef;
159       }
160     }
161   }
162   PetscCall(DMStagVecRestoreArray(dm, vecLocal, &arr));
163   PetscCall(DMCreateGlobalVector(dm, &vecGlobal));
164   PetscCall(DMLocalToGlobal(dm, vecLocal, INSERT_VALUES, vecGlobal));
165   PetscCall(VecDuplicate(vecLocal, &vecLocalCheck));
166   PetscCall(VecSet(vecLocalCheck, -1.0));
167   PetscCall(DMGlobalToLocal(dm, vecGlobal, INSERT_VALUES, vecLocalCheck));
168   PetscCall(DMStagVecGetArrayRead(dm, vecLocalCheck, &arr));
169   for (i = startx; i < startx + nx + nExtrax; ++i) {
170     for (c = 0; c < dof0; ++c) {
171       const PetscScalar valRef = TEST_FUNCTION(i, 0, 0, idxLeft, c);
172       const PetscScalar val    = arr[i][idxLeft + c];
173       PetscCall(CompareValues(i, j, k, c, val, valRef));
174     }
175     if (i < startx + nx) {
176       for (c = 0; c < dof1; ++c) {
177         const PetscScalar valRef = TEST_FUNCTION(i, 0, 0, idxElement, c);
178         const PetscScalar val    = arr[i][idxElement + c];
179         PetscCall(CompareValues(i, j, k, c, val, valRef));
180       }
181     } else {
182       for (c = 0; c < dof1; ++c) {
183         const PetscScalar valRef = -1.0;
184         const PetscScalar val    = arr[i][idxElement + c];
185         PetscCall(CompareValues(i, j, k, c, val, valRef));
186       }
187     }
188   }
189   PetscCall(DMStagVecRestoreArrayRead(dm, vecLocalCheck, &arr));
190   PetscCall(VecDestroy(&vecLocal));
191   PetscCall(VecDestroy(&vecLocalCheck));
192   PetscCall(VecDestroy(&vecGlobal));
193   PetscFunctionReturn(PETSC_SUCCESS);
194 }
195 
Test2_2d(DM dm)196 static PetscErrorCode Test2_2d(DM dm)
197 {
198   Vec            vecLocal, vecLocalCheck, vecGlobal;
199   PetscInt       i, j, startx, starty, nx, ny, nExtrax, nExtray, dof0, dof1, dof2, c, idxLeft, idxDown, idxDownLeft, idxElement;
200   PetscScalar ***arr;
201   const PetscInt k = -1;
202 
203   PetscFunctionBeginUser;
204   PetscCall(DMCreateLocalVector(dm, &vecLocal));
205   PetscCall(VecSet(vecLocal, -1.0));
206   PetscCall(DMStagGetCorners(dm, &startx, &starty, NULL, &nx, &ny, NULL, &nExtrax, &nExtray, NULL));
207   PetscCall(DMStagGetDOF(dm, &dof0, &dof1, &dof2, NULL));
208   PetscCall(DMStagVecGetArray(dm, vecLocal, &arr));
209   if (dof0 > 0) PetscCall(DMStagGetLocationSlot(dm, DMSTAG_DOWN_LEFT, 0, &idxDownLeft));
210   if (dof1 > 0) {
211     PetscCall(DMStagGetLocationSlot(dm, DMSTAG_LEFT, 0, &idxLeft));
212     PetscCall(DMStagGetLocationSlot(dm, DMSTAG_DOWN, 0, &idxDown));
213   }
214   if (dof2 > 0) PetscCall(DMStagGetLocationSlot(dm, DMSTAG_ELEMENT, 0, &idxElement));
215   for (j = starty; j < starty + ny + nExtray; ++j) {
216     for (i = startx; i < startx + nx + nExtrax; ++i) {
217       for (c = 0; c < dof0; ++c) {
218         const PetscScalar valRef   = TEST_FUNCTION(i, j, 0, idxDownLeft, c);
219         arr[j][i][idxDownLeft + c] = valRef;
220       }
221       if (j < starty + ny) {
222         for (c = 0; c < dof1; ++c) {
223           const PetscScalar valRef = TEST_FUNCTION(i, j, 0, idxLeft, c);
224           arr[j][i][idxLeft + c]   = valRef;
225         }
226       }
227       if (i < startx + nx) {
228         for (c = 0; c < dof1; ++c) {
229           const PetscScalar valRef = TEST_FUNCTION(i, j, 0, idxDown, c);
230           arr[j][i][idxDown + c]   = valRef;
231         }
232       }
233       if (i < startx + nx && j < starty + ny) {
234         for (c = 0; c < dof2; ++c) {
235           const PetscScalar valRef  = TEST_FUNCTION(i, j, 0, idxElement, c);
236           arr[j][i][idxElement + c] = valRef;
237         }
238       }
239     }
240   }
241   PetscCall(DMStagVecRestoreArray(dm, vecLocal, &arr));
242   PetscCall(DMCreateGlobalVector(dm, &vecGlobal));
243   PetscCall(DMLocalToGlobal(dm, vecLocal, INSERT_VALUES, vecGlobal));
244   PetscCall(VecDuplicate(vecLocal, &vecLocalCheck));
245   PetscCall(VecSet(vecLocalCheck, -1.0));
246   PetscCall(DMGlobalToLocal(dm, vecGlobal, INSERT_VALUES, vecLocalCheck));
247   PetscCall(DMStagVecGetArrayRead(dm, vecLocalCheck, &arr));
248   for (j = starty; j < starty + ny + nExtray; ++j) {
249     for (i = startx; i < startx + nx + nExtrax; ++i) {
250       for (c = 0; c < dof0; ++c) {
251         const PetscScalar valRef = TEST_FUNCTION(i, j, 0, idxDownLeft, c);
252         const PetscScalar val    = arr[j][i][idxDownLeft + c];
253         PetscCall(CompareValues(i, j, k, c, val, valRef));
254       }
255       if (j < starty + ny) {
256         for (c = 0; c < dof1; ++c) {
257           const PetscScalar valRef = TEST_FUNCTION(i, j, 0, idxLeft, c);
258           const PetscScalar val    = arr[j][i][idxLeft + c];
259           PetscCall(CompareValues(i, j, k, c, val, valRef));
260         }
261       } else {
262         for (c = 0; c < dof1; ++c) {
263           const PetscScalar valRef = -1.0;
264           const PetscScalar val    = arr[j][i][idxLeft + c];
265           PetscCall(CompareValues(i, j, k, c, val, valRef));
266         }
267       }
268       if (i < startx + nx) {
269         for (c = 0; c < dof1; ++c) {
270           const PetscScalar valRef = TEST_FUNCTION(i, j, 0, idxDown, c);
271           const PetscScalar val    = arr[j][i][idxDown + c];
272           PetscCall(CompareValues(i, j, k, c, val, valRef));
273         }
274       } else {
275         for (c = 0; c < dof1; ++c) {
276           const PetscScalar valRef = -1.0;
277           const PetscScalar val    = arr[j][i][idxDown + c];
278           PetscCall(CompareValues(i, j, k, c, val, valRef));
279         }
280       }
281       if (i < startx + nx && j < starty + ny) {
282         for (c = 0; c < dof2; ++c) {
283           const PetscScalar valRef = TEST_FUNCTION(i, j, 0, idxElement, c);
284           const PetscScalar val    = arr[j][i][idxElement + c];
285           PetscCall(CompareValues(i, j, k, c, val, valRef));
286         }
287       } else {
288         for (c = 0; c < dof2; ++c) {
289           const PetscScalar valRef = -1.0;
290           const PetscScalar val    = arr[j][i][idxElement + c];
291           PetscCall(CompareValues(i, j, k, c, val, valRef));
292         }
293       }
294     }
295   }
296   PetscCall(DMStagVecRestoreArrayRead(dm, vecLocalCheck, &arr));
297   PetscCall(VecDestroy(&vecLocal));
298   PetscCall(VecDestroy(&vecLocalCheck));
299   PetscCall(VecDestroy(&vecGlobal));
300   PetscFunctionReturn(PETSC_SUCCESS);
301 }
302 
Test2_3d(DM dm)303 static PetscErrorCode Test2_3d(DM dm)
304 {
305   Vec             vecLocal, vecLocalCheck, vecGlobal;
306   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;
307   PetscScalar ****arr;
308 
309   PetscFunctionBeginUser;
310   PetscCall(DMCreateLocalVector(dm, &vecLocal));
311   PetscCall(VecSet(vecLocal, -1.0));
312   PetscCall(DMStagGetCorners(dm, &startx, &starty, &startz, &nx, &ny, &nz, &nExtrax, &nExtray, &nExtraz));
313   PetscCall(DMStagGetDOF(dm, &dof0, &dof1, &dof2, &dof3));
314   PetscCall(DMStagVecGetArray(dm, vecLocal, &arr));
315   if (dof0 > 0) PetscCall(DMStagGetLocationSlot(dm, DMSTAG_BACK_DOWN_LEFT, 0, &idxBackDownLeft));
316   if (dof1 > 0) {
317     PetscCall(DMStagGetLocationSlot(dm, DMSTAG_BACK_LEFT, 0, &idxBackLeft));
318     PetscCall(DMStagGetLocationSlot(dm, DMSTAG_BACK_DOWN, 0, &idxBackDown));
319     PetscCall(DMStagGetLocationSlot(dm, DMSTAG_DOWN_LEFT, 0, &idxDownLeft));
320   }
321   if (dof2 > 0) {
322     PetscCall(DMStagGetLocationSlot(dm, DMSTAG_LEFT, 0, &idxLeft));
323     PetscCall(DMStagGetLocationSlot(dm, DMSTAG_DOWN, 0, &idxDown));
324     PetscCall(DMStagGetLocationSlot(dm, DMSTAG_BACK, 0, &idxBack));
325   }
326   if (dof3 > 0) PetscCall(DMStagGetLocationSlot(dm, DMSTAG_ELEMENT, 0, &idxElement));
327   for (k = startz; k < startz + nz + nExtraz; ++k) {
328     for (j = starty; j < starty + ny + nExtray; ++j) {
329       for (i = startx; i < startx + nx + nExtrax; ++i) {
330         for (c = 0; c < dof0; ++c) {
331           const PetscScalar valRef          = TEST_FUNCTION(i, j, k, idxBackDownLeft, c);
332           arr[k][j][i][idxBackDownLeft + c] = valRef;
333         }
334         if (k < startz + nz) {
335           for (c = 0; c < dof1; ++c) {
336             const PetscScalar valRef      = TEST_FUNCTION(i, j, k, idxDownLeft, c);
337             arr[k][j][i][idxDownLeft + c] = valRef;
338           }
339         }
340         if (j < starty + ny) {
341           for (c = 0; c < dof1; ++c) {
342             const PetscScalar valRef      = TEST_FUNCTION(i, j, k, idxBackLeft, c);
343             arr[k][j][i][idxBackLeft + c] = valRef;
344           }
345         }
346         if (i < startx + nx) {
347           for (c = 0; c < dof1; ++c) {
348             const PetscScalar valRef      = TEST_FUNCTION(i, j, k, idxBackDown, c);
349             arr[k][j][i][idxBackDown + c] = valRef;
350           }
351         }
352         if (j < starty + ny && k < startz + nz) {
353           for (c = 0; c < dof2; ++c) {
354             const PetscScalar valRef  = TEST_FUNCTION(i, j, k, idxLeft, c);
355             arr[k][j][i][idxLeft + c] = valRef;
356           }
357         }
358         if (i < startx + nx && k < startz + nz) {
359           for (c = 0; c < dof2; ++c) {
360             const PetscScalar valRef  = TEST_FUNCTION(i, j, k, idxDown, c);
361             arr[k][j][i][idxDown + c] = valRef;
362           }
363         }
364         if (i < startx + nx && j < starty + ny) {
365           for (c = 0; c < dof2; ++c) {
366             const PetscScalar valRef  = TEST_FUNCTION(i, j, k, idxBack, c);
367             arr[k][j][i][idxBack + c] = valRef;
368           }
369         }
370         if (i < startx + nx && j < starty + ny && k < startz + nz) {
371           for (c = 0; c < dof3; ++c) {
372             const PetscScalar valRef     = TEST_FUNCTION(i, j, k, idxElement, c);
373             arr[k][j][i][idxElement + c] = valRef;
374           }
375         }
376       }
377     }
378   }
379   PetscCall(DMStagVecRestoreArray(dm, vecLocal, &arr));
380   PetscCall(DMCreateGlobalVector(dm, &vecGlobal));
381   PetscCall(DMLocalToGlobal(dm, vecLocal, INSERT_VALUES, vecGlobal));
382   PetscCall(VecDuplicate(vecLocal, &vecLocalCheck));
383   PetscCall(VecSet(vecLocalCheck, -1.0));
384   PetscCall(DMGlobalToLocal(dm, vecGlobal, INSERT_VALUES, vecLocalCheck));
385   PetscCall(DMStagVecGetArrayRead(dm, vecLocalCheck, &arr));
386   for (k = startz; k < startz + nz + nExtraz; ++k) {
387     for (j = starty; j < starty + ny + nExtray; ++j) {
388       for (i = startx; i < startx + nx + nExtrax; ++i) {
389         for (c = 0; c < dof0; ++c) {
390           const PetscScalar valRef = TEST_FUNCTION(i, j, k, idxBackDownLeft, c);
391           const PetscScalar val    = arr[k][j][i][idxBackDownLeft + c];
392           PetscCall(CompareValues(i, j, k, c, val, valRef));
393         }
394         if (k < startz + nz) {
395           for (c = 0; c < dof1; ++c) {
396             const PetscScalar valRef = TEST_FUNCTION(i, j, k, idxDownLeft, c);
397             const PetscScalar val    = arr[k][j][i][idxDownLeft + c];
398             PetscCall(CompareValues(i, j, k, c, val, valRef));
399           }
400         } else {
401           for (c = 0; c < dof1; ++c) {
402             const PetscScalar valRef = -1.0;
403             const PetscScalar val    = arr[k][j][i][idxDownLeft + c];
404             PetscCall(CompareValues(i, j, k, c, val, valRef));
405           }
406         }
407         if (j < starty + ny) {
408           for (c = 0; c < dof1; ++c) {
409             const PetscScalar valRef = TEST_FUNCTION(i, j, k, idxBackLeft, c);
410             const PetscScalar val    = arr[k][j][i][idxBackLeft + c];
411             PetscCall(CompareValues(i, j, k, c, val, valRef));
412           }
413         } else {
414           for (c = 0; c < dof1; ++c) {
415             const PetscScalar valRef = -1.0;
416             const PetscScalar val    = arr[k][j][i][idxBackLeft + c];
417             PetscCall(CompareValues(i, j, k, c, val, valRef));
418           }
419         }
420         if (i < startx + nx) {
421           for (c = 0; c < dof1; ++c) {
422             const PetscScalar valRef = TEST_FUNCTION(i, j, k, idxBackDown, c);
423             const PetscScalar val    = arr[k][j][i][idxBackDown + c];
424             PetscCall(CompareValues(i, j, k, c, val, valRef));
425           }
426         } else {
427           for (c = 0; c < dof1; ++c) {
428             const PetscScalar valRef = -1.0;
429             const PetscScalar val    = arr[k][j][i][idxBackDown + c];
430             PetscCall(CompareValues(i, j, k, c, val, valRef));
431           }
432         }
433         if (j < starty + ny && k < startz + nz) {
434           for (c = 0; c < dof2; ++c) {
435             const PetscScalar valRef = TEST_FUNCTION(i, j, k, idxLeft, c);
436             const PetscScalar val    = arr[k][j][i][idxLeft + c];
437             PetscCall(CompareValues(i, j, k, c, val, valRef));
438           }
439         } else {
440           for (c = 0; c < dof2; ++c) {
441             const PetscScalar valRef = -1.0;
442             const PetscScalar val    = arr[k][j][i][idxLeft + c];
443             PetscCall(CompareValues(i, j, k, c, val, valRef));
444           }
445         }
446         if (i < startx + nx && k < startz + nz) {
447           for (c = 0; c < dof2; ++c) {
448             const PetscScalar valRef = TEST_FUNCTION(i, j, k, idxDown, c);
449             const PetscScalar val    = arr[k][j][i][idxDown + c];
450             PetscCall(CompareValues(i, j, k, c, val, valRef));
451           }
452         } else {
453           for (c = 0; c < dof2; ++c) {
454             const PetscScalar valRef = -1.0;
455             const PetscScalar val    = arr[k][j][i][idxDown + c];
456             PetscCall(CompareValues(i, j, k, c, val, valRef));
457           }
458         }
459         if (i < startx + nx && j < starty + ny) {
460           for (c = 0; c < dof2; ++c) {
461             const PetscScalar valRef = TEST_FUNCTION(i, j, k, idxBack, c);
462             const PetscScalar val    = arr[k][j][i][idxBack + c];
463             PetscCall(CompareValues(i, j, k, c, val, valRef));
464           }
465         } else {
466           for (c = 0; c < dof2; ++c) {
467             const PetscScalar valRef = -1.0;
468             const PetscScalar val    = arr[k][j][i][idxBack + c];
469             PetscCall(CompareValues(i, j, k, c, val, valRef));
470           }
471         }
472         if (i < startx + nx && j < starty + ny && k < startz + nz) {
473           for (c = 0; c < dof3; ++c) {
474             const PetscScalar valRef = TEST_FUNCTION(i, j, k, idxElement, c);
475             const PetscScalar val    = arr[k][j][i][idxElement + c];
476             PetscCall(CompareValues(i, j, k, c, val, valRef));
477           }
478         } else {
479           for (c = 0; c < dof3; ++c) {
480             const PetscScalar valRef = -1.0;
481             const PetscScalar val    = arr[k][j][i][idxElement + c];
482             PetscCall(CompareValues(i, j, k, c, val, valRef));
483           }
484         }
485       }
486     }
487   }
488   PetscCall(DMStagVecRestoreArrayRead(dm, vecLocalCheck, &arr));
489   PetscCall(VecDestroy(&vecLocal));
490   PetscCall(VecDestroy(&vecLocalCheck));
491   PetscCall(VecDestroy(&vecGlobal));
492   PetscFunctionReturn(PETSC_SUCCESS);
493 }
494 #undef TEST_FUNCTION
495 
496 /*TEST
497 
498    testset:
499       suffix: periodic_1d_seq
500       nsize: 1
501       args: -dm_view -dim 1 -stag_grid_x 4 -stag_boundary_type_x periodic -stag_stencil_width {{0 1 2}separate output}
502 
503    testset:
504       suffix: ghosted_1d_seq
505       nsize: 1
506       args: -dm_view -dim 1 -stag_grid_x 4 -stag_boundary_type_x ghosted -stag_stencil_width {{0 1 2}separate output}
507 
508    testset:
509       suffix: none_1d_seq
510       nsize: 1
511       args: -dm_view -dim 1 -stag_grid_x 4 -stag_boundary_type_x ghosted -stag_stencil_width {{0 1 2}separate output}
512 
513    testset:
514       suffix: periodic_1d_par
515       nsize: 3
516       args: -dm_view -dim 1 -setsizes -stag_boundary_type_x periodic -stag_stencil_width {{0 1 2}separate output}
517 
518    testset:
519       suffix: ghosted_1d_par
520       nsize: 3
521       args: -dm_view -dim 1 -setsizes -stag_boundary_type_x ghosted -stag_stencil_width {{0 1 2}separate output}
522 
523    testset:
524       suffix: none_1d_par
525       nsize: 3
526       args: -dm_view -dim 1 -setsizes -stag_boundary_type_x ghosted -stag_stencil_width {{0 1 2}separate output}
527 
528    testset:
529       suffix: periodic_periodic_2d_seq
530       nsize: 1
531       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}
532 
533    testset:
534       suffix: periodic_ghosted_2d_seq
535       nsize: 1
536       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}
537 
538    testset:
539       suffix: none_none_2d_seq
540       nsize: 1
541       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}
542 
543    testset:
544       suffix: none_ghosted_2d_seq
545       nsize: 1
546       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}
547 
548    testset:
549       suffix: none_periodic_2d_seq
550       nsize: 1
551       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}
552 
553    testset:
554       suffix: periodic_periodic_2d_par
555       nsize: 6
556       args: -dm_view -dim 2 -setsizes -stag_boundary_type_x periodic -stag_boundary_type_y periodic -stag_stencil_width {{0 1 2}separate output}
557 
558    testset:
559       suffix: periodic_ghosted_2d_par
560       nsize: 6
561       args: -dm_view -dim 2 -setsizes -stag_boundary_type_x periodic -stag_boundary_type_y ghosted -stag_stencil_width {{0 1 2}separate output}
562 
563    testset:
564       suffix: none_none_2d_par
565       nsize: 6
566       args: -dm_view -dim 2 -setsizes -stag_boundary_type_x none -stag_boundary_type_y none -stag_stencil_width {{0 1 2}separate output}
567 
568    testset:
569       suffix: none_ghosted_2d_par
570       nsize: 6
571       args: -dm_view -dim 2 -setsizes -stag_boundary_type_x none -stag_boundary_type_y ghosted -stag_stencil_width {{0 1 2}separate output}
572 
573    testset:
574       suffix: none_periodic_2d_par
575       nsize: 6
576       args: -dm_view -dim 2 -setsizes -stag_boundary_type_x none -stag_boundary_type_y periodic -stag_stencil_width {{0 1 2}separate output}
577 
578    testset:
579       suffix: periodic_periodic_periodic_3d_seq
580       nsize: 1
581       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}
582 
583    testset:
584       suffix: periodic_ghosted_periodic_3d_seq
585       nsize: 1
586       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}
587 
588    testset:
589       suffix: none_periodic_ghosted_3d_seq
590       nsize: 1
591       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}
592 
593    testset:
594       suffix: none_none_none_3d_seq
595       nsize: 1
596       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}
597 
598    testset:
599       suffix: periodic_periodic_periodic_3d_par
600       nsize: 12
601       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}
602 
603    testset:
604       suffix: periodic_ghosted_ghosted_3d_par
605       nsize: 12
606       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}
607 
608    testset:
609       suffix: ghosted_periodic_periodic_3d_par
610       nsize: 12
611       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}
612 
613    testset:
614       suffix: none_none_none_3d_par
615       nsize: 12
616       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}
617 
618    test:
619       suffix: periodic_none_none_3d_skinny_seq
620       nsize: 1
621       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
622 
623    test:
624       suffix: periodic_none_none_3d_skinny_par
625       nsize: 4
626       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
627 
628 TEST*/
629