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