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