1 static char help[] = "Spot test DMStag->DMDA routines in 3d\n\n";
2 #include <petscdm.h>
3 #include <petscdmstag.h>
4
main(int argc,char ** argv)5 int main(int argc, char **argv)
6 {
7 DM dm;
8 Vec vec;
9
10 PetscFunctionBeginUser;
11 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
12 PetscCall(DMStagCreate3d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, 4, 4, 4, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, 3, 3, 3, 3, DMSTAG_STENCIL_STAR, 1, NULL, NULL, NULL, &dm));
13 PetscCall(DMSetFromOptions(dm));
14 PetscCall(DMSetUp(dm));
15 PetscCall(DMStagSetUniformCoordinatesProduct(dm, 0.0, 10.0, 0.0, 10.0, 0.0, 10.0));
16
17 PetscCall(DMCreateGlobalVector(dm, &vec));
18 PetscCall(VecSet(vec, 1.234));
19
20 /* All element values */
21 {
22 DM da;
23 Vec vecda;
24 PetscCall(DMStagVecSplitToDMDA(dm, vec, DMSTAG_ELEMENT, -3, &da, &vecda));
25 PetscCall(DMDestroy(&da));
26 PetscCall(VecDestroy(&vecda));
27 }
28
29 /* Pad element values */
30 {
31 DM da;
32 Vec vecda;
33 PetscCall(DMStagVecSplitToDMDA(dm, vec, DMSTAG_ELEMENT, -5, &da, &vecda));
34 PetscCall(DMDestroy(&da));
35 PetscCall(VecDestroy(&vecda));
36 }
37
38 /* 2 element values */
39 {
40 DM da;
41 Vec vecda;
42 PetscCall(DMStagVecSplitToDMDA(dm, vec, DMSTAG_ELEMENT, -2, &da, &vecda));
43 PetscCall(DMDestroy(&da));
44 PetscCall(VecDestroy(&vecda));
45 }
46
47 /* One corner value */
48 {
49 DM da;
50 Vec vecda;
51 PetscCall(DMStagVecSplitToDMDA(dm, vec, DMSTAG_FRONT_DOWN_LEFT, 2, &da, &vecda));
52 PetscCall(DMDestroy(&da));
53 PetscCall(VecDestroy(&vecda));
54 }
55
56 /* One edge value */
57 {
58 DM da;
59 Vec vecda;
60 PetscCall(DMStagVecSplitToDMDA(dm, vec, DMSTAG_BACK_RIGHT, 1, &da, &vecda));
61 PetscCall(DMDestroy(&da));
62 PetscCall(VecDestroy(&vecda));
63 }
64
65 /* One face value */
66 {
67 DM da;
68 Vec vecda;
69 PetscCall(DMStagVecSplitToDMDA(dm, vec, DMSTAG_DOWN, 0, &da, &vecda));
70 PetscCall(DMDestroy(&da));
71 PetscCall(VecDestroy(&vecda));
72 }
73
74 PetscCall(VecDestroy(&vec));
75 PetscCall(DMDestroy(&dm));
76 PetscCall(PetscFinalize());
77 return 0;
78 }
79
80 /*TEST
81
82 test:
83 suffix: 1
84 nsize: 12
85 args: -stag_ranks_x 2 -stag_ranks_y 3 -stag_ranks_z 2
86 output_file: output/empty.out
87
88 TEST*/
89