1 #define PETSCDM_DLL 2 #include <petsc/private/dmpleximpl.h> /*I "petscdmplex.h" I*/ 3 4 /*@C 5 DMPlexCreateFluentFromFile - Create a DMPlex mesh from a Fluent mesh file 6 7 + comm - The MPI communicator 8 . filename - Name of the Fluent mesh 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(), DMPlexCreateFluent(), DMPlexCreate() 17 @*/ 18 PetscErrorCode DMPlexCreateFluentFromFile(MPI_Comm comm, const char filename[], PetscBool interpolate, DM *dm) 19 { 20 PetscViewer viewer; 21 PetscErrorCode ierr; 22 23 PetscFunctionBegin; 24 /* Create file viewer and build plex */ 25 ierr = PetscViewerCreate(comm, &viewer);CHKERRQ(ierr); 26 ierr = PetscViewerSetType(viewer, PETSCVIEWERASCII);CHKERRQ(ierr); 27 ierr = PetscViewerFileSetMode(viewer, FILE_MODE_READ);CHKERRQ(ierr); 28 ierr = PetscViewerFileSetName(viewer, filename);CHKERRQ(ierr); 29 ierr = DMPlexCreateFluent(comm, viewer, interpolate, dm);CHKERRQ(ierr); 30 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 31 PetscFunctionReturn(0); 32 } 33 34 PetscErrorCode DMPlexCreateFluent_ReadString(PetscViewer viewer, char *buffer, char delim) 35 { 36 PetscInt ret, i = 0; 37 PetscErrorCode ierr; 38 39 PetscFunctionBegin; 40 do {ierr = PetscViewerRead(viewer, &(buffer[i++]), 1, &ret, PETSC_CHAR);CHKERRQ(ierr);} 41 while (ret > 0 && buffer[i-1] != '\0' && buffer[i-1] != delim); 42 buffer[i] = '\0'; 43 PetscFunctionReturn(0); 44 } 45 46 PetscErrorCode DMPlexCreateFluent_ReadValues(PetscViewer viewer, void *data, PetscInt count, PetscDataType dtype, PetscBool binary) 47 { 48 int fdes=0; 49 FILE *file; 50 PetscInt i; 51 PetscErrorCode ierr; 52 53 PetscFunctionBegin; 54 if (binary) { 55 /* Extract raw file descriptor to read binary block */ 56 ierr = PetscViewerASCIIGetPointer(viewer, &file);CHKERRQ(ierr); 57 fflush(file); fdes = fileno(file); 58 } 59 60 if (!binary && dtype == PETSC_INT) { 61 char cbuf[256]; 62 int ibuf, snum; 63 /* Parse hexadecimal ascii integers */ 64 for (i = 0; i < count; i++) { 65 ierr = PetscViewerRead(viewer, cbuf, 1, NULL, PETSC_STRING);CHKERRQ(ierr); 66 snum = sscanf(cbuf, "%x", &ibuf); 67 if (snum != 1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Fluent file"); 68 ((PetscInt*)data)[i] = (PetscInt)ibuf; 69 } 70 } else if (binary && dtype == PETSC_INT) { 71 /* Always read 32-bit ints and cast to PetscInt */ 72 int *ibuf; 73 ierr = PetscMalloc1(count, &ibuf);CHKERRQ(ierr); 74 ierr = PetscBinaryRead(fdes, ibuf, count, PETSC_ENUM);CHKERRQ(ierr); 75 ierr = PetscByteSwap(ibuf, PETSC_ENUM, count);CHKERRQ(ierr); 76 for (i = 0; i < count; i++) ((PetscInt*)data)[i] = (PetscInt)(ibuf[i]); 77 ierr = PetscFree(ibuf);CHKERRQ(ierr); 78 79 } else if (binary && dtype == PETSC_SCALAR) { 80 float *fbuf; 81 /* Always read 32-bit floats and cast to PetscScalar */ 82 ierr = PetscMalloc1(count, &fbuf);CHKERRQ(ierr); 83 ierr = PetscBinaryRead(fdes, fbuf, count, PETSC_FLOAT);CHKERRQ(ierr); 84 ierr = PetscByteSwap(fbuf, PETSC_FLOAT, count);CHKERRQ(ierr); 85 for (i = 0; i < count; i++) ((PetscScalar*)data)[i] = (PetscScalar)(fbuf[i]); 86 ierr = PetscFree(fbuf);CHKERRQ(ierr); 87 } else { 88 ierr = PetscViewerASCIIRead(viewer, data, count, NULL, dtype);CHKERRQ(ierr); 89 } 90 PetscFunctionReturn(0); 91 } 92 93 PetscErrorCode DMPlexCreateFluent_ReadSection(PetscViewer viewer, FluentSection *s) 94 { 95 char buffer[PETSC_MAX_PATH_LEN]; 96 int snum; 97 PetscErrorCode ierr; 98 99 PetscFunctionBegin; 100 /* Fast-forward to next section and derive its index */ 101 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, '(');CHKERRQ(ierr); 102 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, ' ');CHKERRQ(ierr); 103 snum = sscanf(buffer, "%d", &(s->index)); 104 /* If we can't match an index return -1 to signal end-of-file */ 105 if (snum < 1) {s->index = -1; PetscFunctionReturn(0);} 106 107 if (s->index == 0) { /* Comment */ 108 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, ')');CHKERRQ(ierr); 109 110 } else if (s->index == 2) { /* Dimension */ 111 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, ')');CHKERRQ(ierr); 112 snum = sscanf(buffer, "%d", &(s->nd)); 113 if (snum != 1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Fluent file"); 114 115 } else if (s->index == 10 || s->index == 2010) { /* Vertices */ 116 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, ')');CHKERRQ(ierr); 117 snum = sscanf(buffer, "(%x %x %x %d %d)", &(s->zoneID), &(s->first), &(s->last), &(s->type), &(s->nd)); 118 if (snum != 5) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Fluent file"); 119 if (s->zoneID > 0) { 120 PetscInt numCoords = s->last - s->first + 1; 121 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, '(');CHKERRQ(ierr); 122 ierr = PetscMalloc1(s->nd*numCoords, (PetscScalar**)&s->data);CHKERRQ(ierr); 123 ierr = DMPlexCreateFluent_ReadValues(viewer, s->data, s->nd*numCoords, PETSC_SCALAR, s->index==2010 ? PETSC_TRUE : PETSC_FALSE);CHKERRQ(ierr); 124 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, ')');CHKERRQ(ierr); 125 } 126 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, ')');CHKERRQ(ierr); 127 128 } else if (s->index == 12 || s->index == 2012) { /* Cells */ 129 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, ')');CHKERRQ(ierr); 130 snum = sscanf(buffer, "(%x", &(s->zoneID)); 131 if (snum != 1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Fluent file"); 132 if (s->zoneID == 0) { /* Header section */ 133 snum = sscanf(buffer, "(%x %x %x %d)", &(s->zoneID), &(s->first), &(s->last), &(s->nd)); 134 if (snum != 4) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Fluent file"); 135 } else { /* Data section */ 136 snum = sscanf(buffer, "(%x %x %x %d %d)", &(s->zoneID), &(s->first), &(s->last), &(s->type), &(s->nd)); 137 if (snum != 5) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Fluent file"); 138 if (s->nd == 0) { 139 /* Read cell type definitions for mixed cells */ 140 PetscInt numCells = s->last - s->first + 1; 141 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, '(');CHKERRQ(ierr); 142 ierr = PetscMalloc1(numCells, (PetscInt**)&s->data);CHKERRQ(ierr); 143 ierr = DMPlexCreateFluent_ReadValues(viewer, s->data, numCells, PETSC_INT, s->index==2012 ? PETSC_TRUE : PETSC_FALSE);CHKERRQ(ierr); 144 ierr = PetscFree(s->data);CHKERRQ(ierr); 145 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, ')');CHKERRQ(ierr); 146 } 147 } 148 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, ')');CHKERRQ(ierr); 149 150 } else if (s->index == 13 || s->index == 2013) { /* Faces */ 151 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, ')');CHKERRQ(ierr); 152 snum = sscanf(buffer, "(%x", &(s->zoneID)); 153 if (snum != 1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Fluent file"); 154 if (s->zoneID == 0) { /* Header section */ 155 snum = sscanf(buffer, "(%x %x %x %d)", &(s->zoneID), &(s->first), &(s->last), &(s->nd)); 156 if (snum != 4) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Fluent file"); 157 } else { /* Data section */ 158 PetscInt f, numEntries, numFaces, numFaceVert; 159 snum = sscanf(buffer, "(%x %x %x %d %d)", &(s->zoneID), &(s->first), &(s->last), &(s->type), &(s->nd)); 160 if (snum != 5) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Fluent file"); 161 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, '(');CHKERRQ(ierr); 162 switch (s->nd) { 163 case 0: numEntries = PETSC_DETERMINE; break; 164 case 2: numEntries = 2 + 2; break; /* linear */ 165 case 3: numEntries = 2 + 3; break; /* triangular */ 166 case 4: numEntries = 2 + 4; break; /* quadrilateral */ 167 default: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Unknown face type in Fluent file"); 168 } 169 numFaces = s->last-s->first + 1; 170 if (numEntries != PETSC_DETERMINE) { 171 /* Allocate space only if we already know the size of the block */ 172 ierr = PetscMalloc1(numEntries*numFaces, (PetscInt**)&s->data);CHKERRQ(ierr); 173 } 174 for (f = 0; f < numFaces; f++) { 175 if (s->nd == 0) { 176 /* Determine the size of the block for "mixed" facets */ 177 ierr = DMPlexCreateFluent_ReadValues(viewer, &numFaceVert, 1, PETSC_INT, s->index==2013 ? PETSC_TRUE : PETSC_FALSE);CHKERRQ(ierr); 178 if (numEntries == PETSC_DETERMINE) { 179 numEntries = numFaceVert + 2; 180 ierr = PetscMalloc1(numEntries*numFaces, (PetscInt**)&s->data);CHKERRQ(ierr); 181 } else { 182 if (numEntries != numFaceVert + 2) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "No support for mixed faces in Fluent files"); 183 } 184 } 185 ierr = DMPlexCreateFluent_ReadValues(viewer, &(((PetscInt*)s->data)[f*numEntries]), numEntries, PETSC_INT, s->index==2013 ? PETSC_TRUE : PETSC_FALSE);CHKERRQ(ierr); 186 } 187 s->nd = numEntries - 2; 188 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, ')');CHKERRQ(ierr); 189 } 190 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, ')');CHKERRQ(ierr); 191 192 } else { /* Unknown section type */ 193 PetscInt depth = 1; 194 do { 195 /* Match parentheses when parsing unknown sections */ 196 do {ierr = PetscViewerRead(viewer, &(buffer[0]), 1, NULL, PETSC_CHAR);CHKERRQ(ierr);} 197 while (buffer[0] != '(' && buffer[0] != ')'); 198 if (buffer[0] == '(') depth++; 199 if (buffer[0] == ')') depth--; 200 } while (depth > 0); 201 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, '\n');CHKERRQ(ierr); 202 } 203 PetscFunctionReturn(0); 204 } 205 206 /*@C 207 DMPlexCreateFluent - Create a DMPlex mesh from a Fluent mesh file. 208 209 Collective on comm 210 211 Input Parameters: 212 + comm - The MPI communicator 213 . viewer - The Viewer associated with a Fluent mesh file 214 - interpolate - Create faces and edges in the mesh 215 216 Output Parameter: 217 . dm - The DM object representing the mesh 218 219 Note: http://aerojet.engr.ucdavis.edu/fluenthelp/html/ug/node1490.htm 220 221 Level: beginner 222 223 .keywords: mesh, fluent, case 224 .seealso: DMPLEX, DMCreate() 225 @*/ 226 PetscErrorCode DMPlexCreateFluent(MPI_Comm comm, PetscViewer viewer, PetscBool interpolate, DM *dm) 227 { 228 PetscMPIInt rank; 229 PetscInt c, f, v, dim = PETSC_DETERMINE, numCells = 0, numVertices = 0, numCellVertices = PETSC_DETERMINE; 230 PetscInt numFaces = PETSC_DETERMINE, numFaceEntries = PETSC_DETERMINE, numFaceVertices = PETSC_DETERMINE; 231 PetscInt *faces = NULL, *cellVertices = NULL, *faceZoneIDs = NULL; 232 PetscInt d, coordSize; 233 PetscScalar *coords, *coordsIn = NULL; 234 PetscSection coordSection; 235 Vec coordinates; 236 PetscErrorCode ierr; 237 238 PetscFunctionBegin; 239 ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 240 241 if (!rank) { 242 FluentSection s; 243 do { 244 ierr = DMPlexCreateFluent_ReadSection(viewer, &s);CHKERRQ(ierr); 245 if (s.index == 2) { /* Dimension */ 246 dim = s.nd; 247 248 } else if (s.index == 10 || s.index == 2010) { /* Vertices */ 249 if (s.zoneID == 0) numVertices = s.last; 250 else { 251 if (coordsIn) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Currently no support for multiple coordinate sets in Fluent files"); 252 coordsIn = (PetscScalar *) s.data; 253 } 254 255 } else if (s.index == 12 || s.index == 2012) { /* Cells */ 256 if (s.zoneID == 0) numCells = s.last; 257 else { 258 switch (s.nd) { 259 case 0: numCellVertices = PETSC_DETERMINE; break; 260 case 1: numCellVertices = 3; break; /* triangular */ 261 case 2: numCellVertices = 4; break; /* tetrahedral */ 262 case 3: numCellVertices = 4; break; /* quadrilateral */ 263 case 4: numCellVertices = 8; break; /* hexahedral */ 264 case 5: numCellVertices = 5; break; /* pyramid */ 265 case 6: numCellVertices = 6; break; /* wedge */ 266 default: numCellVertices = PETSC_DETERMINE; 267 } 268 } 269 270 } else if (s.index == 13 || s.index == 2013) { /* Facets */ 271 if (s.zoneID == 0) { /* Header section */ 272 numFaces = s.last - s.first + 1; 273 if (s.nd == 0 || s.nd == 5) numFaceVertices = PETSC_DETERMINE; 274 else numFaceVertices = s.nd; 275 } else { /* Data section */ 276 if (numFaceVertices != PETSC_DETERMINE && s.nd != numFaceVertices) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Mixed facets in Fluent files are not supported"); 277 if (numFaces < 0) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "No header section for facets in Fluent file"); 278 if (numFaceVertices == PETSC_DETERMINE) numFaceVertices = s.nd; 279 numFaceEntries = numFaceVertices + 2; 280 if (!faces) {ierr = PetscMalloc1(numFaces*numFaceEntries, &faces);CHKERRQ(ierr);} 281 if (!faceZoneIDs) {ierr = PetscMalloc1(numFaces, &faceZoneIDs);CHKERRQ(ierr);} 282 ierr = PetscMemcpy(&(faces[(s.first-1)*numFaceEntries]), s.data, (s.last-s.first+1)*numFaceEntries*sizeof(PetscInt));CHKERRQ(ierr); 283 /* Record the zoneID for each face set */ 284 for (f = s.first -1; f < s.last; f++) faceZoneIDs[f] = s.zoneID; 285 ierr = PetscFree(s.data);CHKERRQ(ierr); 286 } 287 } 288 } while (s.index >= 0); 289 } 290 ierr = MPI_Bcast(&dim, 1, MPIU_INT, 0, comm);CHKERRQ(ierr); 291 if (dim < 0) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Fluent file does not include dimension"); 292 293 /* Allocate cell-vertex mesh */ 294 ierr = DMCreate(comm, dm);CHKERRQ(ierr); 295 ierr = DMSetType(*dm, DMPLEX);CHKERRQ(ierr); 296 ierr = DMSetDimension(*dm, dim);CHKERRQ(ierr); 297 ierr = DMPlexSetChart(*dm, 0, numCells + numVertices);CHKERRQ(ierr); 298 if (!rank) { 299 if (numCells < 0) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Unknown number of cells in Fluent file"); 300 /* If no cell type was given we assume simplices */ 301 if (numCellVertices == PETSC_DETERMINE) numCellVertices = numFaceVertices + 1; 302 for (c = 0; c < numCells; ++c) {ierr = DMPlexSetConeSize(*dm, c, numCellVertices);CHKERRQ(ierr);} 303 } 304 ierr = DMSetUp(*dm);CHKERRQ(ierr); 305 306 if (!rank && faces) { 307 /* Derive cell-vertex list from face-vertex and face-cell maps */ 308 ierr = PetscMalloc1(numCells*numCellVertices, &cellVertices);CHKERRQ(ierr); 309 for (c = 0; c < numCells*numCellVertices; c++) cellVertices[c] = -1; 310 for (f = 0; f < numFaces; f++) { 311 PetscInt *cell; 312 const PetscInt cl = faces[f*numFaceEntries + numFaceVertices]; 313 const PetscInt cr = faces[f*numFaceEntries + numFaceVertices + 1]; 314 const PetscInt *face = &(faces[f*numFaceEntries]); 315 316 if (cl > 0) { 317 cell = &(cellVertices[(cl-1) * numCellVertices]); 318 for (v = 0; v < numFaceVertices; v++) { 319 PetscBool found = PETSC_FALSE; 320 for (c = 0; c < numCellVertices; c++) { 321 if (cell[c] < 0) break; 322 if (cell[c] == face[v]-1 + numCells) {found = PETSC_TRUE; break;} 323 } 324 if (!found) cell[c] = face[v]-1 + numCells; 325 } 326 } 327 if (cr > 0) { 328 cell = &(cellVertices[(cr-1) * numCellVertices]); 329 for (v = 0; v < numFaceVertices; v++) { 330 PetscBool found = PETSC_FALSE; 331 for (c = 0; c < numCellVertices; c++) { 332 if (cell[c] < 0) break; 333 if (cell[c] == face[v]-1 + numCells) {found = PETSC_TRUE; break;} 334 } 335 if (!found) cell[c] = face[v]-1 + numCells; 336 } 337 } 338 } 339 for (c = 0; c < numCells; c++) { 340 ierr = DMPlexSetCone(*dm, c, &(cellVertices[c*numCellVertices]));CHKERRQ(ierr); 341 } 342 } 343 ierr = DMPlexSymmetrize(*dm);CHKERRQ(ierr); 344 ierr = DMPlexStratify(*dm);CHKERRQ(ierr); 345 if (interpolate) { 346 DM idm = NULL; 347 348 ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 349 ierr = DMDestroy(dm);CHKERRQ(ierr); 350 *dm = idm; 351 } 352 353 if (!rank && faces) { 354 PetscInt fi, joinSize, meetSize, *fverts, cells[2]; 355 const PetscInt *join, *meet; 356 ierr = PetscMalloc1(numFaceVertices, &fverts);CHKERRQ(ierr); 357 /* Mark facets by finding the full join of all adjacent vertices */ 358 for (f = 0; f < numFaces; f++) { 359 const PetscInt cl = faces[f*numFaceEntries + numFaceVertices] - 1; 360 const PetscInt cr = faces[f*numFaceEntries + numFaceVertices + 1] - 1; 361 if (cl > 0 && cr > 0) { 362 /* If we know both adjoining cells we can use a single-level meet */ 363 cells[0] = cl; cells[1] = cr; 364 ierr = DMPlexGetMeet(*dm, 2, cells, &meetSize, &meet);CHKERRQ(ierr); 365 if (meetSize != 1) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Could not determine Plex facet for Fluent face %d", f); 366 ierr = DMSetLabelValue(*dm, "Face Sets", meet[0], faceZoneIDs[f]);CHKERRQ(ierr); 367 ierr = DMPlexRestoreMeet(*dm, numFaceVertices, fverts, &meetSize, &meet);CHKERRQ(ierr); 368 } else { 369 for (fi = 0; fi < numFaceVertices; fi++) fverts[fi] = faces[f*numFaceEntries + fi] + numCells - 1; 370 ierr = DMPlexGetFullJoin(*dm, numFaceVertices, fverts, &joinSize, &join);CHKERRQ(ierr); 371 if (joinSize != 1) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Could not determine Plex facet for Fluent face %d", f); 372 ierr = DMSetLabelValue(*dm, "Face Sets", join[0], faceZoneIDs[f]);CHKERRQ(ierr); 373 ierr = DMPlexRestoreJoin(*dm, numFaceVertices, fverts, &joinSize, &join);CHKERRQ(ierr); 374 } 375 } 376 ierr = PetscFree(fverts);CHKERRQ(ierr); 377 } 378 379 /* Read coordinates */ 380 ierr = DMGetCoordinateSection(*dm, &coordSection);CHKERRQ(ierr); 381 ierr = PetscSectionSetNumFields(coordSection, 1);CHKERRQ(ierr); 382 ierr = PetscSectionSetFieldComponents(coordSection, 0, dim);CHKERRQ(ierr); 383 ierr = PetscSectionSetChart(coordSection, numCells, numCells + numVertices);CHKERRQ(ierr); 384 for (v = numCells; v < numCells+numVertices; ++v) { 385 ierr = PetscSectionSetDof(coordSection, v, dim);CHKERRQ(ierr); 386 ierr = PetscSectionSetFieldDof(coordSection, v, 0, dim);CHKERRQ(ierr); 387 } 388 ierr = PetscSectionSetUp(coordSection);CHKERRQ(ierr); 389 ierr = PetscSectionGetStorageSize(coordSection, &coordSize);CHKERRQ(ierr); 390 ierr = VecCreate(PETSC_COMM_SELF, &coordinates);CHKERRQ(ierr); 391 ierr = PetscObjectSetName((PetscObject) coordinates, "coordinates");CHKERRQ(ierr); 392 ierr = VecSetSizes(coordinates, coordSize, PETSC_DETERMINE);CHKERRQ(ierr); 393 ierr = VecSetType(coordinates, VECSTANDARD);CHKERRQ(ierr); 394 ierr = VecGetArray(coordinates, &coords);CHKERRQ(ierr); 395 if (!rank && coordsIn) { 396 for (v = 0; v < numVertices; ++v) { 397 for (d = 0; d < dim; ++d) { 398 coords[v*dim+d] = coordsIn[v*dim+d]; 399 } 400 } 401 } 402 ierr = VecRestoreArray(coordinates, &coords);CHKERRQ(ierr); 403 ierr = DMSetCoordinatesLocal(*dm, coordinates);CHKERRQ(ierr); 404 ierr = VecDestroy(&coordinates);CHKERRQ(ierr); 405 if (!rank) { 406 if (cellVertices) {ierr = PetscFree(cellVertices);CHKERRQ(ierr);} 407 ierr = PetscFree(faces);CHKERRQ(ierr); 408 ierr = PetscFree(faceZoneIDs);CHKERRQ(ierr); 409 ierr = PetscFree(coordsIn);CHKERRQ(ierr); 410 } 411 PetscFunctionReturn(0); 412 } 413