52 { decompositionType::FACE_CENTRE_TRIS,
"faceCentre" },
53 { decompositionType::FACE_DIAG_TRIS,
"faceDiagonal" },
54 { decompositionType::PYRAMID,
"pyramid" },
60void Foam::tetDecomposer::modifyFace
62 polyTopoChange& meshMod,
73 if (nei == -1 || own < nei)
104void Foam::tetDecomposer::addFace
106 polyTopoChange& meshMod,
110 const label masterPointID,
111 const label masterEdgeID,
112 const label masterFaceID,
119 if (nei == -1 || own < nei)
155Foam::label Foam::tetDecomposer::triIndex(
const label facei,
const label fp)
158 const face&
f = mesh_.faces()[facei];
159 const label fp0 =
max(0, mesh_.tetBasePtIs()[facei]);
169 thisTrii =
f.
size()-3;
173 thisTrii = (fp-fp0-1) % (
f.
size()-2);
192 const bitSet& decomposeCell,
196 cellToPoint_.setSize(mesh_.nCells(), -1);
197 forAll(mesh_.cellCentres(), celli)
199 if (decomposeCell[celli])
202 label masterPointi = mesh_.faces()[mesh_.cells()[celli][0]][0];
204 cellToPoint_[celli] = meshMod.
addPoint
206 mesh_.cellCentres()[celli],
216 bitSet decomposeFace(mesh_.nFaces());
218 for (label facei = 0; facei < mesh_.nInternalFaces(); facei++)
220 label own = mesh_.faceOwner()[facei];
221 label nei = mesh_.faceNeighbour()[facei];
222 if (decomposeCell[own] || decomposeCell[nei])
224 decomposeFace[facei] =
true;
228 boolList neiDecomposeCell(mesh_.nBoundaryFaces());
229 forAll(neiDecomposeCell, bFacei)
231 label facei = mesh_.nInternalFaces()+bFacei;
232 label own = mesh_.faceOwner()[facei];
233 neiDecomposeCell[bFacei] = decomposeCell[own];
239 label facei = mesh_.nInternalFaces();
240 facei < mesh_.nFaces();
244 label own = mesh_.faceOwner()[facei];
245 label bFacei = facei-mesh_.nInternalFaces();
246 if (decomposeCell[own] || neiDecomposeCell[bFacei])
248 decomposeFace[facei] =
true;
255 if (decomposeType == FACE_CENTRE_TRIS)
257 faceToPoint_.
setSize(mesh_.nFaces(), -1);
258 forAll(mesh_.faceCentres(), facei)
260 if (decomposeFace[facei])
263 const label masterPointi = mesh_.faces()[facei][0];
265 faceToPoint_[facei] = meshMod.
addPoint
267 mesh_.faceCentres()[facei],
279 faceOwnerCells_.setSize(mesh_.nFaces());
280 faceNeighbourCells_.setSize(mesh_.nFaces());
282 if (decomposeType == FACE_CENTRE_TRIS)
284 forAll(faceOwnerCells_, facei)
286 const face&
f = mesh_.faces()[facei];
288 faceNeighbourCells_[facei].setSize(
f.
size(), -1);
291 else if (decomposeType == FACE_DIAG_TRIS)
294 (void)mesh_.tetBasePtIs();
296 forAll(faceOwnerCells_, facei)
298 const face&
f = mesh_.faces()[facei];
300 faceNeighbourCells_[facei].setSize(
f.
size()-2, -1);
305 forAll(faceOwnerCells_, facei)
307 faceOwnerCells_[facei].setSize(1, -1);
308 faceNeighbourCells_[facei].setSize(1, -1);
315 forAll(mesh_.cells(), celli)
317 const cell& cFaces = mesh_.cells()[celli];
320 bool modifiedCell =
false;
324 label facei = cFaces[cFacei];
325 const face&
f = mesh_.faces()[facei];
330 (mesh_.faceOwner()[facei] == celli)
331 ? faceOwnerCells_[facei]
332 : faceNeighbourCells_[facei]
335 if (decomposeCell[celli])
337 if (decomposeType == FACE_CENTRE_TRIS)
355 mesh_.cellZones().whichZone(celli)
360 else if (decomposeType == FACE_DIAG_TRIS)
362 for (label triI = 0; triI <
f.
size()-2; triI++)
378 mesh_.cellZones().whichZone(celli)
407 mesh_.cellZones().whichZone(celli)
425 forAll(mesh_.faces(), facei)
427 label own = mesh_.faceOwner()[facei];
428 const labelList& addedOwn = faceOwnerCells_[facei];
429 const labelList& addedNei = faceNeighbourCells_[facei];
430 const face&
f = mesh_.faces()[facei];
433 if (facei >= mesh_.nInternalFaces())
435 patchi = mesh_.boundaryMesh().whichPatch(facei);
438 label zoneI = mesh_.faceZones().whichZone(facei);
439 bool zoneFlip =
false;
442 const faceZone& fz = mesh_.faceZones()[zoneI];
449 if (decomposeType == FACE_CENTRE_TRIS && decomposeFace[facei])
501 if (decomposeCell[own])
503 label newOwn = addedOwn[
f.
rcIndex(fp)];
504 label newNei = addedOwn[fp];
527 facei < mesh_.nInternalFaces()
528 && decomposeCell[mesh_.faceNeighbour()[facei]]
531 label newOwn = addedNei[
f.
rcIndex(fp)];
532 label newNei = addedNei[fp];
536 triangle[2] = cellToPoint_[mesh_.faceNeighbour()[facei]];
554 else if (decomposeType == FACE_DIAG_TRIS && decomposeFace[facei])
556 label fp0 =
max(mesh_.tetBasePtIs()[facei], 0);
559 for (label triI = 0; triI <
f.
size()-2; triI++)
561 label nextTri = triI+1;
562 if (nextTri >=
f.
size()-2)
564 nextTri -=
f.
size()-2;
613 if (triI <
f.
size()-3)
615 if (decomposeCell[own])
617 label newOwn = addedOwn[triI];
618 label newNei = addedOwn[nextTri];
642 facei < mesh_.nInternalFaces()
643 && decomposeCell[mesh_.faceNeighbour()[facei]]
646 label newOwn = addedNei[triI];
647 label newNei = addedNei[nextTri];
652 cellToPoint_[mesh_.faceNeighbour()[facei]];
696 forAll(mesh_.cells(), celli)
698 const cell& cFaces = mesh_.cells()[celli];
704 label facei = cFaces[cFacei];
706 if (decomposeCell[celli])
708 const face&
f = mesh_.faces()[facei];
717 if (edgeFnd == edgeToFace.
end())
724 label otherFacei = edgeFnd();
725 const face& otherF = mesh_.faces()[otherFacei];
731 label otherFp = otherF.
find(
p0);
736 else if (otherF.
prevLabel(otherFp) == p1)
738 otherFp = otherF.
rcIndex(otherFp);
748 if (mesh_.faceOwner()[facei] == celli)
762 label thisTet, otherTet;
764 if (decomposeType == FACE_CENTRE_TRIS)
766 if (mesh_.faceOwner()[facei] == celli)
768 thisTet = faceOwnerCells_[facei][fp];
772 thisTet = faceNeighbourCells_[facei][fp];
775 if (mesh_.faceOwner()[otherFacei] == celli)
778 faceOwnerCells_[otherFacei][otherFp];
783 faceNeighbourCells_[otherFacei][otherFp];
786 else if (decomposeType == FACE_DIAG_TRIS)
788 label thisTriI = triIndex(facei, fp);
789 if (mesh_.faceOwner()[facei] == celli)
791 thisTet = faceOwnerCells_[facei][thisTriI];
795 thisTet = faceNeighbourCells_[facei][thisTriI];
798 label otherTriI = triIndex(otherFacei, otherFp);
799 if (mesh_.faceOwner()[otherFacei] == celli)
802 faceOwnerCells_[otherFacei][otherTriI];
807 faceNeighbourCells_[otherFacei][otherTriI];
812 if (mesh_.faceOwner()[facei] == celli)
814 thisTet = faceOwnerCells_[facei][0];
818 thisTet = faceNeighbourCells_[facei][0];
820 if (mesh_.faceOwner()[otherFacei] == celli)
822 otherTet = faceOwnerCells_[otherFacei][0];
827 faceNeighbourCells_[otherFacei][0];
857 forAll(faceOwnerCells_, facei)
861 forAll(faceNeighbourCells_, facei)
Map from edge (expressed as its endpoints) to value. For easier forward declaration it is currently i...
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
bool insert(const Key &key, const T &obj)
Copy insert a new entry, not overwriting existing entries.
iterator find(const Key &key)
Find and return an iterator set at the hashed entry.
void clear()
Clear all entries from table.
iterator end() noexcept
iterator to signal the end (for any HashTable)
void setSize(const label n)
Alias for resize()
void setSize(const label n, unsigned int val=0u)
Alias for resize()
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
A const_iterator for iterating across on values.
A bitSet stores bits (elements with only two states) in packed internal format and supports a variety...
A cell is defined as a list of faces with extra functionality.
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
A subset of mesh faces organised as a primitive patch.
label whichFace(const label globalCellID) const
Helper function to re-direct to zone::localID(...)
const boolList & flipMap() const noexcept
Return face flip map.
A face is a list of labels corresponding to mesh vertices.
label nextLabel(const label i) const
Next vertex on face.
label prevLabel(const label i) const
Previous vertex on face.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
const labelList & reverseCellMap() const
Reverse cell map.
const labelList & reversePointMap() const
Reverse point map.
void updateMesh()
Update for new mesh topology.
Mesh consisting of general polyhedral cells.
Direct mesh changes based on v1.3 polyTopoChange syntax.
label addPoint(const point &pt, const label masterPointID, const label zoneID, const bool inCell)
Add point. Return new point label.
label addCell(const label masterPointID, const label masterEdgeID, const label masterFaceID, const label masterCellID, const label zoneID)
Add cell. Return new cell label.
Decomposes polyMesh into tets or pyramids.
static const Enum< decompositionType > decompositionTypeNames
void setRefinement(const decompositionType decomposeType, const bitSet &decomposeCell, polyTopoChange &meshMod)
Insert all changes into meshMod to convert the polyMesh into.
A triangle primitive used to calculate face normals and swept volumes.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
const volScalarField & p0
#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.
void inplaceRenumber(const labelUList &oldToNew, IntListType &input)
Inplace renumber the values (not the indices) of a list.
errorManip< error > abort(error &err)
#define forAll(list, i)
Loop across all elements in list.