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: 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: }
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;
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: }
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;
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: }
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;
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
496: /*TEST
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}
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}
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}
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}
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}
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}
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}
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}
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}
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}
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}
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}
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}
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}
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}
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}
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}
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}
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}
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}
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}
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}
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}
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}
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
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
628: TEST*/