57 rotatedFace.setSize(
f.size());
60 rotatedFace[i] =
f[fp];
71 const Map<label>& nodeIdToPoints,
72 DynamicList<label>& verts
76 for (label i = 0; i < nVerts; i++)
86 verts.push_back(verti-1);
98 Map<label>& elemToFoam
101 const label sz = foamToElem.size();
102 foamToElem.resize(sz+nShapes);
105 elemToFoam.reserve(elemToFoam.size()+nShapes);
106 for (label shapei = 0; shapei < nShapes; shapei++)
110 foamToElem[sz+shapei] = elemi;
111 elemToFoam.insert(elemi, sz+shapei);
116 for (label shapei = 0; shapei < nShapes; shapei++)
118 foamToElem[sz+shapei] = sz+shapei;
126 const cellModel& model,
127 DynamicList<label>& verts,
153 if (model.mag(verts,
points) < 0)
155 if (verts.size() == 8)
158 std::swap(verts[0], verts[4]);
159 std::swap(verts[1], verts[5]);
160 std::swap(verts[2], verts[6]);
161 std::swap(verts[3], verts[7]);
163 else if (verts.size() == 4)
166 std::swap(verts[0], verts[1]);
168 else if (verts.size() == 5)
171 std::swap(verts[1], verts[3]);
173 else if (verts.size() == 6)
176 std::swap(verts[0], verts[3]);
177 std::swap(verts[1], verts[4]);
178 std::swap(verts[2], verts[5]);
187 const bool read_node_ids,
188 const bool read_elem_ids,
192 Map<label>& nodeIdToPoints,
197 Map<label>& elemIdToCells,
202 Map<label>& elemIdToFaces
210 DynamicList<label> verts;
217 is.readKeyword(line);
219 while (line.empty() && is.good());
223 if (
split.size() == 0)
228 if (
split[0] ==
"part")
232 else if (
split[0] ==
"node_ids")
236 for (label i = 0; i <
nPoints; i++)
242 else if (
split[0] ==
"coordinates")
248 <<
" starting at line " << is.lineNumber()
249 <<
" position " << is.stdStream().tellg() <<
endl;
262 for (label pointi = 0; pointi <
nPoints; pointi++)
266 for (label pointi = 0; pointi <
nPoints; pointi++)
270 for (label pointi = 0; pointi <
nPoints; pointi++)
272 is.read(
points[pointi].z());
275 else if (
split[0] ==
"tetra4")
282 <<
" position " << is.stdStream().tellg() <<
endl;
284 const label celli =
cells.size();
285 cells.resize(celli+nShapes);
297 for (label shapei = 0; shapei < nShapes; shapei++)
299 readVerts(is, 4, nodeIdToPoints, verts);
302 setHandedness(model, verts,
points);
304 const cellShape cellVerts(model, verts);
305 cells[celli+shapei] = cellVerts.faces();
308 else if (
split[0] ==
"pyramid5")
315 <<
" position " << is.stdStream().tellg() <<
endl;
317 const label celli =
cells.size();
318 cells.resize(celli+nShapes);
330 for (label shapei = 0; shapei < nShapes; shapei++)
332 readVerts(is, 5, nodeIdToPoints, verts);
335 setHandedness(model, verts,
points);
337 const cellShape cellVerts(model, verts);
338 cells[celli+shapei] = cellVerts.faces();
341 else if (
split[0] ==
"penta6")
348 <<
" position " << is.stdStream().tellg() <<
endl;
350 const label celli =
cells.size();
351 cells.resize(celli+nShapes);
363 for (label shapei = 0; shapei < nShapes; shapei++)
365 readVerts(is, 6, nodeIdToPoints, verts);
368 setHandedness(model, verts,
points);
370 const cellShape cellVerts(model, verts);
371 cells[celli+shapei] = cellVerts.faces();
374 else if (
split[0] ==
"hexa8")
381 <<
" position " << is.stdStream().tellg() <<
endl;
383 const label celli =
cells.size();
384 cells.resize(celli+nShapes);
396 for (label shapei = 0; shapei < nShapes; shapei++)
398 readVerts(is, 8, nodeIdToPoints, verts);
401 setHandedness(model, verts,
points);
403 const cellShape cellVerts(model, verts);
404 cells[celli+shapei] = cellVerts.faces();
407 else if (
split[0] ==
"nfaced")
414 <<
" position " << is.stdStream().tellg() <<
endl;
416 const label celli =
cells.size();
417 cells.resize(celli+nShapes);
428 for (label shapei = 0; shapei < nShapes; shapei++)
433 cellFaces.setSize(nFaces);
436 for (label shapei = 0; shapei < nShapes; shapei++)
439 forAll(cellFaces, cellFacei)
443 cellFaces[cellFacei].setSize(nVerts);
447 for (label shapei = 0; shapei < nShapes; shapei++)
450 forAll(cellFaces, cellFacei)
452 face&
f = cellFaces[cellFacei];
453 readVerts(is,
f.size(), nodeIdToPoints, verts);
454 f.labelList::operator=(verts);
458 if (
f[fp] < 0 ||
f[fp] >=
points.size())
462 <<
" indexes outside points:" <<
points.size()
469 else if (
split[0] ==
"tria3")
476 <<
" position " << is.stdStream().tellg() <<
endl;
478 const label facei = faces.size();
489 faces.setSize(facei+nShapes);
491 for (label shapei = 0; shapei < nShapes; shapei++)
493 auto&
f = faces[facei+shapei];
495 readVerts(is,
f.size(), nodeIdToPoints, verts);
496 f.labelList::operator=(verts);
499 else if (
split[0] ==
"quad4")
506 <<
" position " << is.stdStream().tellg() <<
endl;
508 const label facei = faces.size();
519 faces.setSize(facei+nShapes);
521 for (label shapei = 0; shapei < nShapes; shapei++)
523 auto&
f = faces[facei+shapei];
525 readVerts(is,
f.size(), nodeIdToPoints, verts);
526 f.labelList::operator=(verts);
529 else if (
split[0] ==
"nsided")
536 <<
" position " << is.stdStream().tellg() <<
endl;
538 const label facei = faces.size();
549 faces.setSize(facei+nShapes);
551 for (label shapei = 0; shapei < nShapes; shapei++)
553 auto&
f = faces[facei+shapei];
559 for (label shapei = 0; shapei < nShapes; shapei++)
561 auto&
f = faces[facei+shapei];
562 readVerts(is,
f.size(), nodeIdToPoints, verts);
563 f.labelList::operator=(verts);
569 <<
" from line " << line
570 <<
" starting at line " << is.lineNumber()
571 <<
" position " << is.stdStream().tellg() <<
endl;
584 const scalar scaleFactor
587 ensightReadFile is(geometryFile_);
590 is.readBinaryHeader();
594 Info<<
"Ensight : " << header <<
endl;
596 Info<<
"Ensight : " << header <<
endl;
599 bool read_node_ids =
false;
600 bool read_elem_ids =
false;
605 List<string> partNames;
607 DynamicList<label> partIDs;
608 PtrList<pointField> partPoints;
609 PtrList<labelList> partNodeIDs;
610 PtrList<Map<label>> partPointIDs;
614 PtrList<faceListList> partCells;
616 PtrList<labelList> partCellIDs;
617 PtrList<Map<label>> partCellElemIDs;
619 PtrList<faceList> partFaces;
621 PtrList<labelList> partFaceIDs;
622 PtrList<Map<label>> partFaceElemIDs;
631 is.readKeyword(line);
633 while (line.empty() && is.good());
636 if (
split[0] ==
"extents")
647 <<
"Read extents " << boundBox(
min,
max)
650 else if (
split[0] ==
"node")
654 if (op ==
"given" || op ==
"ignore")
657 read_node_ids =
true;
660 else if (
split[0] ==
"element")
664 if (op ==
"given" || op ==
"ignore")
667 read_elem_ids =
true;
670 else if (
split[0] ==
"part")
672 bool finished =
false;
676 partIDs.emplace_back();
677 is.read(partIDs.back());
679 partNames.emplace_back();
680 is.read(partNames.back());
683 <<
"Reading part " << partIDs.back()
684 <<
" name " << partNames.back()
686 <<
" position " << is.stdStream().tellg() <<
endl;
689 finished = readGoldPart
695 partPoints.emplace_back(),
696 partNodeIDs.emplace_back(),
697 partPointIDs.emplace_back(),
700 partCells.emplace_back(),
701 partCellIDs.emplace_back(),
702 partCellElemIDs.emplace_back(),
705 partFaces.emplace_back(),
706 partFaceIDs.emplace_back(),
707 partFaceElemIDs.emplace_back()
710 partPoints.back() *= scaleFactor;
713 <<
"For part " << partIDs.back()
714 <<
" read cells " << partCells.back().size()
715 <<
" faces " << partFaces.back().size()
791 extents.add(partPoints[parti]);
792 nPoints += partPoints[parti].size();
794 const scalar mergeDist = mergeTol_*extents.mag();
796 Pout<<
"Merging points closer than " << mergeDist
797 <<
" calculated from bounding box " << extents
798 <<
" and mergeTol " << mergeTol_
803 const auto& pPoints = partPoints[parti];
804 auto& pPointMap = pointToMerged[parti];
805 pPointMap.setSize(pPoints.size(), -1);
807 Pout<<
"Matching part " << parti
808 <<
" name " << partNames[parti]
809 <<
" points " << pPoints.size()
810 <<
" to current merged points " << points_.size()
816 pPointMap =
identity(pPoints.size());
832 forAll(from0To1, partPointi)
834 const label pointi = from0To1[partPointi];
837 pPointMap[partPointi] = pointi;
845 const label nOldPoints = points_.size();
846 points_.setSize(nOldPoints+nAdded);
848 forAll(from0To1, partPointi)
850 if (from0To1[partPointi] == -1)
852 const label
newPointi = nOldPoints+nAdded++;
853 points_[
newPointi] = pPoints[partPointi];
861 nodeIds_.setSize(points_.size());
862 forAll(partNodeIDs, parti)
864 const auto& pPointMap = pointToMerged[parti];
865 UIndirectList<label>(nodeIds_, pPointMap) = partNodeIDs[parti];
868 Pout<<
"Merged from " <<
nPoints <<
" down to " << points_.size()
869 <<
" points" <<
endl;
874 labelList cellOffsets(partCells.size()+1);
880 nCells += partCells[parti].size();
881 cellOffsets[parti+1] = nCells;
884 cellFaces_.setSize(nCells);
885 cellTableId_.setSize(nCells);
886 elementIds_.setSize(nCells, -1);
891 const auto&
cells = partCells[parti];
893 SubList<faceList> cellSlice
913 ) = partCellIDs[parti];
916 const auto& pointMap = pointToMerged[parti];
918 for (
auto& cellFaces : cellSlice)
920 for (
auto&
f : cellFaces)
932 cellTable_.setMaterial(parti,
"fluid");
933 cellTable_.setName(parti,
"part" +
Foam::name(partIDs[parti]));
940 HashTable<cellFaceIdentifier, face, face::symmHasher> vertsToCell
951 const auto& cFaces = cellFaces_[celli];
954 const face&
f = rotateFace(cFaces[cFacei], rotatedFace);
956 const auto fFnd = vertsToCell.find(
f);
960 vertsToCell.erase(fFnd);
964 vertsToCell.insert(
f, cellFaceIdentifier(celli, cFacei));
974 if (partFaces[parti].size())
976 Pout<<
"Using part " << parti
977 <<
" name " << partNames[parti]
986 patchTypes_.setSize(
nPatches,
"wall");
988 forAll(patchNames_, patchi)
990 const label parti = patchToPart[patchi];
992 Pout<<
"Matching part " << parti
993 <<
" name " << partNames[parti]
994 <<
" faces " << partFaces[parti].size()
996 <<
" to merged faces " << vertsToCell.size()
997 <<
", merged points " << points_.size()
1002 const auto& pointMap = pointToMerged[parti];
1003 const auto& faces = partFaces[parti];
1005 auto& partCellAndFace = boundaryIds_[patchi];
1006 partCellAndFace.setSize(faces.size());
1008 label patchFacei = 0;
1011 if (faces[facei].empty())
1013 Pout<<
"Ignoring empty face:" << facei <<
endl;
1018 const face newF(pointMap, faces[facei]);
1020 const auto cAndF = vertsToCell.find
1041 partCellAndFace[patchFacei++] = cAndF();
1042 vertsToCell.erase(cAndF);
1045 partCellAndFace.setSize(patchFacei);
1047 patchPhysicalTypes_.setSize(
nPatches,
"wall");
1051 if (vertsToCell.size())
1054 boundaryIds_.emplace_back(vertsToCell.size());
1056 auto& cellAndFaces = boundaryIds_.back();
1060 cellAndFaces[i++] = iter.val();
1064 patchTypes_.push_back(
"empty");
1065 patchNames_.push_back(
"defaultFaces");
1066 patchPhysicalTypes_.push_back(
"empty");
1067 patchStarts_.push_back(0);
1068 patchSizes_.push_back(0);
1070 Pout<<
"Introducing default patch " << patchNames_.size()-1
1071 <<
" name " << patchNames_.back() <<
endl;
1082 const fileName& geomFile,
1083 const objectRegistry& registry,
1084 const scalar mergeTol,
1085 const scalar scaleFactor,
1086 const bool setHandedness
1089 meshReader(geomFile, scaleFactor),
1090 mergeTol_(mergeTol),
1091 setHandedness_(setHandedness)
static word validate(const std::string &s, const bool prefix=false)
Construct validated word (no invalid characters).
List< faceList > faceListList
List of faceList.
Ostream & indent(Ostream &os)
Indent stream.
errorManipArg< error, int > exit(error &err, const int errNo=1)
A face is a list of labels corresponding to mesh vertices.
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Ostream & endl(Ostream &os)
Add newline and flush stream.
List< labelList > labelListList
List of labelList.
Determine correspondence between points. See below.
#define forAll(list, i)
Loop across all elements in list.
static const cellModel & ref(const modelType model)
Look up reference to cellModel by enumeration. Fatal on failure.
List< face > faceList
List of faces.
label findMin(const ListType &input, label start=0)
Linear search for the index of the min element, similar to std::min_element but for lists and returns...
bool matchPoints(const UList< point > &pts0, const UList< point > &pts1, const UList< scalar > &matchDistance, const bool verbose, labelList &from0To1, const point &origin=point::zero)
Determine correspondence between pointFields. Gets passed.
vectorField pointField
pointField is a vectorField.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
labelList identity(const label len, label start=0)
Return an identity map of the given length with (map[i] == i), works like std::iota() but returning a...
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
defineTypeNameAndDebug(combustionModel, 0)
Ostream & decrIndent(Ostream &os)
Decrement the indent level.
Geometric merging of points. See below.
Foam::SubStrings< StringType > splitSpace(const StringType &str)
Split string into sub-strings at whitespace (TAB, NL, VT, FF, CR, SPC)
label lineNumber() const noexcept
Const access to the current stream line number.
static bool split(const std::string &line, std::string &key, std::string &val)
void inplaceRenumber(const labelUList &oldToNew, IntListType &input)
Inplace renumber the values (not the indices) of a list.
vector point
Point is a vector.
#define WarningInFunction
Report a warning using Foam::Warning.
messageStream Info
Information stream (stdout output on master, null elsewhere)
List< label > labelList
A List of labels.
Ostream & incrIndent(Ostream &os)
Increment the indent level.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
forAllConstIters(mixture.phases(), phase)