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 */
19: PetscInitialize(&argc, &argv, (char *)0, help);
20: dim = 2;
21: PetscOptionsGetInt(NULL, NULL, "-dim", &dim, NULL);
22: setSizes = PETSC_FALSE;
23: PetscOptionsGetBool(NULL, NULL, "-setsizes", &setSizes, NULL);
24: useInjective = PETSC_TRUE;
25: PetscOptionsGetBool(NULL, NULL, "-useinjective", &useInjective, NULL);
27: /* Creation (normal) */
28: if (!setSizes) {
29: switch (dim) {
30: case 1:
31: DMStagCreate1d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, 3, 1, 1, DMSTAG_STENCIL_BOX, 1, NULL, &dm);
32: break;
33: case 2:
34: 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: 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: MPI_Comm_size(PETSC_COMM_WORLD, &size);
50: switch (dim) {
51: case 1:
53: DMStagCreate1d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, mx, 1, 1, DMSTAG_STENCIL_BOX, 1, lx, &dm);
54: break;
55: case 2:
57: 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:
61: 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: DMSetFromOptions(dm);
70: DMSetUp(dm);
72: /* Populate Additional Injective Local-to-Global Map */
73: if (useInjective) DMStagPopulateLocalToGlobalInjective(dm);
75: /* Test: Make sure L2G inverts G2L */
76: Test1(dm);
78: /* Test: Make sure that G2L inverts L2G, on its domain */
79: DMGetDimension(dm, &dim);
80: switch (dim) {
81: case 1:
82: Test2_1d(dm);
83: break;
84: case 2:
85: Test2_2d(dm);
86: break;
87: case 3:
88: 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: DMDestroy(&dm);
96: PetscFinalize();
97: return 0;
98: }
100: static PetscErrorCode Test1(DM dm)
101: {
102: Vec vecLocal, vecGlobal, vecGlobalCheck;
103: PetscRandom rctx;
104: PetscBool equal;
107: DMCreateLocalVector(dm, &vecLocal);
108: DMCreateGlobalVector(dm, &vecGlobal);
109: PetscRandomCreate(PETSC_COMM_WORLD, &rctx);
110: VecSetRandom(vecGlobal, rctx);
111: VecSetRandom(vecLocal, rctx); /* garbage */
112: PetscRandomDestroy(&rctx);
113: VecDuplicate(vecGlobal, &vecGlobalCheck);
114: DMGlobalToLocal(dm, vecGlobal, INSERT_VALUES, vecLocal);
115: DMLocalToGlobal(dm, vecLocal, INSERT_VALUES, vecGlobalCheck);
116: VecEqual(vecGlobal, vecGlobalCheck, &equal);
118: VecDestroy(&vecLocal);
119: VecDestroy(&vecGlobal);
120: VecDestroy(&vecGlobalCheck);
121: return 0;
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: {
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: return 0;
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;
145: DMCreateLocalVector(dm, &vecLocal);
146: VecSet(vecLocal, -1.0);
147: DMStagGetCorners(dm, &startx, NULL, NULL, &nx, NULL, NULL, &nExtrax, NULL, NULL);
148: DMStagGetDOF(dm, &dof0, &dof1, NULL, NULL);
149: DMStagVecGetArray(dm, vecLocal, &arr);
150: if (dof0 > 0) DMStagGetLocationSlot(dm, DMSTAG_LEFT, 0, &idxLeft);
151: if (dof1 > 0) 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: DMStagVecRestoreArray(dm, vecLocal, &arr);
165: DMCreateGlobalVector(dm, &vecGlobal);
166: DMLocalToGlobal(dm, vecLocal, INSERT_VALUES, vecGlobal);
167: VecDuplicate(vecLocal, &vecLocalCheck);
168: VecSet(vecLocalCheck, -1.0);
169: DMGlobalToLocal(dm, vecGlobal, INSERT_VALUES, vecLocalCheck);
170: 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: 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: 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: CompareValues(i, j, k, c, val, valRef);
188: }
189: }
190: }
191: DMStagVecRestoreArrayRead(dm, vecLocalCheck, &arr);
192: VecDestroy(&vecLocal);
193: VecDestroy(&vecLocalCheck);
194: VecDestroy(&vecGlobal);
195: return 0;
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;
206: DMCreateLocalVector(dm, &vecLocal);
207: VecSet(vecLocal, -1.0);
208: DMStagGetCorners(dm, &startx, &starty, NULL, &nx, &ny, NULL, &nExtrax, &nExtray, NULL);
209: DMStagGetDOF(dm, &dof0, &dof1, &dof2, NULL);
210: DMStagVecGetArray(dm, vecLocal, &arr);
211: if (dof0 > 0) DMStagGetLocationSlot(dm, DMSTAG_DOWN_LEFT, 0, &idxDownLeft);
212: if (dof1 > 0) {
213: DMStagGetLocationSlot(dm, DMSTAG_LEFT, 0, &idxLeft);
214: DMStagGetLocationSlot(dm, DMSTAG_DOWN, 0, &idxDown);
215: }
216: if (dof2 > 0) 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: DMStagVecRestoreArray(dm, vecLocal, &arr);
244: DMCreateGlobalVector(dm, &vecGlobal);
245: DMLocalToGlobal(dm, vecLocal, INSERT_VALUES, vecGlobal);
246: VecDuplicate(vecLocal, &vecLocalCheck);
247: VecSet(vecLocalCheck, -1.0);
248: DMGlobalToLocal(dm, vecGlobal, INSERT_VALUES, vecLocalCheck);
249: 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: 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: 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: 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: 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: 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: 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: CompareValues(i, j, k, c, val, valRef);
294: }
295: }
296: }
297: }
298: DMStagVecRestoreArrayRead(dm, vecLocalCheck, &arr);
299: VecDestroy(&vecLocal);
300: VecDestroy(&vecLocalCheck);
301: VecDestroy(&vecGlobal);
302: return 0;
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;
312: DMCreateLocalVector(dm, &vecLocal);
313: VecSet(vecLocal, -1.0);
314: DMStagGetCorners(dm, &startx, &starty, &startz, &nx, &ny, &nz, &nExtrax, &nExtray, &nExtraz);
315: DMStagGetDOF(dm, &dof0, &dof1, &dof2, &dof3);
316: DMStagVecGetArray(dm, vecLocal, &arr);
317: if (dof0 > 0) DMStagGetLocationSlot(dm, DMSTAG_BACK_DOWN_LEFT, 0, &idxBackDownLeft);
318: if (dof1 > 0) {
319: DMStagGetLocationSlot(dm, DMSTAG_BACK_LEFT, 0, &idxBackLeft);
320: DMStagGetLocationSlot(dm, DMSTAG_BACK_DOWN, 0, &idxBackDown);
321: DMStagGetLocationSlot(dm, DMSTAG_DOWN_LEFT, 0, &idxDownLeft);
322: }
323: if (dof2 > 0) {
324: DMStagGetLocationSlot(dm, DMSTAG_LEFT, 0, &idxLeft);
325: DMStagGetLocationSlot(dm, DMSTAG_DOWN, 0, &idxDown);
326: DMStagGetLocationSlot(dm, DMSTAG_BACK, 0, &idxBack);
327: }
328: if (dof3 > 0) 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: DMStagVecRestoreArray(dm, vecLocal, &arr);
382: DMCreateGlobalVector(dm, &vecGlobal);
383: DMLocalToGlobal(dm, vecLocal, INSERT_VALUES, vecGlobal);
384: VecDuplicate(vecLocal, &vecLocalCheck);
385: VecSet(vecLocalCheck, -1.0);
386: DMGlobalToLocal(dm, vecGlobal, INSERT_VALUES, vecLocalCheck);
387: 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: 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: 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: 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: 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: 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: 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: 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: 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: 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: 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: 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: 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: 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: 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: CompareValues(i, j, k, c, val, valRef);
485: }
486: }
487: }
488: }
489: }
490: DMStagVecRestoreArrayRead(dm, vecLocalCheck, &arr);
491: VecDestroy(&vecLocal);
492: VecDestroy(&vecLocalCheck);
493: VecDestroy(&vecGlobal);
494: return 0;
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*/