Actual source code: ex19.c
1: static char help[] = "(Partially) test DMStag default interpolation, 2d faces-only.\n\n";
3: #include <petscdm.h>
4: #include <petscdmstag.h>
5: #include <petscksp.h>
7: PetscErrorCode CreateSystem(DM dm, Mat *A, Vec *b);
9: int main(int argc, char **argv)
10: {
11: DM dm, dmCoarse;
12: Mat Ai;
15: PetscInitialize(&argc, &argv, (char *)0, help);
16: DMStagCreate2d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, 2, 4, PETSC_DECIDE, PETSC_DECIDE, 0, 1, 1, DMSTAG_STENCIL_BOX, 1, NULL, NULL, &dm);
17: DMSetFromOptions(dm);
18: DMSetUp(dm);
19: DMCoarsen(dm, MPI_COMM_NULL, &dmCoarse);
20: DMCreateInterpolation(dmCoarse, dm, &Ai, NULL);
22: /* See what happens to a constant value on each sub-grid */
23: {
24: Vec localCoarse, globalCoarse, globalFine, localFine;
25: DMGetGlobalVector(dm, &globalFine);
26: DMGetGlobalVector(dmCoarse, &globalCoarse);
27: DMGetLocalVector(dmCoarse, &localCoarse);
28: DMGetLocalVector(dm, &localFine);
29: VecSet(localCoarse, -1.0);
30: VecSet(localFine, -1.0);
31: {
32: PetscInt i, j, startx, starty, nx, ny, extrax, extray;
33: PetscInt p, vx, vy;
34: PetscScalar ***arr;
35: DMStagGetCorners(dmCoarse, &startx, &starty, NULL, &nx, &ny, NULL, &extrax, &extray, NULL);
36: DMStagVecGetArray(dmCoarse, localCoarse, &arr);
37: DMStagGetLocationSlot(dmCoarse, DMSTAG_LEFT, 0, &vx);
38: DMStagGetLocationSlot(dmCoarse, DMSTAG_DOWN, 0, &vy);
39: DMStagGetLocationSlot(dmCoarse, DMSTAG_ELEMENT, 0, &p);
40: for (j = starty; j < starty + ny + extray; ++j) {
41: for (i = startx; i < startx + nx + extrax; ++i) {
42: arr[j][i][vy] = (i < startx + nx) ? 10.0 : -1;
43: arr[j][i][vx] = (j < starty + ny) ? 20.0 : -1;
44: arr[j][i][p] = (i < startx + nx) && (j < starty + ny) ? 30.0 : -1;
45: }
46: }
47: DMStagVecRestoreArray(dmCoarse, localCoarse, &arr);
48: }
49: DMLocalToGlobal(dmCoarse, localCoarse, INSERT_VALUES, globalCoarse);
50: MatInterpolate(Ai, globalCoarse, globalFine);
51: DMGlobalToLocal(dm, globalFine, INSERT_VALUES, localFine);
52: {
53: PetscInt i, j, startx, starty, nx, ny, extrax, extray;
54: PetscInt p, vx, vy;
55: PetscScalar ***arr;
56: DMStagGetCorners(dm, &startx, &starty, NULL, &nx, &ny, NULL, &extrax, &extray, NULL);
57: DMStagVecGetArrayRead(dm, localFine, &arr);
58: DMStagGetLocationSlot(dm, DMSTAG_LEFT, 0, &vx);
59: DMStagGetLocationSlot(dm, DMSTAG_DOWN, 0, &vy);
60: DMStagGetLocationSlot(dm, DMSTAG_ELEMENT, 0, &p);
61: for (j = starty; j < starty + ny + extray; ++j) {
62: for (i = startx; i < startx + nx + extrax; ++i) {
63: const PetscScalar expected_vy = (i < startx + nx) ? 10.0 : -1;
64: const PetscScalar expected_vx = (j < starty + ny) ? 20.0 : -1;
65: const PetscScalar expected_p = (i < startx + nx) && (j < starty + ny) ? 30.0 : -1;
66: if (arr[j][i][vy] != expected_vy) PetscPrintf(PETSC_COMM_SELF, "wrong %" PetscInt_FMT " %" PetscInt_FMT "\n", i, j);
67: if (arr[j][i][vx] != expected_vx) PetscPrintf(PETSC_COMM_SELF, "wrong %" PetscInt_FMT " %" PetscInt_FMT "\n", i, j);
68: if (arr[j][i][p] != expected_p) PetscPrintf(PETSC_COMM_SELF, "wrong %" PetscInt_FMT " %" PetscInt_FMT "\n", i, j);
69: }
70: }
71: DMStagVecRestoreArrayRead(dm, localFine, &arr);
72: }
73: DMRestoreLocalVector(dmCoarse, &localCoarse);
74: DMRestoreLocalVector(dm, &localFine);
75: DMRestoreGlobalVector(dmCoarse, &globalCoarse);
76: DMRestoreGlobalVector(dm, &globalFine);
77: }
79: MatDestroy(&Ai);
80: DMDestroy(&dm);
81: DMDestroy(&dmCoarse);
82: PetscFinalize();
83: return 0;
84: }
86: /*TEST
88: test:
89: suffix: 1
90: nsize: 1
91: args:
93: test:
94: suffix: 2
95: nsize: 4
96: args: -stag_grid_x 8 -stag_grid_y 4
98: TEST*/