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 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 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 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 testset: 194 suffix: B2 # 2D ball 195 requires: define(PETSC_GMSH_EXE) 196 args: -dir ${wPETSC_DIR}/share/petsc/datafiles/meshes 197 args: -msh {{B2tri B2qua}} 198 args: -dim 2 -integral 3.141592653589793 # pi 199 args: -order {{2 3 4 5 6 7 8 9}} -tol 0.05 200 201 testset: 202 suffix: B2_bnd # 2D ball boundary 203 requires: define(PETSC_GMSH_EXE) 204 args: -dir ${wPETSC_DIR}/share/petsc/datafiles/meshes 205 args: -dm_plex_gmsh_spacedim 2 206 args: -msh {{B2tri B2qua}} 207 args: -dim 1 -integral 6.283185307179586 # 2*pi 208 args: -order {{2 3 4 5 6 7 8 9}} -tol 0.05 209 210 testset: 211 suffix: B3 # 3D ball 212 requires: define(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 requires: define(PETSC_GMSH_EXE) 221 args: -dir ${wPETSC_DIR}/share/petsc/datafiles/meshes 222 args: -dm_plex_gmsh_spacedim 3 223 args: -msh {{B3tet B3hex}} 224 args: -dim 2 -integral 12.566370614359172 # 4*pi 225 args: -order {{2 3 4 5 6 7 8 9}} -tol 0.25 226 227 testset: 228 suffix: B_lin # linear discretizations 229 requires: define(PETSC_GMSH_EXE) 230 args: -dir ${wPETSC_DIR}/share/petsc/datafiles/meshes 231 args: -dm_plex_gmsh_highorder true 232 args: -dm_plex_gmsh_project true 233 args: -dm_plex_gmsh_project_petscspace_degree {{1 2 3}separate_output} 234 args: -dm_plex_gmsh_fe_view 235 args: -dm_plex_gmsh_project_fe_view 236 args: -order 1 -tol 1e-4 237 test: 238 suffix: dim-1 239 args: -dm_plex_gmsh_spacedim 2 240 args: -msh {{B2tri B2qua}separate_output} 241 args: -dim 1 -integral 5.656854249492381 # 4*sqrt(2) 242 test: 243 suffix: dim-2 244 args: -dm_plex_gmsh_spacedim 2 245 args: -msh {{B2tri B2qua}separate_output} 246 args: -dim 2 -integral 2.000000000000000 # 2 247 test: 248 suffix: dim-2_msh-B3tet 249 args: -dm_plex_gmsh_spacedim 3 250 args: -msh B3tet -dim 2 -integral 9.914478 251 test: 252 suffix: dim-2_msh-B3hex 253 args: -dm_plex_gmsh_spacedim 3 254 args: -msh B3hex -dim 2 -integral 8.000000 255 test: 256 suffix: dim-3_msh-B3tet 257 args: -dm_plex_gmsh_spacedim 3 258 args: -msh B3tet -dim 3 -integral 2.666649 259 test: 260 suffix: dim-3_msh-B3hex 261 args: -dm_plex_gmsh_spacedim 3 262 args: -msh B3hex -dim 3 -integral 1.539600 263 264 TEST*/ 265