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