1 static char help[] = "Test DMStagVecSplitToDMDA()\n\n";
2
3 #include <petscdm.h>
4 #include <petscdmstag.h>
5
main(int argc,char ** argv)6 int main(int argc, char **argv)
7 {
8 DM dm;
9 Vec x;
10 PetscBool coords;
11 PetscInt dim, dof[4];
12 PetscInt n_loc[4];
13 DMStagStencilLocation loc[4][3];
14
15 PetscFunctionBeginUser;
16 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
17 dim = 2;
18 PetscCall(PetscOptionsGetInt(NULL, NULL, "-dim", &dim, NULL));
19 coords = PETSC_TRUE;
20 PetscCall(PetscOptionsGetBool(NULL, NULL, "-coords", &coords, NULL));
21
22 // Create DMStag and setup set of locations to test
23 switch (dim) {
24 case 1:
25 PetscCall(DMStagCreate1d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, 3, 1, 2, DMSTAG_STENCIL_BOX, 1, NULL, &dm));
26 n_loc[0] = 1;
27 loc[0][0] = DMSTAG_LEFT;
28 n_loc[1] = 1;
29 loc[1][0] = DMSTAG_ELEMENT;
30 n_loc[2] = 0;
31 n_loc[3] = 0;
32 break;
33 case 2:
34 PetscCall(DMStagCreate2d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, 3, 2, PETSC_DECIDE, PETSC_DECIDE, 2, 1, 3, DMSTAG_STENCIL_BOX, 1, NULL, NULL, &dm));
35 n_loc[0] = 1;
36 loc[0][0] = DMSTAG_DOWN_LEFT;
37 n_loc[1] = 2;
38 loc[1][0] = DMSTAG_LEFT;
39 loc[1][1] = DMSTAG_DOWN;
40 n_loc[2] = 1;
41 loc[2][0] = DMSTAG_ELEMENT;
42 n_loc[3] = 0;
43 break;
44 case 3:
45 PetscCall(DMStagCreate3d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, 3, 2, 2, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, 2, 3, 1, 2, DMSTAG_STENCIL_BOX, 1, NULL, NULL, NULL, &dm));
46 n_loc[0] = 1;
47 loc[0][0] = DMSTAG_BACK_DOWN_LEFT;
48 n_loc[1] = 3;
49 loc[1][0] = DMSTAG_DOWN_LEFT;
50 loc[1][1] = DMSTAG_BACK_LEFT;
51 loc[1][2] = DMSTAG_BACK_DOWN;
52 n_loc[2] = 3;
53 loc[2][0] = DMSTAG_LEFT;
54 loc[2][1] = DMSTAG_DOWN;
55 loc[2][2] = DMSTAG_BACK;
56 n_loc[3] = 1;
57 loc[3][0] = DMSTAG_ELEMENT;
58 break;
59 default:
60 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_SUP, "No support for dimension %" PetscInt_FMT, dim);
61 }
62
63 PetscCall(DMSetFromOptions(dm));
64 PetscCall(DMSetUp(dm));
65 if (coords) PetscCall(DMStagSetUniformCoordinatesProduct(dm, -1.0, 1.0, -2.0, 2.0, -3.0, 3.0));
66
67 PetscCall(DMCreateGlobalVector(dm, &x));
68 PetscCall(VecSet(x, 1.2345));
69
70 PetscCall(DMStagGetDOF(dm, &dof[0], &dof[1], &dof[2], &dof[3]));
71 for (PetscInt stratum = 0; stratum < dim + 1; ++stratum) {
72 for (PetscInt i_loc = 0; i_loc < n_loc[stratum]; ++i_loc) {
73 // Extract 3 components, padding or truncating
74 {
75 DM da;
76 Vec x_da;
77
78 PetscCall(DMStagVecSplitToDMDA(dm, x, loc[stratum][i_loc], -3, &da, &x_da));
79 PetscCall(DMView(da, PETSC_VIEWER_STDOUT_WORLD));
80 PetscCall(VecView(x_da, PETSC_VIEWER_STDOUT_WORLD));
81 PetscCall(DMDestroy(&da));
82 PetscCall(VecDestroy(&x_da));
83 }
84
85 // Extract individual components
86 for (PetscInt c = 0; c < dof[stratum]; ++c) {
87 DM da;
88 Vec x_da;
89
90 PetscCall(DMStagVecSplitToDMDA(dm, x, loc[stratum][i_loc], c, &da, &x_da));
91 PetscCall(DMView(da, PETSC_VIEWER_STDOUT_WORLD));
92 PetscCall(VecView(x_da, PETSC_VIEWER_STDOUT_WORLD));
93 PetscCall(DMDestroy(&da));
94 PetscCall(VecDestroy(&x_da));
95 }
96 }
97 }
98
99 PetscCall(VecDestroy(&x));
100 PetscCall(DMDestroy(&dm));
101 PetscCall(PetscFinalize());
102 return 0;
103 }
104
105 /*TEST
106
107 test:
108 requires: !complex
109 args: -dim {{1 2 3}separate output} -coords {{true false}separate output}
110
111 TEST*/
112