1 #define PETSCDM_DLL 2 #include <petsc-private/dmpleximpl.h> /*I "petscdmplex.h" I*/ 3 4 #undef __FUNCT__ 5 #define __FUNCT__ "DMPlexCreateFluentFromFile" 6 /*@C 7 DMPlexCreateFluentFromFile - Create a DMPlex mesh from a Fluent mesh file 8 9 + comm - The MPI communicator 10 . filename - Name of the Fluent mesh 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(), DMPlexCreateFluent(), DMPlexCreate() 19 @*/ 20 PetscErrorCode DMPlexCreateFluentFromFile(MPI_Comm comm, const char filename[], PetscBool interpolate, DM *dm) 21 { 22 PetscViewer viewer; 23 PetscErrorCode ierr; 24 25 PetscFunctionBegin; 26 /* Create file viewer and build plex */ 27 ierr = PetscViewerCreate(comm, &viewer);CHKERRQ(ierr); 28 ierr = PetscViewerSetType(viewer, PETSCVIEWERASCII);CHKERRQ(ierr); 29 ierr = PetscViewerFileSetMode(viewer, FILE_MODE_READ);CHKERRQ(ierr); 30 ierr = PetscViewerFileSetName(viewer, filename);CHKERRQ(ierr); 31 ierr = DMPlexCreateFluent(comm, viewer, interpolate, dm);CHKERRQ(ierr); 32 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 33 PetscFunctionReturn(0); 34 } 35 36 #undef __FUNCT__ 37 #define __FUNCT__ "DMPlexCreateFluent_ReadString" 38 PetscErrorCode DMPlexCreateFluent_ReadString(PetscViewer viewer, char *buffer, char delim) 39 { 40 PetscInt i = 0; 41 PetscErrorCode ierr; 42 43 PetscFunctionBegin; 44 do {ierr = PetscViewerRead(viewer, &(buffer[i++]), 1, PETSC_CHAR);CHKERRQ(ierr);} 45 while (buffer[i-1] != '\0' && buffer[i-1] != delim); 46 buffer[i] = '\0'; 47 PetscFunctionReturn(0); 48 } 49 50 #undef __FUNCT__ 51 #define __FUNCT__ "DMPlexCreateFluent_ReadSection" 52 PetscErrorCode DMPlexCreateFluent_ReadSection(PetscViewer viewer, FluentSection *s) 53 { 54 char buffer[PETSC_MAX_PATH_LEN]; 55 int snum; 56 PetscErrorCode ierr; 57 58 PetscFunctionBegin; 59 /* Fast-forward to next section and derive its index */ 60 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, '(');CHKERRQ(ierr); 61 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, ' ');CHKERRQ(ierr); 62 snum = sscanf(buffer, "%d", &(s->index)); 63 /* If we can't match an index return -1 to signal end-of-file */ 64 if (snum < 1) {s->index = -1; PetscFunctionReturn(0);} 65 66 if (s->index == 0) { /* Comment */ 67 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, ')');CHKERRQ(ierr); 68 69 } else if (s->index == 2) { /* Dimension */ 70 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, ')');CHKERRQ(ierr); 71 snum = sscanf(buffer, "%d", &(s->nd)); 72 if (snum != 1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Fluent file"); 73 74 } else if (s->index == 10) { /* Vertices */ 75 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, ')');CHKERRQ(ierr); 76 snum = sscanf(buffer, "(%x %x %x %d %d)", &(s->zoneID), &(s->first), &(s->last), &(s->type), &(s->nd)); 77 if (snum != 5) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Fluent file"); 78 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, ')');CHKERRQ(ierr); 79 80 } else if (s->index == 12) { /* Cells */ 81 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, ')');CHKERRQ(ierr); 82 snum = sscanf(buffer, "(%x", &(s->zoneID)); 83 if (snum != 1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Fluent file"); 84 if (s->zoneID == 0) { /* Header section */ 85 snum = sscanf(buffer, "(%x %x %x %d)", &(s->zoneID), &(s->first), &(s->last), &(s->nd)); 86 if (snum != 4) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Fluent file"); 87 } else { /* Data section */ 88 snum = sscanf(buffer, "(%x %x %x %d %d)", &(s->zoneID), &(s->first), &(s->last), &(s->type), &(s->nd)); 89 if (snum != 5) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Fluent file"); 90 } 91 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, ')');CHKERRQ(ierr); 92 93 } else if (s->index == 13) { /* Faces */ 94 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, ')');CHKERRQ(ierr); 95 snum = sscanf(buffer, "(%x", &(s->zoneID)); 96 if (snum != 1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Fluent file"); 97 if (s->zoneID == 0) { /* Header section */ 98 snum = sscanf(buffer, "(%x %x %x %d)", &(s->zoneID), &(s->first), &(s->last), &(s->nd)); 99 if (snum != 4) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Fluent file"); 100 } else { /* Data section */ 101 PetscInt f, e, numEntries, numFaces; 102 snum = sscanf(buffer, "(%x %x %x %d %d)", &(s->zoneID), &(s->first), &(s->last), &(s->type), &(s->nd)); 103 if (snum != 5) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Fluent file"); 104 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, '(');CHKERRQ(ierr); 105 switch (s->nd) { 106 case 0: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Mixed faces in Fluent files are not supported"); 107 case 2: numEntries = 2 + 2; break; /* linear */ 108 case 3: numEntries = 2 + 3; break; /* triangular */ 109 case 4: numEntries = 2 + 4; break; /* quadrilateral */ 110 default: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Unknown face type in Fluent file"); 111 } 112 numFaces = s->last-s->first + 1; 113 ierr = PetscMalloc1(numEntries*numFaces, (PetscInt**)&s->data);CHKERRQ(ierr); 114 for (f = 0; f < numFaces; f++) { 115 for (e = 0; e < numEntries; e++) { 116 ierr = PetscViewerRead(viewer, buffer, 1, PETSC_STRING);CHKERRQ(ierr); 117 snum = sscanf(buffer, "%x", &(((PetscInt*)s->data)[f*numEntries + e])); 118 if (snum != 1) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "File is not a valid Fluent file"); 119 } 120 } 121 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, ')');CHKERRQ(ierr); 122 } 123 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, ')');CHKERRQ(ierr); 124 125 } else { /* Unknown section type */ 126 PetscInt depth = 1; 127 do { 128 /* Match parentheses when parsing unknown sections */ 129 do {ierr = PetscViewerRead(viewer, &(buffer[0]), 1, PETSC_CHAR);CHKERRQ(ierr);} 130 while (buffer[0] != '(' && buffer[0] != ')'); 131 if (buffer[0] == '(') depth++; 132 if (buffer[0] == ')') depth--; 133 } while (depth > 0); 134 ierr = DMPlexCreateFluent_ReadString(viewer, buffer, '\n');CHKERRQ(ierr); 135 } 136 PetscFunctionReturn(0); 137 } 138 139 #undef __FUNCT__ 140 #define __FUNCT__ "DMPlexCreateFluent" 141 /*@C 142 DMPlexCreateFluent - Create a DMPlex mesh from a Fluent mesh file. 143 144 Collective on comm 145 146 Input Parameters: 147 + comm - The MPI communicator 148 . viewer - The Viewer associated with a Fluent mesh 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://aerojet.engr.ucdavis.edu/fluenthelp/html/ug/node1490.htm 155 156 Level: beginner 157 158 .keywords: mesh, fluent, case 159 .seealso: DMPLEX, DMCreate() 160 @*/ 161 PetscErrorCode DMPlexCreateFluent(MPI_Comm comm, PetscViewer viewer, PetscBool interpolate, DM *dm) 162 { 163 PetscMPIInt rank; 164 PetscInt c, f, v, dim = PETSC_DETERMINE, numCells = 0, numVertices = 0, numCellVertices = PETSC_DETERMINE; 165 PetscInt numFaces = PETSC_DETERMINE, numFaceEntries = PETSC_DETERMINE, numFaceVertices = PETSC_DETERMINE; 166 PetscInt *faces = NULL, *cellVertices; 167 PetscErrorCode ierr; 168 169 PetscFunctionBegin; 170 ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 171 172 if (!rank) { 173 FluentSection s; 174 do { 175 ierr = DMPlexCreateFluent_ReadSection(viewer, &s);CHKERRQ(ierr); 176 if (s.index == 2) { /* Dimension */ 177 dim = s.nd; 178 179 } else if (s.index == 10) { /* Vertices */ 180 if (s.zoneID == 0) numVertices = s.last; 181 182 } else if (s.index == 12) { /* Cells */ 183 if (s.zoneID == 0) numCells = s.last; 184 else { 185 switch (s.nd) { 186 case 0: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Mixed elements in Fluent files are not supported"); 187 case 1: numCellVertices = 3; break; /* triangular */ 188 case 2: numCellVertices = 4; break; /* tetrahedral */ 189 case 3: numCellVertices = 4; break; /* quadrilateral */ 190 case 4: numCellVertices = 8; break; /* hexahedral */ 191 case 5: numCellVertices = 5; break; /* pyramid */ 192 case 6: numCellVertices = 6; break; /* wedge */ 193 default: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Unknown cell element-type in Fluent file"); 194 } 195 } 196 197 } else if (s.index == 13) { /* Facets */ 198 if (s.zoneID == 0) { /* Header section */ 199 numFaces = s.last - s.first + 1; 200 } else { /* Data section */ 201 if (s.nd == 0 || (numFaceVertices != PETSC_DETERMINE && s.nd != numFaceVertices)) { 202 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Mixed facets in Fluent files are currently not supported"); 203 } 204 if (numFaces == PETSC_DETERMINE) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "No header section for facets in Fluent file"); 205 if (numFaceVertices == PETSC_DETERMINE) numFaceVertices = s.nd; 206 numFaceEntries = numFaceVertices + 2; 207 if (!faces) {ierr = PetscMalloc1(numFaces*numFaceEntries, &faces);CHKERRQ(ierr);} 208 ierr = PetscMemcpy(&(faces[(s.first-1)*numFaceEntries]), s.data, (s.last-s.first+1)*numFaceEntries*sizeof(PetscInt));CHKERRQ(ierr); 209 ierr = PetscFree(s.data);CHKERRQ(ierr); 210 } 211 } 212 } while (s.index >= 0); 213 } 214 ierr = MPI_Bcast(&dim, 1, MPIU_INT, 0, comm);CHKERRQ(ierr); 215 if (dim < 0) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Fluent file does not include dimension"); 216 217 /* Allocate cell-vertex mesh */ 218 ierr = DMCreate(comm, dm);CHKERRQ(ierr); 219 ierr = DMSetType(*dm, DMPLEX);CHKERRQ(ierr); 220 ierr = DMSetDimension(*dm, dim);CHKERRQ(ierr); 221 ierr = DMPlexSetChart(*dm, 0, numCells + numVertices);CHKERRQ(ierr); 222 if (!rank) { 223 if (numCellVertices == PETSC_DETERMINE) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Cell type not specified in Fluent file"); 224 for (c = 0; c < numCells; ++c) {ierr = DMPlexSetConeSize(*dm, c, numCellVertices);CHKERRQ(ierr);} 225 } 226 ierr = DMSetUp(*dm);CHKERRQ(ierr); 227 228 if (!rank) { 229 /* Derive cell-vertex list from face-vertex and face-cell maps */ 230 if (numCells == PETSC_DETERMINE || numCellVertices == PETSC_DETERMINE) { 231 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Insufficent cell header information in Fluent file"); 232 } 233 ierr = PetscMalloc1(numCells*numCellVertices, &cellVertices);CHKERRQ(ierr); 234 for (c = 0; c < numCells*numCellVertices; c++) cellVertices[c] = -1; 235 for (f = 0; f < numFaces; f++) { 236 PetscInt *cell; 237 const PetscInt cl = faces[f*numFaceEntries + numFaceVertices]; 238 const PetscInt cr = faces[f*numFaceEntries + numFaceVertices + 1]; 239 const PetscInt *face = &(faces[f*numFaceEntries]); 240 241 if (cl > 0) { 242 cell = &(cellVertices[(cl-1) * numCellVertices]); 243 for (v = 0; v < numFaceVertices; v++) { 244 PetscBool found = PETSC_FALSE; 245 for (c = 0; c < numCellVertices; c++) { 246 if (cell[c] < 0) break; 247 if (cell[c] == face[v]-1 + numCells) {found = PETSC_TRUE; break;} 248 } 249 if (!found) cell[c] = face[v]-1 + numCells; 250 } 251 } 252 if (cr > 0) { 253 cell = &(cellVertices[(cr-1) * numCellVertices]); 254 for (v = 0; v < numFaceVertices; v++) { 255 PetscBool found = PETSC_FALSE; 256 for (c = 0; c < numCellVertices; c++) { 257 if (cell[c] < 0) break; 258 if (cell[c] == face[v]-1 + numCells) {found = PETSC_TRUE; break;} 259 } 260 if (!found) cell[c] = face[v]-1 + numCells; 261 } 262 } 263 } 264 for (c = 0; c < numCells; c++) { 265 ierr = DMPlexSetCone(*dm, c, &(cellVertices[c*numCellVertices]));CHKERRQ(ierr); 266 } 267 } 268 ierr = DMPlexSymmetrize(*dm);CHKERRQ(ierr); 269 ierr = DMPlexStratify(*dm);CHKERRQ(ierr); 270 if (interpolate) { 271 DM idm = NULL; 272 273 ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 274 ierr = DMDestroy(dm);CHKERRQ(ierr); 275 *dm = idm; 276 } 277 278 if (!rank) { 279 ierr = PetscFree(cellVertices);CHKERRQ(ierr); 280 ierr = PetscFree(faces);CHKERRQ(ierr); 281 } 282 PetscFunctionReturn(0); 283 } 284