1 #define PETSCDM_DLL 2 #include <petsc/private/dmpleximpl.h> /*I "petscdmplex.h" I*/ 3 4 /*@C 5 DMPlexCreateGmshFromFile - Create a DMPlex mesh from a Gmsh file 6 7 + comm - The MPI communicator 8 . filename - Name of the Gmsh file 9 - interpolate - Create faces and edges in the mesh 10 11 Output Parameter: 12 . dm - The DM object representing the mesh 13 14 Level: beginner 15 16 .seealso: DMPlexCreateFromFile(), DMPlexCreateGmsh(), DMPlexCreate() 17 @*/ 18 PetscErrorCode DMPlexCreateGmshFromFile(MPI_Comm comm, const char filename[], PetscBool interpolate, DM *dm) 19 { 20 PetscViewer viewer, vheader; 21 PetscMPIInt rank; 22 PetscViewerType vtype; 23 char line[PETSC_MAX_PATH_LEN]; 24 int snum; 25 PetscBool match; 26 int fT; 27 PetscInt fileType; 28 float version; 29 PetscErrorCode ierr; 30 31 PetscFunctionBegin; 32 ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 33 /* Determine Gmsh file type (ASCII or binary) from file header */ 34 ierr = PetscViewerCreate(comm, &vheader);CHKERRQ(ierr); 35 ierr = PetscViewerSetType(vheader, PETSCVIEWERASCII);CHKERRQ(ierr); 36 ierr = PetscViewerFileSetMode(vheader, FILE_MODE_READ);CHKERRQ(ierr); 37 ierr = PetscViewerFileSetName(vheader, filename);CHKERRQ(ierr); 38 if (!rank) { 39 /* Read only the first two lines of the Gmsh file */ 40 ierr = PetscViewerRead(vheader, line, 1, NULL, PETSC_STRING);CHKERRQ(ierr); 41 ierr = PetscStrncmp(line, "$MeshFormat", PETSC_MAX_PATH_LEN, &match);CHKERRQ(ierr); 42 if (!match) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Gmsh file"); 43 ierr = PetscViewerRead(vheader, line, 2, NULL, PETSC_STRING);CHKERRQ(ierr); 44 snum = sscanf(line, "%f %d", &version, &fT); 45 fileType = (PetscInt) fT; 46 if (snum != 2) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Unable to parse Gmsh file header: %s", line); 47 if (version < 2.0) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Gmsh file must be at least version 2.0"); 48 } 49 ierr = MPI_Bcast(&fileType, 1, MPIU_INT, 0, comm);CHKERRQ(ierr); 50 /* Create appropriate viewer and build plex */ 51 if (fileType == 0) vtype = PETSCVIEWERASCII; 52 else vtype = PETSCVIEWERBINARY; 53 ierr = PetscViewerCreate(comm, &viewer);CHKERRQ(ierr); 54 ierr = PetscViewerSetType(viewer, vtype);CHKERRQ(ierr); 55 ierr = PetscViewerFileSetMode(viewer, FILE_MODE_READ);CHKERRQ(ierr); 56 ierr = PetscViewerFileSetName(viewer, filename);CHKERRQ(ierr); 57 ierr = DMPlexCreateGmsh(comm, viewer, interpolate, dm);CHKERRQ(ierr); 58 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 59 ierr = PetscViewerDestroy(&vheader);CHKERRQ(ierr); 60 PetscFunctionReturn(0); 61 } 62 63 static PetscErrorCode DMPlexCreateGmsh_ReadElement(PetscViewer viewer, PetscInt numCells, PetscBool binary, PetscBool byteSwap, GmshElement **gmsh_elems) 64 { 65 PetscInt c, p; 66 GmshElement *elements; 67 int i, cellType, dim, numNodes, numElem, numTags; 68 int ibuf[16]; 69 PetscErrorCode ierr; 70 71 PetscFunctionBegin; 72 ierr = PetscMalloc1(numCells, &elements);CHKERRQ(ierr); 73 for (c = 0; c < numCells;) { 74 ierr = PetscViewerRead(viewer, &ibuf, 3, NULL, PETSC_ENUM);CHKERRQ(ierr); 75 if (byteSwap) ierr = PetscByteSwap(&ibuf, PETSC_ENUM, 3);CHKERRQ(ierr); 76 if (binary) { 77 cellType = ibuf[0]; 78 numElem = ibuf[1]; 79 numTags = ibuf[2]; 80 } else { 81 elements[c].id = ibuf[0]; 82 cellType = ibuf[1]; 83 numTags = ibuf[2]; 84 numElem = 1; 85 } 86 switch (cellType) { 87 case 1: /* 2-node line */ 88 dim = 1; 89 numNodes = 2; 90 break; 91 case 2: /* 3-node triangle */ 92 dim = 2; 93 numNodes = 3; 94 break; 95 case 3: /* 4-node quadrangle */ 96 dim = 2; 97 numNodes = 4; 98 break; 99 case 4: /* 4-node tetrahedron */ 100 dim = 3; 101 numNodes = 4; 102 break; 103 case 5: /* 8-node hexahedron */ 104 dim = 3; 105 numNodes = 8; 106 break; 107 case 15: /* 1-node vertex */ 108 dim = 0; 109 numNodes = 1; 110 break; 111 default: 112 SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Unsupported Gmsh element type %d", cellType); 113 } 114 if (binary) { 115 const PetscInt nint = numNodes + numTags + 1; 116 for (i = 0; i < numElem; ++i, ++c) { 117 /* Loop over inner binary element block */ 118 elements[c].dim = dim; 119 elements[c].numNodes = numNodes; 120 elements[c].numTags = numTags; 121 122 ierr = PetscViewerRead(viewer, &ibuf, nint, NULL, PETSC_ENUM);CHKERRQ(ierr); 123 if (byteSwap) ierr = PetscByteSwap( &ibuf, PETSC_ENUM, nint);CHKERRQ(ierr); 124 elements[c].id = ibuf[0]; 125 for (p = 0; p < numTags; p++) elements[c].tags[p] = ibuf[1 + p]; 126 for (p = 0; p < numNodes; p++) elements[c].nodes[p] = ibuf[1 + numTags + p]; 127 } 128 } else { 129 elements[c].dim = dim; 130 elements[c].numNodes = numNodes; 131 elements[c].numTags = numTags; 132 ierr = PetscViewerRead(viewer, elements[c].tags, elements[c].numTags, NULL, PETSC_ENUM);CHKERRQ(ierr); 133 ierr = PetscViewerRead(viewer, elements[c].nodes, elements[c].numNodes, NULL, PETSC_ENUM);CHKERRQ(ierr); 134 c++; 135 } 136 } 137 *gmsh_elems = elements; 138 PetscFunctionReturn(0); 139 } 140 141 /*@ 142 DMPlexCreateGmsh - Create a DMPlex mesh from a Gmsh file viewer 143 144 Collective on comm 145 146 Input Parameters: 147 + comm - The MPI communicator 148 . viewer - The Viewer associated with a Gmsh file 149 - interpolate - Create faces and edges in the mesh 150 151 Output Parameter: 152 . dm - The DM object representing the mesh 153 154 Note: http://www.geuz.org/gmsh/doc/texinfo/#MSH-ASCII-file-format 155 and http://www.geuz.org/gmsh/doc/texinfo/#MSH-binary-file-format 156 157 Level: beginner 158 159 .keywords: mesh,Gmsh 160 .seealso: DMPLEX, DMCreate() 161 @*/ 162 PetscErrorCode DMPlexCreateGmsh(MPI_Comm comm, PetscViewer viewer, PetscBool interpolate, DM *dm) 163 { 164 PetscViewerType vtype; 165 GmshElement *gmsh_elem; 166 PetscSection coordSection; 167 Vec coordinates; 168 PetscScalar *coords, *coordsIn = NULL; 169 PetscInt dim = 0, coordSize, c, v, d, r, cell; 170 int i, numVertices = 0, numCells = 0, trueNumCells = 0, numRegions = 0, snum; 171 PetscMPIInt num_proc, rank; 172 char line[PETSC_MAX_PATH_LEN]; 173 PetscBool match, binary, bswap = PETSC_FALSE; 174 PetscErrorCode ierr; 175 176 PetscFunctionBegin; 177 ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 178 ierr = MPI_Comm_size(comm, &num_proc);CHKERRQ(ierr); 179 ierr = DMCreate(comm, dm);CHKERRQ(ierr); 180 ierr = DMSetType(*dm, DMPLEX);CHKERRQ(ierr); 181 ierr = PetscLogEventBegin(DMPLEX_CreateGmsh,*dm,0,0,0);CHKERRQ(ierr); 182 ierr = PetscViewerGetType(viewer, &vtype);CHKERRQ(ierr); 183 ierr = PetscStrcmp(vtype, PETSCVIEWERBINARY, &binary);CHKERRQ(ierr); 184 if (!rank || binary) { 185 PetscBool match; 186 int fileType, dataSize; 187 float version; 188 189 /* Read header */ 190 ierr = PetscViewerRead(viewer, line, 1, NULL, PETSC_STRING);CHKERRQ(ierr); 191 ierr = PetscStrncmp(line, "$MeshFormat", PETSC_MAX_PATH_LEN, &match);CHKERRQ(ierr); 192 if (!match) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Gmsh file"); 193 ierr = PetscViewerRead(viewer, line, 3, NULL, PETSC_STRING);CHKERRQ(ierr); 194 snum = sscanf(line, "%f %d %d", &version, &fileType, &dataSize); 195 if (snum != 3) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Unable to parse Gmsh file header: %s", line); 196 if (version < 2.0) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Gmsh file must be at least version 2.0"); 197 if (dataSize != sizeof(double)) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Data size %d is not valid for a Gmsh file", dataSize); 198 if (binary) { 199 int checkInt; 200 ierr = PetscViewerRead(viewer, &checkInt, 1, NULL, PETSC_ENUM);CHKERRQ(ierr); 201 if (checkInt != 1) { 202 ierr = PetscByteSwap(&checkInt, PETSC_ENUM, 1);CHKERRQ(ierr); 203 if (checkInt == 1) bswap = PETSC_TRUE; 204 else SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File type %d is not a valid Gmsh binary file", fileType); 205 } 206 } else if (fileType) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File type %d is not a valid Gmsh ASCII file", fileType); 207 ierr = PetscViewerRead(viewer, line, 1, NULL, PETSC_STRING);CHKERRQ(ierr); 208 ierr = PetscStrncmp(line, "$EndMeshFormat", PETSC_MAX_PATH_LEN, &match);CHKERRQ(ierr); 209 if (!match) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Gmsh file"); 210 /* OPTIONAL Read physical names */ 211 ierr = PetscViewerRead(viewer, line, 1, NULL, PETSC_STRING);CHKERRQ(ierr); 212 ierr = PetscStrncmp(line, "$PhysicalNames", PETSC_MAX_PATH_LEN, &match);CHKERRQ(ierr); 213 if (match) { 214 ierr = PetscViewerRead(viewer, line, 1, NULL, PETSC_STRING);CHKERRQ(ierr); 215 snum = sscanf(line, "%d", &numRegions); 216 if (snum != 1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Gmsh file"); 217 for (r = 0; r < numRegions; ++r) {ierr = PetscViewerRead(viewer, line, 3, NULL, PETSC_STRING);CHKERRQ(ierr);} 218 ierr = PetscViewerRead(viewer, line, 1, NULL, PETSC_STRING);CHKERRQ(ierr); 219 ierr = PetscStrncmp(line, "$EndPhysicalNames", PETSC_MAX_PATH_LEN, &match);CHKERRQ(ierr); 220 if (!match) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Gmsh file"); 221 /* Initial read for vertex section */ 222 ierr = PetscViewerRead(viewer, line, 1, NULL, PETSC_STRING);CHKERRQ(ierr); 223 } 224 /* Read vertices */ 225 ierr = PetscStrncmp(line, "$Nodes", PETSC_MAX_PATH_LEN, &match);CHKERRQ(ierr); 226 if (!match) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Gmsh file"); 227 ierr = PetscViewerRead(viewer, line, 1, NULL, PETSC_STRING);CHKERRQ(ierr); 228 snum = sscanf(line, "%d", &numVertices); 229 if (snum != 1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Gmsh file"); 230 ierr = PetscMalloc1(numVertices*3, &coordsIn);CHKERRQ(ierr); 231 if (binary) { 232 size_t doubleSize, intSize; 233 PetscInt elementSize; 234 char *buffer; 235 PetscScalar *baseptr; 236 ierr = PetscDataTypeGetSize(PETSC_ENUM, &intSize);CHKERRQ(ierr); 237 ierr = PetscDataTypeGetSize(PETSC_DOUBLE, &doubleSize);CHKERRQ(ierr); 238 elementSize = (intSize + 3*doubleSize); 239 ierr = PetscMalloc1(elementSize*numVertices, &buffer);CHKERRQ(ierr); 240 ierr = PetscViewerRead(viewer, buffer, elementSize*numVertices, NULL, PETSC_CHAR);CHKERRQ(ierr); 241 if (bswap) ierr = PetscByteSwap(buffer, PETSC_CHAR, elementSize*numVertices);CHKERRQ(ierr); 242 for (v = 0; v < numVertices; ++v) { 243 baseptr = ((PetscScalar*)(buffer+v*elementSize+intSize)); 244 coordsIn[v*3+0] = baseptr[0]; 245 coordsIn[v*3+1] = baseptr[1]; 246 coordsIn[v*3+2] = baseptr[2]; 247 } 248 ierr = PetscFree(buffer);CHKERRQ(ierr); 249 } else { 250 for (v = 0; v < numVertices; ++v) { 251 ierr = PetscViewerRead(viewer, &i, 1, NULL, PETSC_ENUM);CHKERRQ(ierr); 252 ierr = PetscViewerRead(viewer, &(coordsIn[v*3]), 3, NULL, PETSC_DOUBLE);CHKERRQ(ierr); 253 if (i != (int)v+1) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid node number %d should be %d", i, v+1); 254 } 255 } 256 ierr = PetscViewerRead(viewer, line, 1, NULL, PETSC_STRING);CHKERRQ(ierr);; 257 ierr = PetscStrncmp(line, "$EndNodes", PETSC_MAX_PATH_LEN, &match);CHKERRQ(ierr); 258 if (!match) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Gmsh file"); 259 /* Read cells */ 260 ierr = PetscViewerRead(viewer, line, 1, NULL, PETSC_STRING);CHKERRQ(ierr);; 261 ierr = PetscStrncmp(line, "$Elements", PETSC_MAX_PATH_LEN, &match);CHKERRQ(ierr); 262 if (!match) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Gmsh file"); 263 ierr = PetscViewerRead(viewer, line, 1, NULL, PETSC_STRING);CHKERRQ(ierr);; 264 snum = sscanf(line, "%d", &numCells); 265 if (snum != 1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Gmsh file"); 266 } 267 268 if (!rank || binary) { 269 /* Gmsh elements can be of any dimension/co-dimension, so we need to traverse the 270 file contents multiple times to figure out the true number of cells and facets 271 in the given mesh. To make this more efficient we read the file contents only 272 once and store them in memory, while determining the true number of cells. */ 273 ierr = DMPlexCreateGmsh_ReadElement(viewer, numCells, binary, bswap, &gmsh_elem);CHKERRQ(ierr); 274 for (trueNumCells=0, c = 0; c < numCells; ++c) { 275 if (gmsh_elem[c].dim > dim) {dim = gmsh_elem[c].dim; trueNumCells = 0;} 276 if (gmsh_elem[c].dim == dim) trueNumCells++; 277 } 278 ierr = PetscViewerRead(viewer, line, 1, NULL, PETSC_STRING);CHKERRQ(ierr);; 279 ierr = PetscStrncmp(line, "$EndElements", 12, &match);CHKERRQ(ierr); 280 if (!match) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Gmsh file"); 281 } 282 /* For binary we read on all ranks, but only build the plex on rank 0 */ 283 if (binary && rank) {trueNumCells = 0; numVertices = 0;}; 284 /* Allocate the cell-vertex mesh */ 285 ierr = DMPlexSetChart(*dm, 0, trueNumCells+numVertices);CHKERRQ(ierr); 286 if (!rank) { 287 for (cell = 0, c = 0; c < numCells; ++c) { 288 if (gmsh_elem[c].dim == dim) { 289 ierr = DMPlexSetConeSize(*dm, cell, gmsh_elem[c].numNodes);CHKERRQ(ierr); 290 cell++; 291 } 292 } 293 } 294 ierr = DMSetUp(*dm);CHKERRQ(ierr); 295 /* Add cell-vertex connections */ 296 if (!rank) { 297 PetscInt pcone[8], corner; 298 for (cell = 0, c = 0; c < numCells; ++c) { 299 if (gmsh_elem[c].dim == dim) { 300 for (corner = 0; corner < gmsh_elem[c].numNodes; ++corner) { 301 pcone[corner] = gmsh_elem[c].nodes[corner] + trueNumCells-1; 302 } 303 if (dim == 3) { 304 /* Tetrahedra are inverted */ 305 if (gmsh_elem[c].numNodes == 4) { 306 PetscInt tmp = pcone[0]; 307 pcone[0] = pcone[1]; 308 pcone[1] = tmp; 309 } 310 /* Hexahedra are inverted */ 311 if (gmsh_elem[c].numNodes == 8) { 312 PetscInt tmp = pcone[1]; 313 pcone[1] = pcone[3]; 314 pcone[3] = tmp; 315 } 316 } 317 ierr = DMPlexSetCone(*dm, cell, pcone);CHKERRQ(ierr); 318 cell++; 319 } 320 } 321 } 322 ierr = MPI_Bcast(&dim, 1, MPIU_INT, 0, comm);CHKERRQ(ierr); 323 ierr = DMSetDimension(*dm, dim);CHKERRQ(ierr); 324 ierr = DMPlexSymmetrize(*dm);CHKERRQ(ierr); 325 ierr = DMPlexStratify(*dm);CHKERRQ(ierr); 326 if (interpolate) { 327 DM idm = NULL; 328 329 ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 330 ierr = DMDestroy(dm);CHKERRQ(ierr); 331 *dm = idm; 332 } 333 334 if (!rank) { 335 /* Apply boundary IDs by finding the relevant facets with vertex joins */ 336 PetscInt pcone[8], corner, vStart, vEnd; 337 338 ierr = DMPlexGetDepthStratum(*dm, 0, &vStart, &vEnd);CHKERRQ(ierr); 339 for (c = 0; c < numCells; ++c) { 340 if (gmsh_elem[c].dim == dim-1) { 341 PetscInt joinSize; 342 const PetscInt *join; 343 for (corner = 0; corner < gmsh_elem[c].numNodes; ++corner) { 344 pcone[corner] = gmsh_elem[c].nodes[corner] + vStart - 1; 345 } 346 ierr = DMPlexGetFullJoin(*dm, gmsh_elem[c].numNodes, (const PetscInt *) pcone, &joinSize, &join);CHKERRQ(ierr); 347 if (joinSize != 1) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Could not determine Plex facet for element %d", gmsh_elem[c].id); 348 ierr = DMSetLabelValue(*dm, "Face Sets", join[0], gmsh_elem[c].tags[0]);CHKERRQ(ierr); 349 ierr = DMPlexRestoreJoin(*dm, gmsh_elem[c].numNodes, (const PetscInt *) pcone, &joinSize, &join);CHKERRQ(ierr); 350 } 351 } 352 353 /* Create cell sets */ 354 for (cell = 0, c = 0; c < numCells; ++c) { 355 if (gmsh_elem[c].dim == dim) { 356 if (gmsh_elem[c].numTags > 0) { 357 ierr = DMSetLabelValue(*dm, "Cell Sets", cell, gmsh_elem[c].tags[0]);CHKERRQ(ierr); 358 cell++; 359 } 360 } 361 } 362 363 /* Create vertex sets */ 364 for (c = 0; c < numCells; ++c) { 365 if (gmsh_elem[c].dim == 0) { 366 if (gmsh_elem[c].numTags > 0) { 367 ierr = DMSetLabelValue(*dm, "Vertex Sets", gmsh_elem[c].nodes[0] + vStart - 1, gmsh_elem[c].tags[0]);CHKERRQ(ierr); 368 } 369 } 370 } 371 } 372 373 /* Read coordinates */ 374 ierr = DMGetCoordinateSection(*dm, &coordSection);CHKERRQ(ierr); 375 ierr = PetscSectionSetNumFields(coordSection, 1);CHKERRQ(ierr); 376 ierr = PetscSectionSetFieldComponents(coordSection, 0, dim);CHKERRQ(ierr); 377 ierr = PetscSectionSetChart(coordSection, trueNumCells, trueNumCells + numVertices);CHKERRQ(ierr); 378 for (v = trueNumCells; v < trueNumCells+numVertices; ++v) { 379 ierr = PetscSectionSetDof(coordSection, v, dim);CHKERRQ(ierr); 380 ierr = PetscSectionSetFieldDof(coordSection, v, 0, dim);CHKERRQ(ierr); 381 } 382 ierr = PetscSectionSetUp(coordSection);CHKERRQ(ierr); 383 ierr = PetscSectionGetStorageSize(coordSection, &coordSize);CHKERRQ(ierr); 384 ierr = VecCreate(PETSC_COMM_SELF, &coordinates);CHKERRQ(ierr); 385 ierr = PetscObjectSetName((PetscObject) coordinates, "coordinates");CHKERRQ(ierr); 386 ierr = VecSetSizes(coordinates, coordSize, PETSC_DETERMINE);CHKERRQ(ierr); 387 ierr = VecSetBlockSize(coordinates, dim);CHKERRQ(ierr); 388 ierr = VecSetType(coordinates, VECSTANDARD);CHKERRQ(ierr); 389 ierr = VecGetArray(coordinates, &coords);CHKERRQ(ierr); 390 if (!rank) { 391 for (v = 0; v < numVertices; ++v) { 392 for (d = 0; d < dim; ++d) { 393 coords[v*dim+d] = coordsIn[v*3+d]; 394 } 395 } 396 } 397 ierr = VecRestoreArray(coordinates, &coords);CHKERRQ(ierr); 398 ierr = PetscFree(coordsIn);CHKERRQ(ierr); 399 ierr = DMSetCoordinatesLocal(*dm, coordinates);CHKERRQ(ierr); 400 ierr = VecDestroy(&coordinates);CHKERRQ(ierr); 401 /* Clean up intermediate storage */ 402 if (!rank || binary) ierr = PetscFree(gmsh_elem);CHKERRQ(ierr); 403 ierr = PetscLogEventEnd(DMPLEX_CreateGmsh,*dm,0,0,0);CHKERRQ(ierr); 404 PetscFunctionReturn(0); 405 } 406