1 static const char help[] = "Test DMCreateInjection() for mapping coordinates in 3D";
2
3 #include <petscvec.h>
4 #include <petscmat.h>
5 #include <petscdm.h>
6 #include <petscdmda.h>
7
test1_DAInjection3d(PetscInt mx,PetscInt my,PetscInt mz)8 PetscErrorCode test1_DAInjection3d(PetscInt mx, PetscInt my, PetscInt mz)
9 {
10 DM dac, daf;
11 PetscViewer vv;
12 Vec ac, af;
13 PetscInt periodicity;
14 DMBoundaryType bx, by, bz;
15
16 PetscFunctionBeginUser;
17 bx = DM_BOUNDARY_NONE;
18 by = DM_BOUNDARY_NONE;
19 bz = DM_BOUNDARY_NONE;
20
21 periodicity = 0;
22
23 PetscCall(PetscOptionsGetInt(NULL, NULL, "-periodic", &periodicity, NULL));
24 if (periodicity == 1) {
25 bx = DM_BOUNDARY_PERIODIC;
26 } else if (periodicity == 2) {
27 by = DM_BOUNDARY_PERIODIC;
28 } else if (periodicity == 3) {
29 bz = DM_BOUNDARY_PERIODIC;
30 }
31
32 PetscCall(DMDACreate3d(PETSC_COMM_WORLD, bx, by, bz, DMDA_STENCIL_BOX, mx + 1, my + 1, mz + 1, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, 1, /* 1 dof */
33 1, /* stencil = 1 */ NULL, NULL, NULL, &daf));
34 PetscCall(DMSetFromOptions(daf));
35 PetscCall(DMSetUp(daf));
36
37 PetscCall(DMCoarsen(daf, MPI_COMM_NULL, &dac));
38
39 PetscCall(DMDASetUniformCoordinates(dac, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0));
40 PetscCall(DMDASetUniformCoordinates(daf, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0));
41
42 {
43 DM cdaf, cdac;
44 Vec coordsc, coordsf, coordsf2;
45 Mat inject;
46 VecScatter vscat;
47 Mat interp;
48 PetscReal norm;
49
50 PetscCall(DMGetCoordinateDM(dac, &cdac));
51 PetscCall(DMGetCoordinateDM(daf, &cdaf));
52
53 PetscCall(DMGetCoordinates(dac, &coordsc));
54 PetscCall(DMGetCoordinates(daf, &coordsf));
55
56 PetscCall(DMCreateInjection(cdac, cdaf, &inject));
57 PetscCall(MatScatterGetVecScatter(inject, &vscat));
58 PetscCall(VecScatterBegin(vscat, coordsf, coordsc, INSERT_VALUES, SCATTER_FORWARD));
59 PetscCall(VecScatterEnd(vscat, coordsf, coordsc, INSERT_VALUES, SCATTER_FORWARD));
60 PetscCall(MatDestroy(&inject));
61
62 PetscCall(DMCreateInterpolation(cdac, cdaf, &interp, NULL));
63 PetscCall(VecDuplicate(coordsf, &coordsf2));
64 PetscCall(MatInterpolate(interp, coordsc, coordsf2));
65 PetscCall(VecAXPY(coordsf2, -1.0, coordsf));
66 PetscCall(VecNorm(coordsf2, NORM_MAX, &norm));
67 /* The fine coordinates are only reproduced in certain cases */
68 if (!bx && !by && !bz && norm > PETSC_SQRT_MACHINE_EPSILON) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Norm %g\n", (double)norm));
69 PetscCall(VecDestroy(&coordsf2));
70 PetscCall(MatDestroy(&interp));
71 }
72
73 if (0) {
74 PetscCall(DMCreateGlobalVector(dac, &ac));
75 PetscCall(VecZeroEntries(ac));
76
77 PetscCall(DMCreateGlobalVector(daf, &af));
78 PetscCall(VecZeroEntries(af));
79
80 PetscCall(PetscViewerASCIIOpen(PETSC_COMM_WORLD, "dac_7.vtu", &vv));
81 PetscCall(VecView(ac, vv));
82 PetscCall(PetscViewerDestroy(&vv));
83
84 PetscCall(PetscViewerASCIIOpen(PETSC_COMM_WORLD, "daf_7.vtu", &vv));
85 PetscCall(VecView(af, vv));
86 PetscCall(PetscViewerDestroy(&vv));
87 PetscCall(VecDestroy(&ac));
88 PetscCall(VecDestroy(&af));
89 }
90 PetscCall(DMDestroy(&dac));
91 PetscCall(DMDestroy(&daf));
92 PetscFunctionReturn(PETSC_SUCCESS);
93 }
94
main(int argc,char ** argv)95 int main(int argc, char **argv)
96 {
97 PetscInt mx, my, mz;
98
99 PetscFunctionBeginUser;
100 PetscCall(PetscInitialize(&argc, &argv, 0, help));
101 mx = 2;
102 my = 2;
103 mz = 2;
104 PetscCall(PetscOptionsGetInt(NULL, NULL, "-mx", &mx, 0));
105 PetscCall(PetscOptionsGetInt(NULL, NULL, "-my", &my, 0));
106 PetscCall(PetscOptionsGetInt(NULL, NULL, "-mz", &mz, 0));
107 PetscCall(test1_DAInjection3d(mx, my, mz));
108 PetscCall(PetscFinalize());
109 return 0;
110 }
111
112 /*TEST
113
114 test:
115 nsize: 5
116 args: -mx 30 -my 30 -mz 30 -periodic 0 -da_processors_x 5
117 output_file: output/empty.out
118
119 test:
120 suffix: 2
121 nsize: 5
122 args: -mx 29 -my 30 -mz 30 -periodic 1 -da_processors_x 5
123 output_file: output/empty.out
124
125 test:
126 suffix: 3
127 nsize: 5
128 args: -mx 30 -my 29 -mz 30 -periodic 2 -da_processors_x 5
129 output_file: output/empty.out
130
131 test:
132 suffix: 4
133 nsize: 5
134 args: -mx 30 -my 30 -mz 29 -periodic 3 -da_processors_x 5
135 output_file: output/empty.out
136
137 test:
138 suffix: 5
139 nsize: 5
140 args: -mx 30 -my 30 -mz 30 -periodic 0 -da_processors_y 5
141 output_file: output/empty.out
142
143 test:
144 suffix: 6
145 nsize: 5
146 args: -mx 29 -my 30 -mz 30 -periodic 1 -da_processors_y 5
147 output_file: output/empty.out
148
149 test:
150 suffix: 7
151 nsize: 5
152 args: -mx 30 -my 29 -mz 30 -periodic 2 -da_processors_y 5
153 output_file: output/empty.out
154
155 test:
156 suffix: 8
157 nsize: 5
158 args: -mx 30 -my 30 -mz 29 -periodic 3 -da_processors_y 5
159 output_file: output/empty.out
160
161 test:
162 suffix: 9
163 nsize: 5
164 args: -mx 30 -my 30 -mz 30 -periodic 0 -da_processors_z 5
165 output_file: output/empty.out
166
167 test:
168 suffix: 10
169 nsize: 5
170 args: -mx 29 -my 30 -mz 30 -periodic 1 -da_processors_z 5
171 output_file: output/empty.out
172
173 test:
174 suffix: 11
175 nsize: 5
176 args: -mx 30 -my 29 -mz 30 -periodic 2 -da_processors_z 5
177 output_file: output/empty.out
178
179 test:
180 suffix: 12
181 nsize: 5
182 args: -mx 30 -my 30 -mz 29 -periodic 3 -da_processors_z 5
183 output_file: output/empty.out
184
185 test:
186 suffix: 13
187 nsize: 5
188 args: -mx 30 -my 30 -mz 30 -periodic 0
189 output_file: output/empty.out
190
191 test:
192 suffix: 14
193 nsize: 5
194 args: -mx 29 -my 30 -mz 30 -periodic 1
195 output_file: output/empty.out
196
197 test:
198 suffix: 15
199 nsize: 5
200 args: -mx 30 -my 29 -mz 30 -periodic 2
201 output_file: output/empty.out
202
203 test:
204 suffix: 16
205 nsize: 5
206 args: -mx 30 -my 30 -mz 29 -periodic 3
207 output_file: output/empty.out
208
209 TEST*/
210