1 static char help[] = "Tests for high order geometry\n\n";
2
3 #include <petscdmplex.h>
4 #include <petscds.h>
5
6 typedef struct {
7 PetscReal volume; /* Analytical volume of the mesh */
8 PetscReal tol; /* Tolerance for volume check */
9 } AppCtx;
10
ProcessOptions(MPI_Comm comm,AppCtx * options)11 PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
12 {
13 PetscFunctionBegin;
14 options->volume = -1.0;
15 options->tol = PETSC_SMALL;
16
17 PetscOptionsBegin(comm, "", "Meshing Interpolation Test Options", "DMPLEX");
18 PetscCall(PetscOptionsReal("-volume", "The analytical volume of the mesh", "ex33.c", options->volume, &options->volume, NULL));
19 PetscCall(PetscOptionsReal("-tol", "The tolerance for the volume check", "ex33.c", options->tol, &options->tol, NULL));
20 PetscOptionsEnd();
21 PetscFunctionReturn(PETSC_SUCCESS);
22 }
23
CreateMesh(MPI_Comm comm,AppCtx * ctx,DM * dm)24 PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *ctx, DM *dm)
25 {
26 PetscFunctionBegin;
27 PetscCall(DMCreate(comm, dm));
28 PetscCall(DMSetType(*dm, DMPLEX));
29 PetscCall(DMSetFromOptions(*dm));
30 PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view"));
31 PetscFunctionReturn(PETSC_SUCCESS);
32 }
33
volume(PetscInt dim,PetscInt Nf,PetscInt NfAux,const PetscInt uOff[],const PetscInt uOff_x[],const PetscScalar u[],const PetscScalar u_t[],const PetscScalar u_x[],const PetscInt aOff[],const PetscInt aOff_x[],const PetscScalar a[],const PetscScalar a_t[],const PetscScalar a_x[],PetscReal t,const PetscReal x[],PetscInt numConstants,const PetscScalar constants[],PetscScalar vol[])34 static void volume(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar vol[])
35 {
36 vol[0] = 1.;
37 }
38
CreateDiscretization(DM dm,AppCtx * ctx)39 static PetscErrorCode CreateDiscretization(DM dm, AppCtx *ctx)
40 {
41 PetscDS ds;
42 PetscFE fe;
43 DMPolytopeType ct;
44 PetscInt dim, cStart;
45
46 PetscFunctionBeginUser;
47 PetscCall(DMGetDimension(dm, &dim));
48 PetscCall(DMPlexGetHeightStratum(dm, 0, &cStart, NULL));
49 PetscCall(DMPlexGetCellType(dm, cStart, &ct));
50 PetscCall(PetscFECreateByCell(PETSC_COMM_SELF, dim, 1, ct, NULL, PETSC_DETERMINE, &fe));
51 PetscCall(PetscFESetName(fe, "scalar"));
52 PetscCall(DMAddField(dm, NULL, (PetscObject)fe));
53 PetscCall(PetscFEDestroy(&fe));
54 PetscCall(DMCreateDS(dm));
55 PetscCall(DMGetDS(dm, &ds));
56 PetscCall(PetscDSSetObjective(ds, 0, volume));
57 PetscFunctionReturn(PETSC_SUCCESS);
58 }
59
CheckVolume(DM dm,AppCtx * ctx)60 static PetscErrorCode CheckVolume(DM dm, AppCtx *ctx)
61 {
62 Vec u;
63 PetscScalar result;
64 PetscReal vol, tol = ctx->tol;
65
66 PetscFunctionBeginUser;
67 PetscCall(DMGetGlobalVector(dm, &u));
68 PetscCall(DMPlexComputeIntegralFEM(dm, u, &result, ctx));
69 vol = PetscRealPart(result);
70 PetscCall(DMRestoreGlobalVector(dm, &u));
71 PetscCall(PetscPrintf(PetscObjectComm((PetscObject)dm), "Volume: %g\n", (double)vol));
72 PetscCheck(ctx->volume <= 0.0 || PetscAbsReal(ctx->volume - vol) <= tol, PetscObjectComm((PetscObject)dm), PETSC_ERR_PLIB, "Calculated volume %g != %g actual volume (error %g > %g tol)", (double)vol, (double)ctx->volume, (double)PetscAbsReal(ctx->volume - vol), (double)tol);
73 PetscFunctionReturn(PETSC_SUCCESS);
74 }
75
main(int argc,char ** argv)76 int main(int argc, char **argv)
77 {
78 DM dm;
79 AppCtx user;
80
81 PetscFunctionBeginUser;
82 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
83 PetscCall(ProcessOptions(PETSC_COMM_WORLD, &user));
84 PetscCall(CreateMesh(PETSC_COMM_WORLD, &user, &dm));
85 PetscCall(CreateDiscretization(dm, &user));
86 PetscCall(CheckVolume(dm, &user));
87 PetscCall(DMDestroy(&dm));
88 PetscCall(PetscFinalize());
89 return 0;
90 }
91
92 /*TEST
93
94 testset:
95 args: -dm_plex_simplex 0 -dm_plex_box_faces 1,1 -dm_plex_box_lower -1.,-1. -dm_plex_box_upper 1.,1. -volume 4.
96
97 test:
98 suffix: square_0
99 args: -dm_coord_petscspace_degree 1
100
101 test:
102 suffix: square_1
103 args: -dm_coord_petscspace_degree 2
104
105 test:
106 suffix: square_2
107 args: -dm_refine 1 -dm_coord_petscspace_degree 1
108
109 test:
110 suffix: square_3
111 args: -dm_refine 1 -dm_coord_petscspace_degree 2
112
113 testset:
114 args: -dm_plex_dim 3 -dm_plex_simplex 0 -dm_plex_box_faces 1,1,1 -dm_plex_box_lower -1.,-1.,-1. -dm_plex_box_upper 1.,1.,1. -volume 8.
115
116 test:
117 suffix: cube_0
118 args: -dm_coord_petscspace_degree 1
119
120 test:
121 suffix: cube_1
122 args: -dm_coord_petscspace_degree 2
123
124 test:
125 suffix: cube_2
126 args: -dm_refine 1 -dm_coord_petscspace_degree 1
127
128 test:
129 suffix: cube_3
130 args: -dm_refine 1 -dm_coord_petscspace_degree 2
131
132 testset:
133 args: -dm_plex_simplex 0 -dm_plex_box_faces 1,1 -dm_plex_box_lower -1.,-1. -dm_plex_box_upper 1.,1. \
134 -volume 4. -dm_coord_remap -dm_coord_map shear -dm_coord_map_params 0.0,0.0,3.0
135
136 test:
137 suffix: shear_0
138 args: -dm_coord_petscspace_degree 1
139
140 test:
141 suffix: shear_1
142 args: -dm_coord_petscspace_degree 2
143
144 test:
145 suffix: shear_2
146 args: -dm_refine 1 -dm_coord_petscspace_degree 1
147
148 test:
149 suffix: shear_3
150 args: -dm_refine 1 -dm_coord_petscspace_degree 2
151
152 testset:
153 args: -dm_plex_dim 3 -dm_plex_simplex 0 -dm_plex_box_faces 1,1,1 -dm_plex_box_lower -1.,-1.,-1. -dm_plex_box_upper 1.,1.,1. \
154 -volume 8. -dm_coord_remap -dm_coord_map shear -dm_coord_map_params 0.0,0.0,3.0,4.0
155
156 test:
157 suffix: shear_4
158 args: -dm_coord_petscspace_degree 1
159
160 test:
161 suffix: shear_5
162 args: -dm_coord_petscspace_degree 2
163
164 test:
165 suffix: shear_6
166 args: -dm_refine 1 -dm_coord_petscspace_degree 1
167
168 test:
169 suffix: shear_7
170 args: -dm_refine 1 -dm_coord_petscspace_degree 2
171
172 testset:
173 args: -dm_plex_simplex 0 -dm_plex_box_faces 1,1 -dm_plex_box_lower 1.,-1. -dm_plex_box_upper 3.,1. \
174 -dm_coord_petscspace_degree 1 -volume 8. -dm_coord_remap -dm_coord_map flare
175
176 test:
177 suffix: flare_0
178 args:
179
180 test:
181 suffix: flare_1
182 args: -dm_refine 2
183
184 testset:
185 # Area: 3/4 \pi = 2.3562
186 args: -dm_plex_simplex 0 -dm_plex_box_faces 1,1 -petscfe_default_quadrature_order 2 \
187 -volume 2.35619449019235 -dm_coord_remap -dm_coord_map annulus
188
189 test:
190 # Area: (a+b)/2 h = 3/\sqrt{2} (sqrt{2} - 1/\sqrt{2}) = 3/2
191 suffix: annulus_0
192 requires: double
193 args: -dm_coord_petscspace_degree 1 -volume 1.5
194
195 test:
196 suffix: annulus_1
197 requires: double
198 args: -dm_refine 3 -dm_coord_petscspace_degree 1 -tol .016
199
200 test:
201 suffix: annulus_2
202 requires: double
203 args: -dm_refine 3 -dm_coord_petscspace_degree 2 -tol .0038
204
205 test:
206 suffix: annulus_3
207 requires: double
208 args: -dm_refine 3 -dm_coord_petscspace_degree 3 -tol 2.2e-6
209
210 test:
211 suffix: annulus_4
212 requires: double
213 args: -dm_refine 2 -dm_coord_petscspace_degree 2 -tol .00012
214
215 test:
216 suffix: annulus_5
217 requires: double
218 args: -dm_refine 2 -dm_coord_petscspace_degree 3 -petscfe_default_quadrature_order 3 -tol 1.2e-7
219
220 testset:
221 # Volume: 4/3 \pi (8 - 1)/2 = 14/3 \pi = 14.66076571675238
222 args: -dm_plex_dim 3 -dm_plex_simplex 0 -dm_plex_box_faces 1,1,1 -dm_plex_box_lower -1.,-1.,-1. -dm_plex_box_upper 1.,1.,1. \
223 -volume 14.66076571675238 -dm_coord_remap -dm_coord_map shell
224
225 test:
226 suffix: shell_0
227 requires: double
228 args: -dm_refine 1 -dm_coord_petscspace_degree 1 -petscfe_default_quadrature_order 1 -volume 5.633164922 -tol 1.0e-7
229
230 test:
231 suffix: shell_1
232 requires: double
233 args: -dm_refine 2 -dm_coord_petscspace_degree 1 -petscfe_default_quadrature_order 1 -tol 3.1
234
235 test:
236 suffix: shell_2
237 requires: double
238 args: -dm_refine 2 -dm_coord_petscspace_degree 2 -petscfe_default_quadrature_order 2 -tol .1
239
240 test:
241 suffix: shell_3
242 requires: double
243 args: -dm_refine 2 -dm_coord_petscspace_degree 3 -petscfe_default_quadrature_order 3 -tol .02
244
245 test:
246 suffix: shell_4
247 requires: double
248 args: -dm_refine 2 -dm_coord_petscspace_degree 4 -petscfe_default_quadrature_order 4 -tol .006
249
250 test:
251 # Volume: 1.0
252 suffix: gmsh_q2
253 requires: double
254 args: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/quads-q2.msh -dm_plex_gmsh_project -volume 1.0 -tol 1e-6
255
256 test:
257 # Volume: 1.0
258 suffix: gmsh_q3
259 requires: double
260 nsize: {{1 2}}
261 args: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/quads-q3.msh -dm_plex_gmsh_project -volume 1.0 -tol 1e-6
262
263 test:
264 # Volume: 1.0
265 suffix: gmsh_3d_q2
266 requires: double
267 args: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/cube_q2.msh -dm_plex_gmsh_project -volume 1.0 -tol 1e-6
268
269 test:
270 # Volume: 1.0
271 suffix: gmsh_3d_q3
272 requires: double
273 nsize: {{1 2}}
274 args: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/cube_q3.msh -dm_plex_gmsh_project -volume 1.0 -tol 1e-6
275
276 TEST*/
277