73 list.append(std::move(tri));
81 const auto* cpp = isA<coupledPolyPatch>(pp);
82 return cpp && cpp->parallel() && !cpp->separated();
88static bool isCollocatedPatch(
const coupledPolyPatch& pp)
90 if (isA<processorPolyPatch>(pp) || isA<cyclicPolyPatch>(pp))
96 <<
"Unhandled coupledPolyPatch type " << pp.type() <<
nl
110 const coupledPolyPatch& pp
114 bitSet collocated(pp.size());
116 if (isA<processorPolyPatch>(pp) || isA<cyclicPolyPatch>(pp))
127 <<
"Unhandled coupledPolyPatch type " << pp.type()
134void Foam::isoSurfacePoint::syncUnseparatedPoints
137 const point& nullValue
142 const polyBoundaryMesh&
patches = mesh_.boundaryMesh();
146 const labelList& procPatches = mesh_.globalData().processorPatches();
149 for (
const label patchi : procPatches)
151 const polyPatch& pp =
patches[patchi];
152 const auto& procPatch = refCast<const processorPolyPatch>(pp);
156 const labelList& meshPts = procPatch.meshPoints();
157 const labelList& nbrPts = procPatch.neighbPoints();
163 const label nbrPointi = nbrPts[pointi];
164 patchInfo[nbrPointi] = pointValues[meshPts[pointi]];
170 procPatch.neighbProcNo()
177 for (
const label patchi : procPatches)
179 const polyPatch& pp =
patches[patchi];
180 const auto& procPatch = refCast<const processorPolyPatch>(pp);
191 procPatch.neighbProcNo()
193 fromNbr >> nbrPatchInfo;
196 const labelList& meshPts = procPatch.meshPoints();
200 const label meshPointi = meshPts[pointi];
203 pointValues[meshPointi],
212 for (
const polyPatch& pp :
patches)
214 const cyclicPolyPatch* cpp = isA<cyclicPolyPatch>(pp);
220 const auto& cycPatch = *cpp;
221 const auto& nbrPatch = cycPatch.neighbPatch();
223 const edgeList& coupledPoints = cycPatch.coupledPoints();
224 const labelList& meshPts = cycPatch.meshPoints();
225 const labelList& nbrMeshPoints = nbrPatch.meshPoints();
232 const edge&
e = coupledPoints[i];
233 half0Values[i] = pointValues[meshPts[
e[0]]];
234 half1Values[i] = pointValues[nbrMeshPoints[
e[1]]];
239 const edge&
e = coupledPoints[i];
240 const label
p0 = meshPts[
e[0]];
241 const label p1 = nbrMeshPoints[
e[1]];
243 minEqOp<point>()(pointValues[
p0], half1Values[i]);
244 minEqOp<point>()(pointValues[p1], half0Values[i]);
250 const globalMeshData& pd = mesh_.globalData();
252 if (pd.nGlobalPoints() > 0)
255 pointField sharedPts(pd.nGlobalPoints(), nullValue);
257 forAll(pd.sharedPointLabels(), i)
259 const label meshPointi = pd.sharedPointLabels()[i];
261 sharedPts[pd.sharedPointAddr()[i]] = pointValues[meshPointi];
269 forAll(pd.sharedPointLabels(), i)
271 const label meshPointi = pd.sharedPointLabels()[i];
272 pointValues[meshPointi] = sharedPts[pd.sharedPointAddr()[i]];
278Foam::scalar Foam::isoSurfacePoint::isoFraction
284 const scalar d = s1-s0;
295bool Foam::isoSurfacePoint::isEdgeOfFaceCut
305 const bool fpLower = (pVals[
f[fp]] < iso_);
310 || fpLower != neiLower
311 || fpLower != (pVals[
f[
f.
fcIndex(fp)]] < iso_)
322void Foam::isoSurfacePoint::getNeighbour
333 const labelList& own = mesh_.faceOwner();
334 const labelList& nei = mesh_.faceNeighbour();
336 if (mesh_.isInternalFace(facei))
338 const label nbr = (own[facei] == celli ? nei[facei] : own[facei]);
339 nbrValue = cVals[nbr];
340 nbrPoint = meshC[nbr];
344 const label bFacei = facei-mesh_.nInternalFaces();
345 const label patchi = boundaryRegion[bFacei];
346 const polyPatch& pp = mesh_.boundaryMesh()[patchi];
347 const label patchFacei = facei-pp.start();
349 nbrValue = cVals.boundaryField()[patchi][patchFacei];
350 nbrPoint = meshC.boundaryField()[patchi][patchFacei];
355void Foam::isoSurfacePoint::calcCutTypes
363 const polyBoundaryMesh&
patches = mesh_.boundaryMesh();
364 const labelList& own = mesh_.faceOwner();
365 const labelList& nei = mesh_.faceNeighbour();
367 faceCutType_.
resize(mesh_.nFaces());
368 faceCutType_ = cutType::NOTCUT;
370 for (label facei = 0; facei < mesh_.nInternalFaces(); ++facei)
372 const scalar ownValue = cVals[own[facei]];
373 const bool ownLower = (ownValue < iso_);
388 const bool neiLower = (nbrValue < iso_);
390 if (ownLower != neiLower)
392 faceCutType_[facei] = cutType::CUT;
398 const face&
f = mesh_.faces()[facei];
400 if (isEdgeOfFaceCut(pVals,
f, ownLower, neiLower))
402 faceCutType_[facei] = cutType::CUT;
407 for (
const polyPatch& pp :
patches)
409 for (
const label facei : pp.range())
411 const scalar ownValue = cVals[own[facei]];
412 const bool ownLower = (ownValue < iso_);
427 const bool neiLower = (nbrValue < iso_);
429 if (ownLower != neiLower)
431 faceCutType_[facei] = cutType::CUT;
437 const face&
f = mesh_.faces()[facei];
439 if (isEdgeOfFaceCut(pVals,
f, ownLower, neiLower))
441 faceCutType_[facei] = cutType::CUT;
448 cellCutType_.
resize(mesh_.nCells());
449 cellCutType_ = cutType::NOTCUT;
454 for (label facei = 0; facei < mesh_.nInternalFaces(); ++facei)
456 if (faceCutType_[facei] == cutType::NOTCUT)
461 if (cellCutType_[own[facei]] == cutType::NOTCUT)
463 cellCutType_[own[facei]] = cutType::CUT;
466 if (cellCutType_[nei[facei]] == cutType::NOTCUT)
468 cellCutType_[nei[facei]] = cutType::CUT;
476 for (label facei = mesh_.nInternalFaces(); facei < mesh_.nFaces(); ++facei)
478 if (faceCutType_[facei] == cutType::NOTCUT)
483 if (cellCutType_[own[facei]] == cutType::NOTCUT)
485 cellCutType_[own[facei]] = cutType::CUT;
492 Pout<<
"isoSurfacePoint : candidate cut cells "
493 << nCutCells_ <<
" / " << mesh_.nCells() <<
endl;
498Foam::point Foam::isoSurfacePoint::calcCentre(
const triSurface&
s)
506 sum +=
s[i].centre(
s.points());
512void Foam::isoSurfacePoint::calcSnappedCc
519 DynamicList<point>& snappedPoints,
526 snappedCc.
setSize(mesh_.nCells());
530 DynamicList<point, 64> localTriPoints(64);
532 forAll(mesh_.cells(), celli)
534 if (cellCutType_[celli] == cutType::CUT)
536 const scalar cVal = cVals[celli];
538 localTriPoints.clear();
545 for (
const label facei : mesh_.cells()[celli])
547 const face&
f = mesh_.faces()[facei];
563 FixedList<scalar, 3>
s;
564 FixedList<point, 3> pt;
567 s[2] = isoFraction(cVal, nbrValue);
568 pt[2] = (1.0-
s[2])*cc[celli] +
s[2]*nbrPoint;
574 s[0] = isoFraction(cVal, pVals[
p0]);
575 pt[0] = (1.0-
s[0])*cc[celli] +
s[0]*pts[
p0];
579 s[1] = isoFraction(cVal, pVals[p1]);
580 pt[1] = (1.0-
s[1])*cc[celli] +
s[1]*pts[p1];
584 (
s[0] >= 0.0 &&
s[0] <= 0.5)
585 && (
s[1] >= 0.0 &&
s[1] <= 0.5)
586 && (
s[2] >= 0.0 &&
s[2] <= 0.5)
589 localTriPoints.append(pt[0]);
590 localTriPoints.append(pt[1]);
591 localTriPoints.append(pt[2]);
598 if (
s[i] >= 0.0 &&
s[i] <= 0.5)
600 otherPointSum += pt[i];
608 if (localTriPoints.size() == 0)
614 snappedCc[celli] = snappedPoints.size();
615 snappedPoints.append(otherPointSum/nOther);
622 else if (localTriPoints.size() == 3)
625 snappedCc[celli] = snappedPoints.size();
626 snappedPoints.append(
sum(localTriPoints)/3);
627 localTriPoints.clear();
652 label
nZones = surf.markZones
660 snappedCc[celli] = snappedPoints.size();
661 snappedPoints.append(calcCentre(surf));
672void Foam::isoSurfacePoint::calcSnappedPoint
674 const bitSet& isBoundaryPoint,
680 DynamicList<point>& snappedPoints,
690 DynamicList<point, 64> localTriPoints(100);
692 forAll(mesh_.pointFaces(), pointi)
694 if (isBoundaryPoint.test(pointi))
703 for (
const label facei :
pFaces)
705 if (faceCutType_[facei] == cutType::CUT)
718 localTriPoints.clear();
722 for (
const label facei :
pFaces)
726 const face&
f = mesh_.faces()[facei];
727 const label own = mesh_.faceOwner()[facei];
744 FixedList<scalar, 4>
s;
745 FixedList<point, 4> pt;
747 label fp =
f.
find(pointi);
748 s[0] = isoFraction(pVals[pointi], cVals[own]);
749 pt[0] = (1.0-
s[0])*pts[pointi] +
s[0]*cc[own];
751 s[1] = isoFraction(pVals[pointi], nbrValue);
752 pt[1] = (1.0-
s[1])*pts[pointi] +
s[1]*nbrPoint;
755 s[2] = isoFraction(pVals[pointi], pVals[nextPointi]);
756 pt[2] = (1.0-
s[2])*pts[pointi] +
s[2]*pts[nextPointi];
759 s[3] = isoFraction(pVals[pointi], pVals[prevPointi]);
760 pt[3] = (1.0-
s[3])*pts[pointi] +
s[3]*pts[prevPointi];
764 (
s[0] >= 0.0 &&
s[0] <= 0.5)
765 && (
s[1] >= 0.0 &&
s[1] <= 0.5)
766 && (
s[2] >= 0.0 &&
s[2] <= 0.5)
769 localTriPoints.append(pt[0]);
770 localTriPoints.append(pt[1]);
771 localTriPoints.append(pt[2]);
775 (
s[0] >= 0.0 &&
s[0] <= 0.5)
776 && (
s[1] >= 0.0 &&
s[1] <= 0.5)
777 && (
s[3] >= 0.0 &&
s[3] <= 0.5)
780 localTriPoints.append(pt[3]);
781 localTriPoints.append(pt[0]);
782 localTriPoints.append(pt[1]);
788 if (
s[i] >= 0.0 &&
s[i] <= 0.5)
790 otherPointSum += pt[i];
796 if (localTriPoints.size() == 0)
802 collapsedPoint[pointi] = otherPointSum/nOther;
805 else if (localTriPoints.size() == 3)
831 label
nZones = surf.markZones
839 collapsedPoint[pointi] = calcCentre(surf);
846 syncUnseparatedPoints(collapsedPoint,
point::max);
849 snappedPoint.setSize(mesh_.nPoints());
852 forAll(collapsedPoint, pointi)
856 snappedPoint[pointi] = snappedPoints.size();
857 snappedPoints.append(collapsedPoint[pointi]);
865 const bool checkDuplicates,
866 const List<point>& triPoints,
871 label nTris = triPoints.
size()/3;
873 if ((triPoints.size() % 3) != 0)
876 <<
"Problem: number of points " << triPoints.size()
893 Pout<<
"isoSurfacePoint : merged from " << triPointReverseMap.size()
894 <<
" down to " << newPoints.size() <<
" unique points." <<
endl;
910 <<
"Merged points contain duplicates"
911 <<
" when merging with distance " << mergeDistance_ <<
endl
912 <<
"merged:" << oldToNew.size() <<
" re-merged:"
913 << newNewPoints.size()
919 List<labelledTri> tris;
921 DynamicList<labelledTri> dynTris(nTris);
923 DynamicList<label> newToOldTri(nTris);
925 for (label oldTriI = 0; oldTriI < nTris; oldTriI++)
929 triPointReverseMap[rawPointi],
930 triPointReverseMap[rawPointi+1],
931 triPointReverseMap[rawPointi+2],
938 newToOldTri.append(oldTriI);
943 triMap.transfer(newToOldTri);
944 tris.transfer(dynTris);
959 DynamicList<label> newToOldTri(tris.size());
963 const labelledTri& tri = tris[triI];
972 label nbrTriI =
pFaces[i];
974 if (nbrTriI > triI && (tris[nbrTriI] == tri))
987 label newTriI = newToOldTri.size();
988 newToOldTri.append(triMap[triI]);
989 tris[newTriI] = tris[triI];
993 triMap.transfer(newToOldTri);
994 tris.setSize(triMap.size());
998 Pout<<
"isoSurfacePoint : merged from " << nTris
999 <<
" down to " << tris.size() <<
" unique triangles." <<
endl;
1009 const labelledTri&
f = surf[facei];
1010 const labelList& fFaces = surf.faceFaces()[facei];
1014 label nbrFacei = fFaces[i];
1016 if (nbrFacei <= facei)
1022 const labelledTri& nbrF = surf[nbrFacei];
1028 <<
" triangle " << facei <<
" vertices " <<
f
1029 <<
" fc:" <<
f.centre(surf.points())
1030 <<
" has the same vertices as triangle " << nbrFacei
1031 <<
" vertices " << nbrF
1032 <<
" fc:" << nbrF.centre(surf.points())
1044void Foam::isoSurfacePoint::trimToPlanes
1046 const PtrList<plane>& planes,
1048 DynamicList<point>& newTriPoints
1052 DynamicList<triPoints> insideTrisA;
1053 storeOp insideOpA(insideTrisA);
1056 DynamicList<triPoints> insideTrisB;
1057 storeOp insideOpB(insideTrisB);
1059 triPointRef::dummyOp dop;
1062 insideOpA(triPoints(tri.a(), tri.b(), tri.c()));
1069 const plane& pln = planes[faceI];
1073 insideTrisB.clear();
1074 for (
const triPoints& tri : insideTrisA)
1076 triPointRef(tri).sliceWithPlane(pln, insideOpB, dop);
1081 insideTrisA.clear();
1082 for (
const triPoints& tri : insideTrisB)
1084 triPointRef(tri).sliceWithPlane(pln, insideOpA, dop);
1091 DynamicList<triPoints>& newTris = (useA ? insideTrisA : insideTrisB);
1093 newTriPoints.reserve(newTriPoints.size() + 3*newTris.size());
1096 for (
const triPoints& tri : newTris)
1098 newTriPoints.append(tri[0]);
1099 newTriPoints.append(tri[1]);
1100 newTriPoints.append(tri[2]);
1105void Foam::isoSurfacePoint::trimToBox
1107 const treeBoundBox& bb,
1108 DynamicList<point>& triPoints,
1109 DynamicList<label>& triMap
1114 Pout<<
"isoSurfacePoint : trimming to " << bb <<
endl;
1122 planes.set(faceI,
new plane(bb.faceCentre(faceI), -
n));
1125 label nTris = triPoints.size()/3;
1127 DynamicList<point> newTriPoints(triPoints.size()/16);
1128 triMap.setCapacity(nTris/16);
1131 for (label triI = 0; triI < nTris; triI++)
1133 const point&
p0 = triPoints[vertI++];
1134 const point& p1 = triPoints[vertI++];
1135 const point& p2 = triPoints[vertI++];
1137 label oldNPoints = newTriPoints.
size();
1145 label nCells = (newTriPoints.size()-oldNPoints)/3;
1146 for (label i = 0; i < nCells; i++)
1148 triMap.append(triI);
1154 Pout<<
"isoSurfacePoint : trimmed from " << nTris
1155 <<
" down to " << triMap.size()
1156 <<
" triangles." <<
endl;
1159 triPoints.transfer(newTriPoints);
1163void Foam::isoSurfacePoint::trimToBox
1165 const treeBoundBox& bb,
1166 DynamicList<point>& triPoints,
1167 DynamicList<label>& triMap,
1170 List<FixedList<label, 3>>& interpolatedOldPoints,
1171 List<FixedList<scalar, 3>>& interpolationWeights
1174 const List<point> oldTriPoints(triPoints);
1177 trimToBox(bb, triPoints, triMap);
1184 label sz = oldTriPoints.size()/100;
1185 DynamicList<label> dynInterpolatedPoints(sz);
1186 DynamicList<FixedList<label, 3>> dynInterpolatedOldPoints(sz);
1187 DynamicList<FixedList<scalar, 3>> dynInterpolationWeights(sz);
1190 triPointMap.setSize(triPoints.size());
1193 label oldTriI = triMap[triI];
1196 for (label i = 0; i < 3; i++)
1198 label pointI = 3*triI+i;
1199 const point& pt = triPoints[pointI];
1202 label matchPointI = -1;
1203 for (label j = 0; j < 3; j++)
1205 label oldPointI = 3*oldTriI+j;
1206 if (pt == oldTriPoints[oldPointI])
1208 matchPointI = oldPointI;
1213 triPointMap[pointI] = matchPointI;
1216 if (matchPointI == -1)
1218 dynInterpolatedPoints.append(pointI);
1220 FixedList<label, 3> oldPoints
1222 {3*oldTriI, 3*oldTriI+1, 3*oldTriI+2}
1224 dynInterpolatedOldPoints.append(oldPoints);
1228 FixedList<scalar, 3> weights({bary.a(), bary.b(), bary.c()});
1230 dynInterpolationWeights.append(weights);
1235 interpolatedPoints.transfer(dynInterpolatedPoints);
1236 interpolatedOldPoints.transfer(dynInterpolatedOldPoints);
1237 interpolationWeights.transfer(dynInterpolationWeights);
1243 const triSurface&
s,
1251 ListOps::createWithValue<bool>(
s.size(), newToOldFaces,
true,
false)
1254 newToOldPoints.setSize(
s.points().size());
1255 oldToNewPoints.setSize(
s.points().size());
1256 oldToNewPoints = -1;
1260 forAll(include, oldFacei)
1262 if (include[oldFacei])
1265 const labelledTri& tri =
s[oldFacei];
1269 label oldPointi = tri[fp];
1271 if (oldToNewPoints[oldPointi] == -1)
1273 oldToNewPoints[oldPointi] = pointi;
1274 newToOldPoints[pointi++] = oldPointi;
1279 newToOldPoints.setSize(pointi);
1284 forAll(newToOldPoints, i)
1286 newPoints[i] =
s.points()[newToOldPoints[i]];
1289 List<labelledTri> newTriangles(newToOldFaces.size());
1294 const labelledTri& tri =
s[newToOldFaces[i]];
1296 newTriangles[i][0] = oldToNewPoints[tri[0]];
1297 newTriangles[i][1] = oldToNewPoints[tri[1]];
1298 newTriangles[i][2] = oldToNewPoints[tri[2]];
1299 newTriangles[i].region() = tri.region();
1303 return triSurface(newTriangles,
s.patches(), newPoints,
true);
1321 cellValues.primitiveField(),
1327 mergeDistance_(params.mergeTol()*mesh_.bounds().
mag())
1329 const bool regularise = (params.
filter() != filterType::NONE);
1334 Pout<<
"isoSurfacePoint:" <<
nl
1338 <<
" isoValue : " << iso <<
nl
1339 <<
" filter : " <<
Switch(regularise) <<
nl
1351 cValsPtr_.reset(adaptPatchFields(
cellValues).ptr());
1390 collocatedFaces(refCast<const coupledPolyPatch>(pp))
1395 if (!isCollocated[i])
1401 else if (isA<emptyPolyPatch>(pp))
1409 bfld.
set(patchi,
nullptr);
1446 "isoSurfacePoint.cutType",
1459 forAll(cellCutType_, celli)
1461 debugFld[celli] = cellCutType_[celli];
1464 Pout<<
"Writing cut types:"
1498 Pout<<
"isoSurfacePoint : shifted " << snappedPoints.
size()
1499 <<
" cell centres to intersection." <<
endl;
1502 label nCellSnaps = snappedPoints.
size();
1520 refCast<const coupledPolyPatch>(pp);
1522 bitSet isCollocated(collocatedFaces(cpp));
1526 if (!isCollocated[i])
1530 isBoundaryPoint.
set(
f);
1540 isBoundaryPoint.
set(
f);
1565 Pout<<
"isoSurfacePoint : shifted " << snappedPoints.
size()-nCellSnaps
1566 <<
" vertices to intersection." <<
endl;
1595 Pout<<
"isoSurfacePoint : generated " << triMeshCells.
size()
1597 <<
" unmerged points." <<
endl;
1614 interpolatedPoints_,
1615 interpolatedOldPoints_,
1616 interpolationWeights_
1618 triMeshCells =
labelField(triMeshCells, trimTriMap);
1624 tmpsurf = stitchTriPoints
1634 Pout<<
"isoSurfacePoint : generated " << triMap.
size()
1635 <<
" merged triangles." <<
endl;
1645 labelList newTriPointMergeMap(nOldPoints, -1);
1646 forAll(trimTriPointMap, trimPointI)
1648 label oldPointI = trimTriPointMap[trimPointI];
1651 label pointI = triPointMergeMap_[trimPointI];
1654 newTriPointMergeMap[oldPointI] = pointI;
1658 triPointMergeMap_.
transfer(newTriPointMergeMap);
1670 Pout<<
"isoSurfacePoint : checking " << tmpsurf.
size()
1671 <<
" triangles for validity." <<
endl;
1679 Pout<<
"Dumping surface to " << stlFile <<
endl;
1680 tmpsurf.
write(stlFile);
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
static constexpr label size() noexcept
Return the number of elements in the FixedList.
Generic GeometricBoundaryField class.
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
label size() const noexcept
The number of elements in table.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
fileName objectPath() const
The complete path + object name.
void transfer(List< T > &list)
void setSize(const label n)
Alias for resize()
void resize(const label len)
Adjust allocated size of list.
static void listCombineAllGather(const List< commsStruct > &comms, List< T > &values, const CombineOp &cop, const int tag, const label comm)
After completion all processors have the same data.
const T * set(const label i) const
Specialization of GeometricField which holds slices of given complete fields in a form that they act ...
A List obtained as a section of another List.
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
fileName path() const
Return path.
static word timeName(const scalar t, const int precision=precision_)
label find(const T &val, label pos=0) const
Find index of the first occurrence of the value.
label rcIndex(const label i) const noexcept
void size(const label n)
Older name for setAddressableSize.
label fcIndex(const label i) const noexcept
static bool & parRun() noexcept
Test if this a parallel run.
static constexpr direction size() noexcept
The number of elements in the VectorSpace = Ncmpts.
A bitSet stores bits (elements with only two states) in packed internal format and supports a variety...
void set(const bitSet &bitset)
Set specified bits from another bitset.
static const FixedList< vector, 6 > faceNormals
The unit normal per face.
The boundaryRegion persistent data saved as a Map<dictionary>.
This boundary condition is not designed to be evaluated; it is assmued that the value is assigned via...
The coupledPolyPatch is an abstract base class for patches that couple regions of the computational d...
A face is a list of labels corresponding to mesh vertices.
A class for handling file names.
Mesh data needed to do the Finite Volume discretisation.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
const Time & time() const
Return the top-level database.
Low-level components common to various iso-surface algorithms.
const scalarField & cellValues() const noexcept
The mesh cell values used for creating the iso-surface.
const scalarField & pVals_
Point values.
labelList meshCells_
For every face, the original cell in mesh.
const polyMesh & mesh_
Reference to mesh.
const scalarField & cVals_
Cell values.
Preferences for controlling iso-surface algorithms.
filterType filter() const noexcept
Get current filter type.
const boundBox & getClipBounds() const noexcept
Get optional clipping bounding box.
scalar mergeTol() const noexcept
Get current merge tolerance.
A surface formed by the iso value. After "Regularised Marching Tetrahedra: improved iso-surface extra...
const Time & time() const noexcept
Return time registry.
static const complex max
complex (VGREAT,VGREAT)
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
virtual const faceList & faces() const
Return raw faces.
const fileName & pointsInstance() const
Return the current instance directory for points.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
virtual const pointField & points() const
Return raw points.
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
A patch is a list of labels that address the faces in the global face list.
const List< T >::subList patchSlice(const UList< T > &l) const
Slice list to patch.
virtual bool coupled() const
Return true if this patch is geometrically coupled (i.e. faces and.
label start() const
Return start label of this patch in the polyMesh face list.
label nBoundaryFaces() const noexcept
Number of boundary faces (== nFaces - nInternalFaces)
const vectorField & faceCentres() const
const vectorField & cellCentres() const
label nPoints() const noexcept
Number of mesh points.
label nCells() const noexcept
Number of mesh cells.
virtual bool write(const bool valid=true) const
Write using setting from DB.
transferModelList & transfer()
Transfer.
Standard boundBox with extra functionality for use in octree.
Triangle storage. Null constructable (unfortunately triangle<point, point> is not)
Triangulated surface description with patch information.
void triFaceFaces(List< face > &plainFaceList) const
Create a list of faces from the triFaces.
void write(Ostream &os) const
Write to Ostream in simple OpenFOAM format.
virtual void swapPoints(pointField &pts)
Swap points. Similar to movePoints, but returns the old points.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
const volScalarField & p0
const polyBoundaryMesh & patches
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Convenience macros for instantiating iso-surface interpolate methods.
#define defineIsoSurfaceInterpolateMethods(ThisClass)
Geometric merging of points. See below.
void invertManyToMany(const label len, const UList< InputIntListType > &input, List< OutputIntListType > &output)
Invert many-to-many.
GeometricField< vector, fvPatchField, volMesh > volVectorField
const dimensionSet dimless
Dimensionless.
void inplaceRenumber(const labelUList &oldToNew, IntListType &input)
Inplace renumber the values (not the indices) of a list.
List< label > labelList
A List of labels.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
List< geometricSurfacePatch > geometricSurfacePatchList
A List of geometricSurfacePatch.
vectorField pointField
pointField is a vectorField.
Barycentric2D< scalar > barycentric2D
A scalar version of the templated Barycentric2D.
vector point
Point is a vector.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
label mergePoints(const PointList &points, labelList &pointToUnique, labelList &uniquePoints, const scalar mergeTol=SMALL, const bool verbose=false)
MinMax< label > minMax(const labelHashSet &set)
Find the min/max values of labelHashSet.
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
List< labelList > labelListList
A List of labelList.
triangle< point, const point & > triPointRef
A triangle using referred points.
errorManip< error > abort(error &err)
Field< label > labelField
Specialisation of Field<T> for label.
static bool collocatedPatch(const polyPatch &pp)
List< bool > boolList
A List of bools.
static constexpr const zero Zero
Global zero (0)
List< edge > edgeList
A List of edges.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
List< surfZone > surfZoneList
constexpr char nl
The newline '\n' character (0x0a)
Info<< "Finished reading KIVA file"<< endl;cellShapeList cellShapes(nPoints);labelList cellZoning(nPoints, -1);const cellModel &hex=cellModel::ref(cellModel::HEX);labelList hexLabels(8);label activeCells=0;labelList pointMap(nPoints);forAll(pointMap, i){ pointMap[i]=i;}for(label i=0;i< nPoints;i++){ if(f[i] > 0.0) { hexLabels[0]=i;hexLabels[1]=i1tab[i];hexLabels[2]=i3tab[i1tab[i]];hexLabels[3]=i3tab[i];hexLabels[4]=i8tab[i];hexLabels[5]=i1tab[i8tab[i]];hexLabels[6]=i3tab[i1tab[i8tab[i]]];hexLabels[7]=i3tab[i8tab[i]];cellShapes[activeCells].reset(hex, hexLabels);edgeList edges=cellShapes[activeCells].edges();forAll(edges, ei) { if(edges[ei].mag(points)< SMALL) { label start=pointMap[edges[ei].start()];while(start !=pointMap[start]) { start=pointMap[start];} label end=pointMap[edges[ei].end()];while(end !=pointMap[end]) { end=pointMap[end];} label minLabel=min(start, end);pointMap[start]=pointMap[end]=minLabel;} } cellZoning[activeCells]=idreg[i];activeCells++;}}cellShapes.setSize(activeCells);cellZoning.setSize(activeCells);forAll(cellShapes, celli){ cellShape &cs=cellShapes[celli];forAll(cs, i) { cs[i]=pointMap[cs[i]];} cs.collapse();}label bcIDs[11]={-1, 0, 2, 4, -1, 5, -1, 6, 7, 8, 9};const label nBCs=12;const word *kivaPatchTypes[nBCs]={ &wallPolyPatch::typeName, &wallPolyPatch::typeName, &wallPolyPatch::typeName, &wallPolyPatch::typeName, &symmetryPolyPatch::typeName, &wedgePolyPatch::typeName, &polyPatch::typeName, &polyPatch::typeName, &polyPatch::typeName, &polyPatch::typeName, &symmetryPolyPatch::typeName, &oldCyclicPolyPatch::typeName};enum patchTypeNames{ PISTON, VALVE, LINER, CYLINDERHEAD, AXIS, WEDGE, INFLOW, OUTFLOW, PRESIN, PRESOUT, SYMMETRYPLANE, CYCLIC};const char *kivaPatchNames[nBCs]={ "piston", "valve", "liner", "cylinderHead", "axis", "wedge", "inflow", "outflow", "presin", "presout", "symmetryPlane", "cyclic"};List< SLList< face > > pFaces[nBCs]
#define forAll(list, i)
Loop across all elements in list.
void operator()(triPoints &&tri)
void operator()(const triPoints &tri)
storeOp(DynamicList< triPoints > &tris)
DynamicList< triPoints > & list