Actual source code: ex33.c
1: static char help[] = "Tests for high order geometry\n\n";
3: #include <petscdmplex.h>
4: #include <petscds.h>
6: typedef struct {
7: PetscReal volume; /* Analytical volume of the mesh */
8: PetscReal tol; /* Tolerance for volume check */
9: } AppCtx;
11: PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
12: {
13: PetscFunctionBegin;
14: options->volume = -1.0;
15: options->tol = PETSC_SMALL;
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: }
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: }
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: }
39: static PetscErrorCode CreateDiscretization(DM dm, AppCtx *ctx)
40: {
41: PetscDS ds;
42: PetscFE fe;
43: DMPolytopeType ct;
44: PetscInt dim, cStart;
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: }
60: static PetscErrorCode CheckVolume(DM dm, AppCtx *ctx)
61: {
62: Vec u;
63: PetscScalar result;
64: PetscReal vol, tol = ctx->tol;
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: if (ctx->volume > 0.0 && PetscAbsReal(ctx->volume - vol) > tol) {
73: SETERRQ(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);
74: }
75: PetscFunctionReturn(PETSC_SUCCESS);
76: }
78: int main(int argc, char **argv)
79: {
80: DM dm;
81: AppCtx user;
83: PetscFunctionBeginUser;
84: PetscCall(PetscInitialize(&argc, &argv, NULL, help));
85: PetscCall(ProcessOptions(PETSC_COMM_WORLD, &user));
86: PetscCall(CreateMesh(PETSC_COMM_WORLD, &user, &dm));
87: PetscCall(CreateDiscretization(dm, &user));
88: PetscCall(CheckVolume(dm, &user));
89: PetscCall(DMDestroy(&dm));
90: PetscCall(PetscFinalize());
91: return 0;
92: }
94: /*TEST
96: testset:
97: 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.
99: test:
100: suffix: square_0
101: args: -dm_coord_petscspace_degree 1
103: test:
104: suffix: square_1
105: args: -dm_coord_petscspace_degree 2
107: test:
108: suffix: square_2
109: args: -dm_refine 1 -dm_coord_petscspace_degree 1
111: test:
112: suffix: square_3
113: args: -dm_refine 1 -dm_coord_petscspace_degree 2
115: testset:
116: 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.
118: test:
119: suffix: cube_0
120: args: -dm_coord_petscspace_degree 1
122: test:
123: suffix: cube_1
124: args: -dm_coord_petscspace_degree 2
126: test:
127: suffix: cube_2
128: args: -dm_refine 1 -dm_coord_petscspace_degree 1
130: test:
131: suffix: cube_3
132: args: -dm_refine 1 -dm_coord_petscspace_degree 2
134: testset:
135: args: -dm_plex_simplex 0 -dm_plex_box_faces 1,1 -dm_plex_box_lower -1.,-1. -dm_plex_box_upper 1.,1. \
136: -volume 4. -dm_coord_remap -dm_coord_map shear -dm_coord_map_params 0.0,0.0,3.0
138: test:
139: suffix: shear_0
140: args: -dm_coord_petscspace_degree 1
142: test:
143: suffix: shear_1
144: args: -dm_coord_petscspace_degree 2
146: test:
147: suffix: shear_2
148: args: -dm_refine 1 -dm_coord_petscspace_degree 1
150: test:
151: suffix: shear_3
152: args: -dm_refine 1 -dm_coord_petscspace_degree 2
154: testset:
155: 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. \
156: -volume 8. -dm_coord_remap -dm_coord_map shear -dm_coord_map_params 0.0,0.0,3.0,4.0
158: test:
159: suffix: shear_4
160: args: -dm_coord_petscspace_degree 1
162: test:
163: suffix: shear_5
164: args: -dm_coord_petscspace_degree 2
166: test:
167: suffix: shear_6
168: args: -dm_refine 1 -dm_coord_petscspace_degree 1
170: test:
171: suffix: shear_7
172: args: -dm_refine 1 -dm_coord_petscspace_degree 2
174: testset:
175: args: -dm_plex_simplex 0 -dm_plex_box_faces 1,1 -dm_plex_box_lower 1.,-1. -dm_plex_box_upper 3.,1. \
176: -dm_coord_petscspace_degree 1 -volume 8. -dm_coord_remap -dm_coord_map flare
178: test:
179: suffix: flare_0
180: args:
182: test:
183: suffix: flare_1
184: args: -dm_refine 2
186: testset:
187: # Area: 3/4 \pi = 2.3562
188: args: -dm_plex_simplex 0 -dm_plex_box_faces 1,1 -petscfe_default_quadrature_order 2 \
189: -volume 2.35619449019235 -dm_coord_remap -dm_coord_map annulus
191: test:
192: # Area: (a+b)/2 h = 3/\sqrt{2} (sqrt{2} - 1/\sqrt{2}) = 3/2
193: suffix: annulus_0
194: requires: double
195: args: -dm_coord_petscspace_degree 1 -volume 1.5
197: test:
198: suffix: annulus_1
199: requires: double
200: args: -dm_refine 3 -dm_coord_petscspace_degree 1 -tol .016
202: test:
203: suffix: annulus_2
204: requires: double
205: args: -dm_refine 3 -dm_coord_petscspace_degree 2 -tol .0038
207: test:
208: suffix: annulus_3
209: requires: double
210: args: -dm_refine 3 -dm_coord_petscspace_degree 3 -tol 2.2e-6
212: test:
213: suffix: annulus_4
214: requires: double
215: args: -dm_refine 2 -dm_coord_petscspace_degree 2 -tol .00012
217: test:
218: suffix: annulus_5
219: requires: double
220: args: -dm_refine 2 -dm_coord_petscspace_degree 3 -petscfe_default_quadrature_order 3 -tol 1.2e-7
222: testset:
223: # Volume: 4/3 \pi (8 - 1)/2 = 14/3 \pi = 14.66076571675238
224: 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. \
225: -volume 14.66076571675238 -dm_coord_remap -dm_coord_map shell
227: test:
228: suffix: shell_0
229: requires: double
230: args: -dm_refine 1 -dm_coord_petscspace_degree 1 -petscfe_default_quadrature_order 1 -volume 5.633164922 -tol 1.0e-7
232: test:
233: suffix: shell_1
234: requires: double
235: args: -dm_refine 2 -dm_coord_petscspace_degree 1 -petscfe_default_quadrature_order 1 -tol 3.1
237: test:
238: suffix: shell_2
239: requires: double
240: args: -dm_refine 2 -dm_coord_petscspace_degree 2 -petscfe_default_quadrature_order 2 -tol .1
242: test:
243: suffix: shell_3
244: requires: double
245: args: -dm_refine 2 -dm_coord_petscspace_degree 3 -petscfe_default_quadrature_order 3 -tol .02
247: test:
248: suffix: shell_4
249: requires: double
250: args: -dm_refine 2 -dm_coord_petscspace_degree 4 -petscfe_default_quadrature_order 4 -tol .006
252: test:
253: # Volume: 1.0
254: suffix: gmsh_q2
255: requires: double
256: args: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/quads-q2.msh -dm_plex_gmsh_project -volume 1.0 -tol 1e-6
258: test:
259: # Volume: 1.0
260: suffix: gmsh_q3
261: requires: double
262: nsize: {{1 2}}
263: args: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/quads-q3.msh -dm_plex_gmsh_project -volume 1.0 -tol 1e-6
265: test:
266: # Volume: 1.0
267: suffix: gmsh_3d_q2
268: requires: double
269: args: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/cube_q2.msh -dm_plex_gmsh_project -volume 1.0 -tol 1e-6
271: test:
272: # Volume: 1.0
273: suffix: gmsh_3d_q3
274: requires: double
275: nsize: {{1 2}}
276: args: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/cube_q3.msh -dm_plex_gmsh_project -volume 1.0 -tol 1e-6
278: TEST*/