STARCDMeshWriter.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2020 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "STARCDMeshWriter.H"
27 #include "Time.H"
28 #include "OFstream.H"
29 #include "OSspecific.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 const char* Foam::meshWriters::STARCD::defaultBoundaryName =
34  "Default_Boundary_Region";
35 
37 {
38  { 4, 5, 2, 3, 0, 1 }, // 11 = pro-STAR hex
39  { 0, 1, 4, 5, 2, -1 }, // 12 = pro-STAR prism
40  { 5, 4, 2, 0, -1, -1 }, // 13 = pro-STAR tetra
41  { 0, 4, 3, 5, 2, -1 } // 14 = pro-STAR pyramid
42 };
43 
44 
45 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
46 
47 Foam::label Foam::meshWriters::STARCD::findDefaultBoundary() const
48 {
49  const polyBoundaryMesh& patches = mesh_.boundaryMesh();
50 
51  label id = -1;
52 
53  // find Default_Boundary_Region if it exists
54  forAll(patches, patchi)
55  {
56  if (defaultBoundaryName == patches[patchi].name())
57  {
58  id = patchi;
59  break;
60  }
61  }
62  return id;
63 }
64 
65 
66 void Foam::meshWriters::STARCD::getCellTable()
67 {
68  // read constant/polyMesh/propertyName
69  IOList<label> ioList
70  (
71  IOobject
72  (
73  "cellTableId",
74  mesh_.time().constant(),
76  mesh_,
79  false
80  )
81  );
82 
83  bool useCellZones = false;
85 
86  // get information from constant/polyMesh/cellTableId if possible
87  if (ioList.headerOk())
88  {
89  if (ioList.size() == mesh_.nCells())
90  {
91  cellTableId_.transfer(ioList);
92 
93  if (cellTable_.empty())
94  {
95  Info<< "no cellTable information available" << endl;
96  }
97  }
98  else
99  {
101  << ioList.objectPath() << " has incorrect number of cells "
102  << " - use cellZone information"
103  << endl;
104 
105  ioList.clear();
106  useCellZones = true;
107  }
108  }
109  else
110  {
111  useCellZones = true;
112  }
113 
114 
115  if (useCellZones)
116  {
117  if (cellTable_.empty())
118  {
119  Info<< "created cellTable from cellZones" << endl;
120  cellTable_ = mesh_;
121  }
122 
123  // track if there are unzoned cells
124  label nUnzoned = mesh_.nCells();
125 
126  // get the cellZone <-> cellTable correspondence
127  Info<< "matching cellZones to cellTable" << endl;
128 
129  forAll(mesh_.cellZones(), zoneI)
130  {
131  const cellZone& cZone = mesh_.cellZones()[zoneI];
132  if (cZone.size())
133  {
134  nUnzoned -= cZone.size();
135 
136  label tableId = cellTable_.findIndex(cZone.name());
137  if (tableId < 0)
138  {
139  dictionary dict;
140 
141  dict.add("Label", cZone.name());
142  dict.add("MaterialType", "fluid");
143  tableId = cellTable_.append(dict);
144  }
145 
146  forAll(cZone, i)
147  {
148  cellTableId_[cZone[i]] = tableId;
149  }
150  }
151  }
152 
153  if (nUnzoned)
154  {
155  dictionary dict;
156 
157  dict.add("Label", "__unZonedCells__");
158  dict.add("MaterialType", "fluid");
159  label tableId = cellTable_.append(dict);
160 
161  forAll(cellTableId_, i)
162  {
163  if (cellTableId_[i] < 0)
164  {
165  cellTableId_[i] = tableId;
166  }
167  }
168  }
169  }
170 }
171 
172 
173 void Foam::meshWriters::STARCD::writeHeader(Ostream& os, const char* filetype)
174 {
175  os << "PROSTAR_" << filetype << nl
176  << 4000
177  << " " << 0
178  << " " << 0
179  << " " << 0
180  << " " << 0
181  << " " << 0
182  << " " << 0
183  << " " << 0
184  << endl;
185 }
186 
187 
188 void Foam::meshWriters::STARCD::writePoints(const fileName& prefix) const
189 {
190  OFstream os(prefix + ".vrt");
191  writeHeader(os, "VERTEX");
192 
193  // Set the precision of the points data to 10
194  os.precision(10);
195 
196  // force decimal point for Fortran input
197  os.setf(std::ios::showpoint);
198 
199  const pointField& points = mesh_.points();
200 
201  Info<< "Writing " << os.name() << " : "
202  << points.size() << " points" << endl;
203 
204  forAll(points, ptI)
205  {
206  // convert [m] -> [mm]
207  os
208  << ptI + 1 << " "
209  << scaleFactor_ * points[ptI].x() << " "
210  << scaleFactor_ * points[ptI].y() << " "
211  << scaleFactor_ * points[ptI].z() << nl;
212  }
213  os.flush();
214 
215 }
216 
217 
218 void Foam::meshWriters::STARCD::writeCells(const fileName& prefix) const
219 {
220  OFstream os(prefix + ".cel");
221  writeHeader(os, "CELL");
222 
223  // this is what we seem to need
224  // map foam cellModeller index -> star shape
225  Map<label> shapeLookupIndex;
226  shapeLookupIndex.insert(hexModel->index(), 11);
227  shapeLookupIndex.insert(prismModel->index(), 12);
228  shapeLookupIndex.insert(tetModel->index(), 13);
229  shapeLookupIndex.insert(pyrModel->index(), 14);
230 
231  const cellShapeList& shapes = mesh_.cellShapes();
232  const cellList& cells = mesh_.cells();
233  const faceList& faces = mesh_.faces();
234  const labelList& owner = mesh_.faceOwner();
235 
236  Info<< "Writing " << os.name() << " : "
237  << cells.size() << " cells" << endl;
238 
239  forAll(cells, cellId)
240  {
241  label tableId = cellTableId_[cellId];
242  label materialType = 1; // 1(fluid)
243  if (cellTable_.found(tableId))
244  {
245  const dictionary& dict = cellTable_[tableId];
246  if (dict.found("MaterialType"))
247  {
248  word matType;
249  dict.lookup("MaterialType") >> matType;
250  if (matType == "solid")
251  {
252  materialType = 2;
253  }
254 
255  }
256  }
257 
258  const cellShape& shape = shapes[cellId];
259  label mapIndex = shape.model().index();
260 
261  // a registered primitive type
262  if (shapeLookupIndex.found(mapIndex))
263  {
264  label shapeId = shapeLookupIndex[mapIndex];
265  const labelList& vrtList = shapes[cellId];
266 
267  os << cellId + 1
268  << " " << shapeId
269  << " " << vrtList.size()
270  << " " << tableId
271  << " " << materialType;
272 
273  // primitives have <= 8 vertices, but prevent overrun anyhow
274  // indent following lines for ease of reading
275  label count = 0;
276  forAll(vrtList, i)
277  {
278  if ((count % 8) == 0)
279  {
280  os << nl
281  << " " << cellId + 1;
282  }
283  os << " " << vrtList[i] + 1;
284  count++;
285  }
286  os << endl;
287 
288  }
289  else
290  {
291  label shapeId = 255; // treat as general polyhedral
292  const labelList& cFaces = cells[cellId];
293 
294  // create (beg,end) indices
295  List<label> indices(cFaces.size() + 1);
296  indices[0] = indices.size();
297 
298  label count = indices.size();
299  // determine the total number of vertices
300  forAll(cFaces, facei)
301  {
302  count += faces[cFaces[facei]].size();
303  indices[facei+1] = count;
304  }
305 
306  os << cellId + 1
307  << " " << shapeId
308  << " " << count
309  << " " << tableId
310  << " " << materialType;
311 
312  // write indices - max 8 per line
313  // indent following lines for ease of reading
314  count = 0;
315  forAll(indices, i)
316  {
317  if ((count % 8) == 0)
318  {
319  os << nl
320  << " " << cellId + 1;
321  }
322  os << " " << indices[i];
323  count++;
324  }
325 
326  // write faces - max 8 per line
327  forAll(cFaces, facei)
328  {
329  label meshFace = cFaces[facei];
330  face f;
331 
332  if (owner[meshFace] == cellId)
333  {
334  f = faces[meshFace];
335  }
336  else
337  {
338  f = faces[meshFace].reverseFace();
339  }
340 
341  forAll(f, i)
342  {
343  if ((count % 8) == 0)
344  {
345  os << nl
346  << " " << cellId + 1;
347  }
348 
349  os << " " << f[i] + 1;
350  count++;
351  }
352  }
353 
354  os << endl;
355  }
356  }
357 }
358 
359 
360 void Foam::meshWriters::STARCD::writeBoundary(const fileName& prefix) const
361 {
362  OFstream os(prefix + ".bnd");
363  writeHeader(os, "BOUNDARY");
364 
365  const cellShapeList& shapes = mesh_.cellShapes();
366  const cellList& cells = mesh_.cells();
367  const faceList& faces = mesh_.faces();
368  const labelList& owner = mesh_.faceOwner();
369  const polyBoundaryMesh& patches = mesh_.boundaryMesh();
370 
371  // this is what we seem to need
372  // these MUST correspond to foamToStarFaceAddr
373  //
374  Map<label> faceLookupIndex;
375  faceLookupIndex.insert(hexModel->index(), 0);
376  faceLookupIndex.insert(prismModel->index(), 1);
377  faceLookupIndex.insert(tetModel->index(), 2);
378  faceLookupIndex.insert(pyrModel->index(), 3);
379 
380  Info<< "Writing " << os.name() << " : "
381  << (mesh_.nFaces() - patches[0].start()) << " boundaries" << endl;
382 
383 
384  label defaultId = findDefaultBoundary();
385 
386  //
387  // write boundary faces - skip Default_Boundary_Region entirely
388  //
389  label boundId = 0;
390  forAll(patches, patchi)
391  {
392  label regionId = patchi;
393  if (regionId == defaultId)
394  {
395  continue; // skip - already written
396  }
397  else if (defaultId == -1 || regionId < defaultId)
398  {
399  regionId++;
400  }
401 
402  label patchStart = patches[patchi].start();
403  label patchSize = patches[patchi].size();
404  word bndType = boundaryRegion_.boundaryType(patches[patchi].name());
405 
406  for
407  (
408  label facei = patchStart;
409  facei < (patchStart + patchSize);
410  ++facei
411  )
412  {
413  label cellId = owner[facei];
414  const labelList& cFaces = cells[cellId];
415  const cellShape& shape = shapes[cellId];
416  label cellFaceId = findIndex(cFaces, facei);
417 
418  // Info<< "cell " << cellId + 1 << " face " << facei
419  // << " == " << faces[facei]
420  // << " is index " << cellFaceId << " from " << cFaces;
421 
422  // Unfortunately, the order of faces returned by
423  // primitiveMesh::cells() is not necessarily the same
424  // as defined by primitiveMesh::cellShapes()
425  // Thus, for registered primitive types, do the lookup ourselves.
426  // Finally, the cellModel face number is re-mapped to the
427  // STAR-CD local face number
428 
429  label mapIndex = shape.model().index();
430 
431  // a registered primitive type
432  if (faceLookupIndex.found(mapIndex))
433  {
434  const faceList sFaces = shape.faces();
435  forAll(sFaces, sFacei)
436  {
437  if (faces[facei] == sFaces[sFacei])
438  {
439  cellFaceId = sFacei;
440  break;
441  }
442  }
443 
444  mapIndex = faceLookupIndex[mapIndex];
445  cellFaceId = foamToStarFaceAddr[mapIndex][cellFaceId];
446  }
447  // Info<< endl;
448 
449  boundId++;
450 
451  os
452  << boundId
453  << " " << cellId + 1
454  << " " << cellFaceId + 1
455  << " " << regionId
456  << " " << 0
457  << " " << bndType.c_str()
458  << endl;
459  }
460  }
461 }
462 
463 
464 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
465 
467 (
468  const polyMesh& mesh,
469  const scalar scaleFactor
470 )
471 :
472  meshWriter(mesh, scaleFactor)
473 {
476  getCellTable();
477 }
478 
479 
480 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
481 
483 {}
484 
485 
486 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
487 
488 void Foam::meshWriters::STARCD::rmFiles(const fileName& baseName) const
489 {
490  rm(baseName + ".vrt");
491  rm(baseName + ".cel");
492  rm(baseName + ".bnd");
493  rm(baseName + ".inp");
494 }
495 
496 
498 {
499  fileName baseName(meshName);
500 
501  if (baseName.empty())
502  {
503  baseName = meshWriter::defaultMeshName;
504 
505  if
506  (
507  mesh_.time().timeName() != "0"
508  && mesh_.time().timeName() != mesh_.time().constant()
509  )
510  {
511  baseName += "_" + mesh_.time().timeName();
512  }
513  }
514 
515  rmFiles(baseName);
516  writePoints(baseName);
517  writeCells(baseName);
518 
519  if (writeBoundary_)
520  {
521  writeBoundary(baseName);
522  }
523 
524  return true;
525 }
526 
527 
528 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:434
label append(const dictionary &)
Append to the end, return index.
Definition: cellTable.C:130
static string defaultMeshName
Specify a default mesh name.
Definition: meshWriter.H:126
dictionary dict
static const label foamToStarFaceAddr[4][6]
Face addressing from OpenFOAM faces -> pro-STAR faces.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
A class for handling file names.
Definition: fileName.H:79
virtual bool write(const fileName &meshName=fileName::null) const
Write volume mesh.
const meshCellZones & cellZones() const
Return cell zones.
Definition: polyMesh.H:482
static const cellModel * tetModel
Definition: meshWriter.H:115
static const cellModel * prismModel
Definition: meshWriter.H:117
label nFaces() const
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
Definition: polyMesh.H:312
bool writeBoundary_
Write bnd file.
Definition: meshWriter.H:102
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
const cellShapeList & cellShapes() const
Return cell shapes.
boundaryRegion boundaryRegion_
boundaryRegion persistent data saved as a dictionary
Definition: meshWriter.H:105
bool empty() const
Return true if the hash table is empty.
Definition: HashTableI.H:72
List< face > faceList
Definition: faceListFwd.H:43
label nCells() const
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
label count(const ListType &l, typename ListType::const_reference x)
Count the number of occurrences of a value in a list.
const cellList & cells() const
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:636
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1131
STARCD(const polyMesh &, const scalar scaleFactor=1.0)
Open a file for writing.
Functions used by OpenFOAM that are specific to POSIX compliant operating systems and need to be repl...
void readDict(const objectRegistry &, const word &name="cellTable", const fileName &instance="constant")
Read constant/cellTable.
Definition: cellTable.C:313
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
List< cellShape > cellShapeList
List of cellShapes and PtrList of List of cellShape.
Definition: cellShapeList.H:43
static char meshName[]
Definition: globalFoam.H:7
static const cellModel * pyrModel
Definition: meshWriter.H:116
bool found(const Key &) const
Return true if hashedEntry is found in table.
Definition: HashTable.C:113
const word & constant() const
Return constant name.
Definition: TimePaths.H:123
virtual const labelList & faceOwner() const
Return face owner.
Definition: polyMesh.C:1169
List< label > labelList
A List of labels.
Definition: labelList.H:56
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1156
scalar scaleFactor_
Scaling factor for points (eg, [m] -> [mm])
Definition: meshWriter.H:99
labelList cellTableId_
cellTable IDs for each cell
Definition: meshWriter.H:111
static const char nl
Definition: Ostream.H:260
word boundaryType(const word &name) const
Return BoundaryType corresponding to patch &#39;name&#39;.
const Time & time() const
Return time.
labelList f(nPoints)
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurrence of given element and return index,.
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
label index() const
Return index of model in the model list.
Definition: cellModelI.H:44
void readDict(const objectRegistry &, const word &name="boundaryRegion", const fileName &instance="constant")
Read constant/boundaryRegion.
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
void rmFiles(const fileName &baseName) const
Remove STAR-CD files for the baseName.
void setSize(const label)
Reset size of List.
Definition: List.C:281
label patchi
#define WarningInFunction
Report a warning using Foam::Warning.
virtual ~STARCD()
Destructor.
label cellId
messageStream Info
static const cellModel * hexModel
Definition: meshWriter.H:118
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
bool rm(const fileName &)
Remove a file, returning true if successful otherwise false.
Definition: POSIX.C:1017
meshWriter(const polyMesh &, const scalar scaleFactor=1.0)
Create a writer object.
Definition: meshWriter.C:71
label findIndex(const word &name) const
Return index corresponding to name.
Definition: cellTable.C:206
cellTable cellTable_
cellTable persistent data saved as a dictionary
Definition: meshWriter.H:108
void transfer(List< T > &)
Transfer the contents of the argument List into this list.
Definition: List.C:342
List< cell > cellList
list of cells
Definition: cellList.H:42
const polyMesh & mesh_
Mesh reference.
Definition: meshWriter.H:96