1 static char help[] = "Tests DMPlex Gmsh reader.\n\n";
2
3 #include <petscdmplex.h>
4
5 #if !defined(PETSC_GMSH_EXE)
6 #define PETSC_GMSH_EXE "gmsh"
7 #endif
8
9 #include <petscds.h>
10
one(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 value[])11 static void one(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 value[])
12 {
13 value[0] = (PetscReal)1;
14 }
15
CreateFE(DM dm)16 static PetscErrorCode CreateFE(DM dm)
17 {
18 DM cdm;
19 PetscSpace P;
20 PetscDualSpace Q;
21 DM K;
22 PetscFE fe;
23 DMPolytopeType ptype;
24
25 PetscInt dim, k;
26 PetscBool isSimplex;
27
28 PetscDS ds;
29
30 PetscFunctionBeginUser;
31 PetscCall(DMGetCoordinateDM(dm, &cdm));
32 PetscCall(DMGetField(cdm, 0, NULL, (PetscObject *)&fe));
33 PetscCall(PetscFEGetBasisSpace(fe, &P));
34 PetscCall(PetscFEGetDualSpace(fe, &Q));
35 PetscCall(PetscDualSpaceGetDM(Q, &K));
36 PetscCall(DMGetDimension(K, &dim));
37 PetscCall(PetscSpaceGetDegree(P, &k, NULL));
38 PetscCall(DMPlexGetCellType(K, 0, &ptype));
39 switch (ptype) {
40 case DM_POLYTOPE_QUADRILATERAL:
41 case DM_POLYTOPE_HEXAHEDRON:
42 isSimplex = PETSC_FALSE;
43 break;
44 default:
45 isSimplex = PETSC_TRUE;
46 break;
47 }
48
49 PetscCall(PetscFECreateLagrange(PETSC_COMM_SELF, dim, 1, isSimplex, k, PETSC_DETERMINE, &fe));
50 PetscCall(PetscFESetName(fe, "scalar"));
51 PetscCall(DMAddField(dm, NULL, (PetscObject)fe));
52 PetscCall(PetscFEDestroy(&fe));
53 PetscCall(DMCreateDS(dm));
54
55 PetscCall(DMGetDS(dm, &ds));
56 PetscCall(PetscDSSetObjective(ds, 0, one));
57 PetscFunctionReturn(PETSC_SUCCESS);
58 }
59
CheckIntegral(DM dm,PetscReal integral,PetscReal tol)60 static PetscErrorCode CheckIntegral(DM dm, PetscReal integral, PetscReal tol)
61 {
62 Vec u;
63 PetscReal rval;
64 PetscScalar result;
65
66 PetscFunctionBeginUser;
67 PetscCall(DMGetGlobalVector(dm, &u));
68 PetscCall(DMPlexComputeIntegralFEM(dm, u, &result, NULL));
69 PetscCall(DMRestoreGlobalVector(dm, &u));
70 rval = PetscRealPart(result);
71 if (integral > 0 && PetscAbsReal(integral - rval) > tol) {
72 PetscCall(PetscPrintf(PetscObjectComm((PetscObject)dm), "Calculated value %g != %g actual value (error %g > %g tol)\n", (double)rval, (double)integral, (double)PetscAbsReal(integral - rval), (double)tol));
73 }
74 PetscFunctionReturn(PETSC_SUCCESS);
75 }
76
main(int argc,char ** argv)77 int main(int argc, char **argv)
78 {
79 DM dm;
80 const char *const mshlist[] = {"seg", "tri", "qua", "tet", "wed", "hex", "vtx", "B2tri", "B2qua", "B3tet", "B3hex"};
81 const char *const fmtlist[] = {"msh22", "msh40", "msh41"};
82 PetscInt msh = 5;
83 PetscInt fmt = 2;
84 PetscBool bin = PETSC_TRUE;
85 PetscInt dim = 3;
86 PetscInt order = 2;
87
88 const char cmdtemplate[] = "%s -format %s %s -%d -order %d %s -o %s";
89 char gmsh[PETSC_MAX_PATH_LEN] = PETSC_GMSH_EXE;
90 char tag[PETSC_MAX_PATH_LEN], path[PETSC_MAX_PATH_LEN];
91 char geo[PETSC_MAX_PATH_LEN], geodir[PETSC_MAX_PATH_LEN] = ".";
92 char out[PETSC_MAX_PATH_LEN], outdir[PETSC_MAX_PATH_LEN] = ".";
93 char cmd[PETSC_MAX_PATH_LEN * 4];
94 PetscBool set, flg;
95 FILE *fp;
96
97 PetscFunctionBeginUser;
98 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
99
100 PetscCall(PetscStrncpy(geodir, "${PETSC_DIR}/share/petsc/datafiles/meshes", sizeof(geodir)));
101 PetscCall(PetscOptionsGetenv(PETSC_COMM_SELF, "GMSH", path, sizeof(path), &set));
102 if (set) PetscCall(PetscStrncpy(gmsh, path, sizeof(gmsh)));
103 PetscCall(PetscOptionsGetString(NULL, NULL, "-gmsh", gmsh, sizeof(gmsh), NULL));
104 PetscCall(PetscOptionsGetString(NULL, NULL, "-dir", geodir, sizeof(geodir), NULL));
105 PetscCall(PetscOptionsGetString(NULL, NULL, "-out", outdir, sizeof(outdir), NULL));
106 PetscCall(PetscOptionsGetEList(NULL, NULL, "-msh", mshlist, PETSC_STATIC_ARRAY_LENGTH(mshlist), &msh, NULL));
107 PetscCall(PetscOptionsGetEList(NULL, NULL, "-fmt", fmtlist, PETSC_STATIC_ARRAY_LENGTH(fmtlist), &fmt, NULL));
108 PetscCall(PetscOptionsGetBool(NULL, NULL, "-bin", &bin, NULL));
109 PetscCall(PetscOptionsGetInt(NULL, NULL, "-dim", &dim, NULL));
110 PetscCall(PetscOptionsGetInt(NULL, NULL, "-order", &order, NULL));
111 if (fmt == 1) bin = PETSC_FALSE; /* Recent Gmsh releases cannot generate msh40+binary format*/
112
113 { /* This test requires Gmsh >= 4.2.0 */
114 char space[PETSC_MAX_PATH_LEN];
115 int inum = 0, major = 0, minor = 0, micro = 0;
116 PetscCall(PetscSNPrintf(cmd, sizeof(cmd), "%s -info", gmsh));
117 PetscCall(PetscPOpen(PETSC_COMM_SELF, NULL, cmd, "r", &fp));
118 if (fp) inum = fscanf(fp, "Version %s %d.%d.%d", space, &major, &minor, µ);
119 PetscCall(PetscPClose(PETSC_COMM_SELF, fp));
120 if (inum != 4 || major < 4 || (major == 4 && minor < 2)) {
121 PetscCall(PetscPrintf(PETSC_COMM_SELF, "Gmsh>=4.2.0 not available\n"));
122 goto finish;
123 }
124 }
125
126 PetscCall(PetscSNPrintf(tag, sizeof(tag), "%s-%d-%d-%s%s", mshlist[msh], (int)dim, (int)order, fmtlist[fmt], bin ? "-bin" : ""));
127 PetscCall(PetscSNPrintf(geo, sizeof(geo), "%s/gmsh-%s.geo", geodir, mshlist[msh]));
128 PetscCall(PetscSNPrintf(out, sizeof(out), "%s/mesh-%s.msh", outdir, tag));
129 PetscCall(PetscStrreplace(PETSC_COMM_SELF, geo, path, sizeof(path)));
130 PetscCall(PetscFixFilename(path, geo));
131 PetscCall(PetscStrreplace(PETSC_COMM_SELF, out, path, sizeof(path)));
132 PetscCall(PetscFixFilename(path, out));
133 PetscCall(PetscTestFile(geo, 'r', &flg));
134 PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_USER_INPUT, "File not found: %s", geo);
135
136 PetscCall(PetscSNPrintf(cmd, sizeof(cmd), cmdtemplate, gmsh, fmtlist[fmt], bin ? "-bin" : "", (int)dim, (int)order, geo, out));
137 PetscCall(PetscPOpen(PETSC_COMM_SELF, NULL, cmd, "r", &fp));
138 PetscCall(PetscPClose(PETSC_COMM_SELF, fp));
139
140 PetscCall(DMPlexCreateFromFile(PETSC_COMM_SELF, out, "ex99_plex", PETSC_TRUE, &dm));
141 PetscCall(PetscSNPrintf(tag, sizeof(tag), "mesh-%s", mshlist[msh]));
142 PetscCall(PetscObjectSetName((PetscObject)dm, tag));
143 PetscCall(DMSetFromOptions(dm));
144 PetscCall(DMViewFromOptions(dm, NULL, "-dm_view"));
145 {
146 PetscBool check;
147 PetscReal integral = 0, tol = (PetscReal)1.0e-4;
148 PetscCall(PetscOptionsGetReal(NULL, NULL, "-integral", &integral, &check));
149 PetscCall(PetscOptionsGetReal(NULL, NULL, "-tol", &tol, NULL));
150 if (check) {
151 PetscCall(CreateFE(dm));
152 PetscCall(CheckIntegral(dm, integral, tol));
153 }
154 }
155 PetscCall(DMDestroy(&dm));
156
157 finish:
158 PetscCall(PetscFinalize());
159 return 0;
160 }
161
162 /*TEST
163
164 build:
165 requires: defined(PETSC_HAVE_POPEN)
166
167 test:
168 requires: defined(PETSC_GMSH_EXE)
169 args: -dir ${wPETSC_DIR}/share/petsc/datafiles/meshes
170 args: -msh {{vtx}separate_output}
171 args: -order 1
172 args: -fmt {{msh22 msh40 msh41}} -bin {{0 1}}
173 args: -dm_view ::ascii_info_detail
174 args: -dm_plex_check_all
175 args: -dm_plex_gmsh_highorder false
176 args: -dm_plex_boundary_label marker
177 args: -dm_plex_gmsh_spacedim 3
178
179 test:
180 requires: defined(PETSC_GMSH_EXE)
181 args: -dir ${wPETSC_DIR}/share/petsc/datafiles/meshes
182 args: -msh {{seg tri qua tet wed hex}separate_output}
183 args: -order {{1 2 3 7}}
184 args: -fmt {{msh22 msh40 msh41}} -bin {{0 1}}
185 args: -dm_view ::ascii_info_detail
186 args: -dm_plex_check_all
187 args: -dm_plex_gmsh_highorder false
188 args: -dm_plex_boundary_label marker
189
190 testset:
191 suffix: B2 # 2D ball
192 output_file: output/empty.out
193 requires: defined(PETSC_GMSH_EXE)
194 args: -dir ${wPETSC_DIR}/share/petsc/datafiles/meshes
195 args: -msh {{B2tri B2qua}}
196 args: -dim 2 -integral 3.141592653589793 # pi
197 args: -order {{2 3 4 5 6 7 8 9}} -tol 0.05
198
199 testset:
200 suffix: B2_bnd # 2D ball boundary
201 output_file: output/empty.out
202 requires: defined(PETSC_GMSH_EXE)
203 args: -dir ${wPETSC_DIR}/share/petsc/datafiles/meshes
204 args: -dm_plex_gmsh_spacedim 2
205 args: -msh {{B2tri B2qua}}
206 args: -dim 1 -integral 6.283185307179586 # 2*pi
207 args: -order {{2 3 4 5 6 7 8 9}} -tol 0.05
208
209 testset:
210 suffix: B3 # 3D ball
211 output_file: output/empty.out
212 requires: defined(PETSC_GMSH_EXE)
213 args: -dir ${wPETSC_DIR}/share/petsc/datafiles/meshes
214 args: -msh {{B3tet B3hex}}
215 args: -dim 3 -integral 4.1887902047863905 # 4/3*pi
216 args: -order {{2 3 4 5}} -tol 0.20
217
218 testset:
219 suffix: B3_bnd # 3D ball boundary
220 output_file: output/empty.out
221 requires: defined(PETSC_GMSH_EXE)
222 args: -dir ${wPETSC_DIR}/share/petsc/datafiles/meshes
223 args: -dm_plex_gmsh_spacedim 3
224 args: -msh {{B3tet B3hex}}
225 args: -dim 2 -integral 12.566370614359172 # 4*pi
226 args: -order {{2 3 4 5 6 7 8 9}} -tol 0.25
227
228 testset:
229 suffix: B_lin # linear discretizations
230 requires: defined(PETSC_GMSH_EXE)
231 args: -dir ${wPETSC_DIR}/share/petsc/datafiles/meshes
232 args: -dm_plex_gmsh_highorder true
233 args: -dm_plex_gmsh_project true
234 args: -dm_plex_gmsh_project_petscspace_degree {{1 2 3}separate_output}
235 args: -dm_plex_gmsh_fe_view
236 args: -dm_plex_gmsh_project_fe_view
237 args: -order 1 -tol 1e-4
238 test:
239 suffix: dim-1
240 args: -dm_plex_gmsh_spacedim 2
241 args: -msh {{B2tri B2qua}separate_output}
242 args: -dim 1 -integral 5.656854249492381 # 4*sqrt(2)
243 test:
244 suffix: dim-2
245 args: -dm_plex_gmsh_spacedim 2
246 args: -msh {{B2tri B2qua}separate_output}
247 args: -dim 2 -integral 2.000000000000000 # 2
248 test:
249 suffix: dim-2_msh-B3tet
250 args: -dm_plex_gmsh_spacedim 3
251 args: -msh B3tet -dim 2 -integral 9.914478
252 test:
253 suffix: dim-2_msh-B3hex
254 args: -dm_plex_gmsh_spacedim 3
255 args: -msh B3hex -dim 2 -integral 8.000000
256 test:
257 suffix: dim-3_msh-B3tet
258 args: -dm_plex_gmsh_spacedim 3
259 args: -msh B3tet -dim 3 -integral 2.666649
260 test:
261 suffix: dim-3_msh-B3hex
262 args: -dm_plex_gmsh_spacedim 3
263 args: -msh B3hex -dim 3 -integral 1.539600
264
265 TEST*/
266