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