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) {ierr = PetscViewerRead(viewer, line, 3, NULL, PETSC_STRING);CHKERRQ(ierr);} 145 ierr = PetscViewerRead(viewer, line, 1, NULL, PETSC_STRING);CHKERRQ(ierr); 146 ierr = PetscStrncmp(line, "$EndPhysicalNames", PETSC_MAX_PATH_LEN, &match);CHKERRQ(ierr); 147 if (!match) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Gmsh file"); 148 /* Initial read for vertex section */ 149 ierr = PetscViewerRead(viewer, line, 1, NULL, PETSC_STRING);CHKERRQ(ierr); 150 } 151 /* Read vertices */ 152 ierr = PetscStrncmp(line, "$Nodes", 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 ierr = PetscViewerRead(viewer, line, 1, NULL, PETSC_STRING);CHKERRQ(ierr); 155 snum = sscanf(line, "%d", &numVertices); 156 if (snum != 1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Gmsh file"); 157 ierr = PetscMalloc1(numVertices*3, &coordsIn);CHKERRQ(ierr); 158 if (binary) { 159 size_t doubleSize, intSize; 160 PetscInt elementSize; 161 char *buffer; 162 PetscScalar *baseptr; 163 ierr = PetscDataTypeGetSize(PETSC_ENUM, &intSize);CHKERRQ(ierr); 164 ierr = PetscDataTypeGetSize(PETSC_DOUBLE, &doubleSize);CHKERRQ(ierr); 165 elementSize = (intSize + 3*doubleSize); 166 ierr = PetscMalloc1(elementSize*numVertices, &buffer);CHKERRQ(ierr); 167 ierr = PetscViewerRead(viewer, buffer, elementSize*numVertices, NULL, PETSC_CHAR);CHKERRQ(ierr); 168 if (bswap) ierr = PetscByteSwap(buffer, PETSC_CHAR, elementSize*numVertices);CHKERRQ(ierr); 169 for (v = 0; v < numVertices; ++v) { 170 baseptr = ((PetscScalar*)(buffer+v*elementSize+intSize)); 171 coordsIn[v*3+0] = baseptr[0]; 172 coordsIn[v*3+1] = baseptr[1]; 173 coordsIn[v*3+2] = baseptr[2]; 174 } 175 ierr = PetscFree(buffer);CHKERRQ(ierr); 176 } else { 177 for (v = 0; v < numVertices; ++v) { 178 ierr = PetscViewerRead(viewer, &i, 1, NULL, PETSC_ENUM);CHKERRQ(ierr); 179 ierr = PetscViewerRead(viewer, &(coordsIn[v*3]), 3, NULL, PETSC_DOUBLE);CHKERRQ(ierr); 180 if (i != (int)v+1) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid node number %d should be %d", i, v+1); 181 } 182 } 183 ierr = PetscViewerRead(viewer, line, 1, NULL, PETSC_STRING);CHKERRQ(ierr);; 184 ierr = PetscStrncmp(line, "$EndNodes", PETSC_MAX_PATH_LEN, &match);CHKERRQ(ierr); 185 if (!match) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Gmsh file"); 186 /* Read cells */ 187 ierr = PetscViewerRead(viewer, line, 1, NULL, PETSC_STRING);CHKERRQ(ierr);; 188 ierr = PetscStrncmp(line, "$Elements", PETSC_MAX_PATH_LEN, &match);CHKERRQ(ierr); 189 if (!match) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Gmsh file"); 190 ierr = PetscViewerRead(viewer, line, 1, NULL, PETSC_STRING);CHKERRQ(ierr);; 191 snum = sscanf(line, "%d", &numCells); 192 if (snum != 1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Gmsh file"); 193 } 194 195 if (!rank || binary) { 196 /* Gmsh elements can be of any dimension/co-dimension, so we need to traverse the 197 file contents multiple times to figure out the true number of cells and facets 198 in the given mesh. To make this more efficient we read the file contents only 199 once and store them in memory, while determining the true number of cells. */ 200 ierr = DMPlexCreateGmsh_ReadElement(viewer, numCells, binary, bswap, &gmsh_elem);CHKERRQ(ierr); 201 for (trueNumCells=0, c = 0; c < numCells; ++c) { 202 if (gmsh_elem[c].dim > dim) {dim = gmsh_elem[c].dim; trueNumCells = 0;} 203 if (gmsh_elem[c].dim == dim) trueNumCells++; 204 } 205 ierr = PetscViewerRead(viewer, line, 1, NULL, PETSC_STRING);CHKERRQ(ierr);; 206 ierr = PetscStrncmp(line, "$EndElements", PETSC_MAX_PATH_LEN, &match);CHKERRQ(ierr); 207 if (!match) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Gmsh file"); 208 } 209 /* For binary we read on all ranks, but only build the plex on rank 0 */ 210 if (binary && rank) {trueNumCells = 0; numVertices = 0;}; 211 /* Allocate the cell-vertex mesh */ 212 ierr = DMPlexSetChart(*dm, 0, trueNumCells+numVertices);CHKERRQ(ierr); 213 if (!rank) { 214 for (cell = 0, c = 0; c < numCells; ++c) { 215 if (gmsh_elem[c].dim == dim) { 216 ierr = DMPlexSetConeSize(*dm, cell, gmsh_elem[c].numNodes);CHKERRQ(ierr); 217 cell++; 218 } 219 } 220 } 221 ierr = DMSetUp(*dm);CHKERRQ(ierr); 222 /* Add cell-vertex connections */ 223 if (!rank) { 224 PetscInt pcone[8], corner; 225 for (cell = 0, c = 0; c < numCells; ++c) { 226 if (gmsh_elem[c].dim == dim) { 227 for (corner = 0; corner < gmsh_elem[c].numNodes; ++corner) { 228 pcone[corner] = gmsh_elem[c].nodes[corner] + trueNumCells-1; 229 } 230 if (dim == 3) { 231 /* Tetrahedra are inverted */ 232 if (gmsh_elem[c].numNodes == 4) { 233 PetscInt tmp = pcone[0]; 234 pcone[0] = pcone[1]; 235 pcone[1] = tmp; 236 } 237 /* Hexahedra are inverted */ 238 if (gmsh_elem[c].numNodes == 8) { 239 PetscInt tmp = pcone[1]; 240 pcone[1] = pcone[3]; 241 pcone[3] = tmp; 242 } 243 } 244 ierr = DMPlexSetCone(*dm, cell, pcone);CHKERRQ(ierr); 245 cell++; 246 } 247 } 248 } 249 ierr = MPI_Bcast(&dim, 1, MPIU_INT, 0, comm);CHKERRQ(ierr); 250 ierr = DMSetDimension(*dm, dim);CHKERRQ(ierr); 251 ierr = DMPlexSymmetrize(*dm);CHKERRQ(ierr); 252 ierr = DMPlexStratify(*dm);CHKERRQ(ierr); 253 if (interpolate) { 254 DM idm = NULL; 255 256 ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 257 ierr = DMDestroy(dm);CHKERRQ(ierr); 258 *dm = idm; 259 } 260 261 if (!rank) { 262 /* Apply boundary IDs by finding the relevant facets with vertex joins */ 263 PetscInt pcone[8], corner, vStart, vEnd; 264 265 ierr = DMPlexGetDepthStratum(*dm, 0, &vStart, &vEnd);CHKERRQ(ierr); 266 for (c = 0; c < numCells; ++c) { 267 if (gmsh_elem[c].dim == dim-1) { 268 PetscInt joinSize; 269 const PetscInt *join; 270 for (corner = 0; corner < gmsh_elem[c].numNodes; ++corner) { 271 pcone[corner] = gmsh_elem[c].nodes[corner] + vStart - 1; 272 } 273 ierr = DMPlexGetFullJoin(*dm, gmsh_elem[c].numNodes, (const PetscInt *) pcone, &joinSize, &join);CHKERRQ(ierr); 274 if (joinSize != 1) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Could not determine Plex facet for element %d", gmsh_elem[c].id); 275 ierr = DMSetLabelValue(*dm, "Face Sets", join[0], gmsh_elem[c].tags[0]);CHKERRQ(ierr); 276 ierr = DMPlexRestoreJoin(*dm, gmsh_elem[c].numNodes, (const PetscInt *) pcone, &joinSize, &join);CHKERRQ(ierr); 277 } 278 } 279 280 /* Create cell sets */ 281 for (cell = 0, c = 0; c < numCells; ++c) { 282 if (gmsh_elem[c].dim == dim) { 283 if (gmsh_elem[c].numTags > 0) { 284 ierr = DMSetLabelValue(*dm, "Cell Sets", cell, gmsh_elem[c].tags[0]);CHKERRQ(ierr); 285 cell++; 286 } 287 } 288 } 289 290 /* Create vertex sets */ 291 for (c = 0; c < numCells; ++c) { 292 if (gmsh_elem[c].dim == 0) { 293 if (gmsh_elem[c].numTags > 0) { 294 ierr = DMSetLabelValue(*dm, "Vertex Sets", gmsh_elem[c].nodes[0] + vStart - 1, gmsh_elem[c].tags[0]);CHKERRQ(ierr); 295 } 296 } 297 } 298 } 299 300 /* Read coordinates */ 301 ierr = DMGetCoordinateSection(*dm, &coordSection);CHKERRQ(ierr); 302 ierr = PetscSectionSetNumFields(coordSection, 1);CHKERRQ(ierr); 303 ierr = PetscSectionSetFieldComponents(coordSection, 0, dim);CHKERRQ(ierr); 304 ierr = PetscSectionSetChart(coordSection, trueNumCells, trueNumCells + numVertices);CHKERRQ(ierr); 305 for (v = trueNumCells; v < trueNumCells+numVertices; ++v) { 306 ierr = PetscSectionSetDof(coordSection, v, dim);CHKERRQ(ierr); 307 ierr = PetscSectionSetFieldDof(coordSection, v, 0, dim);CHKERRQ(ierr); 308 } 309 ierr = PetscSectionSetUp(coordSection);CHKERRQ(ierr); 310 ierr = PetscSectionGetStorageSize(coordSection, &coordSize);CHKERRQ(ierr); 311 ierr = VecCreate(PETSC_COMM_SELF, &coordinates);CHKERRQ(ierr); 312 ierr = PetscObjectSetName((PetscObject) coordinates, "coordinates");CHKERRQ(ierr); 313 ierr = VecSetSizes(coordinates, coordSize, PETSC_DETERMINE);CHKERRQ(ierr); 314 ierr = VecSetBlockSize(coordinates, dim);CHKERRQ(ierr); 315 ierr = VecSetType(coordinates, VECSTANDARD);CHKERRQ(ierr); 316 ierr = VecGetArray(coordinates, &coords);CHKERRQ(ierr); 317 if (!rank) { 318 for (v = 0; v < numVertices; ++v) { 319 for (d = 0; d < dim; ++d) { 320 coords[v*dim+d] = coordsIn[v*3+d]; 321 } 322 } 323 } 324 ierr = VecRestoreArray(coordinates, &coords);CHKERRQ(ierr); 325 ierr = PetscFree(coordsIn);CHKERRQ(ierr); 326 ierr = DMSetCoordinatesLocal(*dm, coordinates);CHKERRQ(ierr); 327 ierr = VecDestroy(&coordinates);CHKERRQ(ierr); 328 /* Clean up intermediate storage */ 329 if (!rank || binary) ierr = PetscFree(gmsh_elem);CHKERRQ(ierr); 330 ierr = PetscLogEventEnd(DMPLEX_CreateGmsh,*dm,0,0,0);CHKERRQ(ierr); 331 PetscFunctionReturn(0); 332 } 333 334 #undef __FUNCT__ 335 #define __FUNCT__ "DMPlexCreateGmsh_ReadElement" 336 PetscErrorCode DMPlexCreateGmsh_ReadElement(PetscViewer viewer, PetscInt numCells, PetscBool binary, PetscBool byteSwap, GmshElement **gmsh_elems) 337 { 338 PetscInt c, p; 339 GmshElement *elements; 340 int i, cellType, dim, numNodes, numElem, numTags; 341 int ibuf[16]; 342 PetscErrorCode ierr; 343 344 PetscFunctionBegin; 345 ierr = PetscMalloc1(numCells, &elements);CHKERRQ(ierr); 346 for (c = 0; c < numCells;) { 347 ierr = PetscViewerRead(viewer, &ibuf, 3, NULL, PETSC_ENUM);CHKERRQ(ierr); 348 if (byteSwap) ierr = PetscByteSwap(&ibuf, PETSC_ENUM, 3);CHKERRQ(ierr); 349 if (binary) { 350 cellType = ibuf[0]; 351 numElem = ibuf[1]; 352 numTags = ibuf[2]; 353 } else { 354 elements[c].id = ibuf[0]; 355 cellType = ibuf[1]; 356 numTags = ibuf[2]; 357 numElem = 1; 358 } 359 switch (cellType) { 360 case 1: /* 2-node line */ 361 dim = 1; 362 numNodes = 2; 363 break; 364 case 2: /* 3-node triangle */ 365 dim = 2; 366 numNodes = 3; 367 break; 368 case 3: /* 4-node quadrangle */ 369 dim = 2; 370 numNodes = 4; 371 break; 372 case 4: /* 4-node tetrahedron */ 373 dim = 3; 374 numNodes = 4; 375 break; 376 case 5: /* 8-node hexahedron */ 377 dim = 3; 378 numNodes = 8; 379 break; 380 case 15: /* 1-node vertex */ 381 dim = 0; 382 numNodes = 1; 383 break; 384 default: 385 SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Unsupported Gmsh element type %d", cellType); 386 } 387 if (binary) { 388 const PetscInt nint = numNodes + numTags + 1; 389 for (i = 0; i < numElem; ++i, ++c) { 390 /* Loop over inner binary element block */ 391 elements[c].dim = dim; 392 elements[c].numNodes = numNodes; 393 elements[c].numTags = numTags; 394 395 ierr = PetscViewerRead(viewer, &ibuf, nint, NULL, PETSC_ENUM);CHKERRQ(ierr); 396 if (byteSwap) ierr = PetscByteSwap( &ibuf, PETSC_ENUM, nint);CHKERRQ(ierr); 397 elements[c].id = ibuf[0]; 398 for (p = 0; p < numTags; p++) elements[c].tags[p] = ibuf[1 + p]; 399 for (p = 0; p < numNodes; p++) elements[c].nodes[p] = ibuf[1 + numTags + p]; 400 } 401 } else { 402 elements[c].dim = dim; 403 elements[c].numNodes = numNodes; 404 elements[c].numTags = numTags; 405 ierr = PetscViewerRead(viewer, elements[c].tags, elements[c].numTags, NULL, PETSC_ENUM);CHKERRQ(ierr); 406 ierr = PetscViewerRead(viewer, elements[c].nodes, elements[c].numNodes, NULL, PETSC_ENUM);CHKERRQ(ierr); 407 c++; 408 } 409 } 410 *gmsh_elems = elements; 411 PetscFunctionReturn(0); 412 } 413