diff --git a/applications/utilities/mesh/advanced/collapseEdges/Make/files b/applications/utilities/mesh/advanced/collapseEdges/Make/files
deleted file mode 100644
index d89ca6e737..0000000000
--- a/applications/utilities/mesh/advanced/collapseEdges/Make/files
+++ /dev/null
@@ -1,3 +0,0 @@
-collapseEdges.C
-
-EXE = $(FOAM_APPBIN)/collapseEdges
diff --git a/applications/utilities/mesh/advanced/collapseEdges/Make/options b/applications/utilities/mesh/advanced/collapseEdges/Make/options
deleted file mode 100644
index 7665565e82..0000000000
--- a/applications/utilities/mesh/advanced/collapseEdges/Make/options
+++ /dev/null
@@ -1,10 +0,0 @@
-EXE_INC = \
- -I$(LIB_SRC)/dynamicMesh/lnInclude \
- -I$(LIB_SRC)/meshTools/lnInclude \
- -I$(LIB_SRC)/finiteVolume/lnInclude
-
-EXE_LIBS = \
- -ldynamicMesh \
- -lmeshTools \
- -lgenericPatches \
- -lfiniteVolume
diff --git a/applications/utilities/mesh/advanced/collapseEdges/collapseEdges.C b/applications/utilities/mesh/advanced/collapseEdges/collapseEdges.C
deleted file mode 100644
index da7e4007c9..0000000000
--- a/applications/utilities/mesh/advanced/collapseEdges/collapseEdges.C
+++ /dev/null
@@ -1,258 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration | Website: https://openfoam.org
- \\ / A nd | Copyright (C) 2011-2023 OpenFOAM Foundation
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see "
- << abort(FatalError);
- }
-
-
- // maintain indirectPatchFaces if it is there (default) or force
- // (if collapseFaceSet option provided)
- word faceSetName("indirectPatchFaces");
- IOobject::readOption readFlag = IOobject::READ_IF_PRESENT;
-
- if (args.optionReadIfPresent("collapseFaceSet", faceSetName))
- {
- readFlag = IOobject::MUST_READ;
- }
-
-
-
- labelIOList pointPriority
- (
- IOobject
- (
- "pointPriority",
- runTime.name(),
- runTime,
- IOobject::READ_IF_PRESENT,
- IOobject::AUTO_WRITE
- ),
- labelList(mesh.nPoints(), labelMin)
- );
- forAll(timeDirs, timeI)
- {
- runTime.setTime(timeDirs[timeI], timeI);
-
- Info<< "Time = " << runTime.userTimeName() << endl;
-
- autoPtr meshFilterPtr;
-
- label nBadFaces = 0;
-
- faceSet indirectPatchFaces
- (
- mesh,
- faceSetName,
- readFlag,
- IOobject::AUTO_WRITE
- );
- Info<< "Read faceSet " << indirectPatchFaces.name()
- << " with "
- << returnReduce(indirectPatchFaces.size(), sumOp())
- << " faces" << endl;
-
-
- {
- meshFilterPtr.set
- (
- new polyMeshFilter(mesh, pointPriority, collapseDict)
- );
- polyMeshFilter& meshFilter = meshFilterPtr();
-
- // newMesh will be empty until it is filtered
- const autoPtr& newMesh = meshFilter.filteredMesh();
-
- // Filter small edges only. This reduces the number of faces so that
- // the face filtering is sped up.
- nBadFaces = meshFilter.filterEdges(0);
- {
- polyTopoChange meshMod(newMesh());
-
- meshMod.changeMesh(mesh, false);
-
- polyMeshFilter::copySets(newMesh(), mesh);
- }
-
- pointPriority = meshFilter.pointPriority();
- }
-
- if (collapseFaceSet)
- {
- meshFilterPtr.reset
- (
- new polyMeshFilter(mesh, pointPriority, collapseDict)
- );
- polyMeshFilter& meshFilter = meshFilterPtr();
-
- const autoPtr& newMesh = meshFilter.filteredMesh();
-
- // Filter faces. Pass in the number of bad faces that are present
- // from the previous edge filtering to use as a stopping criterion.
- meshFilter.filter(indirectPatchFaces);
- {
- polyTopoChange meshMod(newMesh);
-
- meshMod.changeMesh(mesh, false);
-
- polyMeshFilter::copySets(newMesh(), mesh);
- }
-
- pointPriority = meshFilter.pointPriority();
- }
-
- if (collapseFaces)
- {
- meshFilterPtr.reset
- (
- new polyMeshFilter(mesh, pointPriority, collapseDict)
- );
- polyMeshFilter& meshFilter = meshFilterPtr();
-
- const autoPtr& newMesh = meshFilter.filteredMesh();
-
- // Filter faces. Pass in the number of bad faces that are present
- // from the previous edge filtering to use as a stopping criterion.
- meshFilter.filter(nBadFaces);
- {
- polyTopoChange meshMod(newMesh);
-
- meshMod.changeMesh(mesh, false);
-
- polyMeshFilter::copySets(newMesh(), mesh);
- }
-
- pointPriority = meshFilter.pointPriority();
- }
-
- // Write resulting mesh
- if (!overwrite)
- {
- runTime++;
- }
- else
- {
- mesh.setInstance(oldInstance);
- }
-
- Info<< nl << "Writing collapsed mesh to time "
- << runTime.name() << nl << endl;
-
- mesh.write();
- pointPriority.write();
- }
-
- Info<< nl << "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
- << " ClockTime = " << runTime.elapsedClockTime() << " s"
- << nl << endl;
-
- Info<< "End\n" << endl;
-
- return 0;
-}
-
-
-// ************************************************************************* //
diff --git a/applications/utilities/mesh/advanced/modifyMesh/Make/files b/applications/utilities/mesh/advanced/modifyMesh/Make/files
deleted file mode 100644
index 170a86fe0f..0000000000
--- a/applications/utilities/mesh/advanced/modifyMesh/Make/files
+++ /dev/null
@@ -1,4 +0,0 @@
-cellSplitter.C
-modifyMesh.C
-
-EXE = $(FOAM_APPBIN)/modifyMesh
diff --git a/applications/utilities/mesh/advanced/modifyMesh/Make/options b/applications/utilities/mesh/advanced/modifyMesh/Make/options
deleted file mode 100644
index 72321d31b5..0000000000
--- a/applications/utilities/mesh/advanced/modifyMesh/Make/options
+++ /dev/null
@@ -1,8 +0,0 @@
-EXE_INC = \
- -I$(LIB_SRC)/meshTools/lnInclude \
- -I$(LIB_SRC)/dynamicMesh/lnInclude
-
-EXE_LIBS = \
- -lmeshTools \
- -lgenericPatches \
- -ldynamicMesh
diff --git a/applications/utilities/mesh/advanced/modifyMesh/cellSplitter.C b/applications/utilities/mesh/advanced/modifyMesh/cellSplitter.C
deleted file mode 100644
index fefbac70ef..0000000000
--- a/applications/utilities/mesh/advanced/modifyMesh/cellSplitter.C
+++ /dev/null
@@ -1,491 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration | Website: https://openfoam.org
- \\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see & cellToCells
-) const
-{
- label oldOwn = mesh_.faceOwner()[facei];
-
- Map::const_iterator fnd = cellToCells.find(oldOwn);
-
- if (fnd == cellToCells.end())
- {
- // Unsplit cell
- return oldOwn;
- }
- else
- {
- // Look up index of face in the cells' faces.
-
- const labelList& newCells = fnd();
-
- const cell& cFaces = mesh_.cells()[oldOwn];
-
- return newCells[findIndex(cFaces, facei)];
- }
-}
-
-
-Foam::label Foam::cellSplitter::newNeighbour
-(
- const label facei,
- const Map& cellToCells
-) const
-{
- label oldNbr = mesh_.faceNeighbour()[facei];
-
- Map::const_iterator fnd = cellToCells.find(oldNbr);
-
- if (fnd == cellToCells.end())
- {
- // Unsplit cell
- return oldNbr;
- }
- else
- {
- // Look up index of face in the cells' faces.
-
- const labelList& newCells = fnd();
-
- const cell& cFaces = mesh_.cells()[oldNbr];
-
- return newCells[findIndex(cFaces, facei)];
- }
-}
-
-
-// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
-
-// Construct from components
-Foam::cellSplitter::cellSplitter(const polyMesh& mesh)
-:
- mesh_(mesh),
- addedPoints_()
-{}
-
-
-// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
-
-Foam::cellSplitter::~cellSplitter()
-{}
-
-
-// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
-
-void Foam::cellSplitter::setRefinement
-(
- const Map& cellToMidPoint,
- polyTopoChange& meshMod
-)
-{
- addedPoints_.clear();
- addedPoints_.resize(cellToMidPoint.size());
-
-
- //
- // Introduce cellToMidPoints.
- //
-
- forAllConstIter(Map, cellToMidPoint, iter)
- {
- label celli = iter.key();
-
- label anchorPoint = mesh_.cellPoints()[celli][0];
-
- label addedPointi =
- meshMod.setAction
- (
- polyAddPoint
- (
- iter(), // point
- anchorPoint, // master point
- -1, // zone for point
- true // supports a cell
- )
- );
- addedPoints_.insert(celli, addedPointi);
-
- // Pout<< "Added point " << addedPointi
- // << iter() << " in cell " << celli << " with centre "
- // << mesh_.cellCentres()[celli] << endl;
- }
-
-
- //
- // Add cells (first one is modified original cell)
- //
-
- Map cellToCells(cellToMidPoint.size());
-
- forAllConstIter(Map, cellToMidPoint, iter)
- {
- label celli = iter.key();
-
- const cell& cFaces = mesh_.cells()[celli];
-
- // Cells created for this cell.
- labelList newCells(cFaces.size());
-
- // First pyramid is the original cell
- newCells[0] = celli;
-
- // Add other pyramids
- for (label i = 1; i < cFaces.size(); i++)
- {
- label addedCelli =
- meshMod.setAction
- (
- polyAddCell
- (
- -1, // master point
- -1, // master edge
- -1, // master face
- celli, // master cell
- -1 // zone
- )
- );
-
- newCells[i] = addedCelli;
- }
-
- cellToCells.insert(celli, newCells);
-
- // Pout<< "Split cell " << celli
- // << " with centre " << mesh_.cellCentres()[celli] << nl
- // << " faces:" << cFaces << nl
- // << " into :" << newCells << endl;
- }
-
-
- //
- // Introduce internal faces. These go from edges of the cell to the mid
- // point.
- //
-
- forAllConstIter(Map, cellToMidPoint, iter)
- {
- label celli = iter.key();
-
- label midPointi = addedPoints_[celli];
-
- const cell& cFaces = mesh_.cells()[celli];
-
- const labelList& cEdges = mesh_.cellEdges()[celli];
-
- forAll(cEdges, i)
- {
- label edgeI = cEdges[i];
- const edge& e = mesh_.edges()[edgeI];
-
- // Get the faces on the cell using the edge
- label face0, face1;
- meshTools::getEdgeFaces(mesh_, celli, edgeI, face0, face1);
-
- // Get the cells on both sides of the face by indexing into cFaces.
- // (since newly created cells are stored in cFaces order)
- const labelList& newCells = cellToCells[celli];
-
- label cell0 = newCells[findIndex(cFaces, face0)];
- label cell1 = newCells[findIndex(cFaces, face1)];
-
- if (cell0 < cell1)
- {
- // Construct face to midpoint that is pointing away from
- // (pyramid split off from) celli
-
- const face& f0 = mesh_.faces()[face0];
-
- label index = findIndex(f0, e[0]);
-
- bool edgeInFaceOrder = (f0[f0.fcIndex(index)] == e[1]);
-
- // Check if celli is the face owner
-
- face newF(3);
- if (edgeInFaceOrder == (mesh_.faceOwner()[face0] == celli))
- {
- // edge used in face order.
- newF[0] = e[1];
- newF[1] = e[0];
- newF[2] = midPointi;
- }
- else
- {
- newF[0] = e[0];
- newF[1] = e[1];
- newF[2] = midPointi;
- }
-
- // Now newF points away from cell0
- meshMod.setAction
- (
- polyAddFace
- (
- newF, // face
- cell0, // owner
- cell1, // neighbour
- -1, // master point
- -1, // master edge
- face0, // master face for addition
- false, // flux flip
- -1, // patch for face
- -1, // zone for face
- false // face zone flip
- )
- );
- }
- else
- {
- // Construct face to midpoint that is pointing away from
- // (pyramid split off from) celli
-
- const face& f1 = mesh_.faces()[face1];
-
- label index = findIndex(f1, e[0]);
-
- bool edgeInFaceOrder = (f1[f1.fcIndex(index)] == e[1]);
-
- // Check if celli is the face owner
-
- face newF(3);
- if (edgeInFaceOrder == (mesh_.faceOwner()[face1] == celli))
- {
- // edge used in face order.
- newF[0] = e[1];
- newF[1] = e[0];
- newF[2] = midPointi;
- }
- else
- {
- newF[0] = e[0];
- newF[1] = e[1];
- newF[2] = midPointi;
- }
-
- // Now newF points away from cell1
- meshMod.setAction
- (
- polyAddFace
- (
- newF, // face
- cell1, // owner
- cell0, // neighbour
- -1, // master point
- -1, // master edge
- face0, // master face for addition
- false, // flux flip
- -1, // patch for face
- -1, // zone for face
- false // face zone flip
- )
- );
- }
- }
- }
-
-
- //
- // Update all existing faces for split owner or neighbour.
- //
-
-
- // Mark off affected face.
- boolList faceUpToDate(mesh_.nFaces(), true);
-
- forAllConstIter(Map, cellToMidPoint, iter)
- {
- label celli = iter.key();
-
- const cell& cFaces = mesh_.cells()[celli];
-
- forAll(cFaces, i)
- {
- label facei = cFaces[i];
-
- faceUpToDate[facei] = false;
- }
- }
-
- forAll(faceUpToDate, facei)
- {
- if (!faceUpToDate[facei])
- {
- const face& f = mesh_.faces()[facei];
-
- if (mesh_.isInternalFace(facei))
- {
- label newOwn = newOwner(facei, cellToCells);
- label newNbr = newNeighbour(facei, cellToCells);
-
- if (newOwn < newNbr)
- {
- meshMod.setAction
- (
- polyModifyFace
- (
- f,
- facei,
- newOwn, // owner
- newNbr, // neighbour
- false, // flux flip
- -1, // patch for face
- false, // remove from zone
- -1, // zone for face
- false // face zone flip
- )
- );
- }
- else
- {
- meshMod.setAction
- (
- polyModifyFace
- (
- f.reverseFace(),
- facei,
- newNbr, // owner
- newOwn, // neighbour
- false, // flux flip
- -1, // patch for face
- false, // remove from zone
- -1, // zone for face
- false // face zone flip
- )
- );
- }
-
- }
- else
- {
- label newOwn = newOwner(facei, cellToCells);
-
- label patchID, zoneID, zoneFlip;
- getFaceInfo(facei, patchID, zoneID, zoneFlip);
-
- meshMod.setAction
- (
- polyModifyFace
- (
- mesh_.faces()[facei],
- facei,
- newOwn, // owner
- -1, // neighbour
- false, // flux flip
- patchID, // patch for face
- false, // remove from zone
- zoneID, // zone for face
- zoneFlip // face zone flip
- )
- );
- }
-
- faceUpToDate[facei] = true;
- }
- }
-}
-
-
-void Foam::cellSplitter::topoChange(const polyTopoChangeMap& map)
-{
- // Create copy since we're deleting entries. Only if both cell and added
- // point get mapped do they get inserted.
- Map newAddedPoints(addedPoints_.size());
-
- forAllConstIter(Map, addedPoints_, iter)
- {
- label oldCelli = iter.key();
-
- label newCelli = map.reverseCellMap()[oldCelli];
-
- label oldPointi = iter();
-
- label newPointi = map.reversePointMap()[oldPointi];
-
- if (newCelli >= 0 && newPointi >= 0)
- {
- newAddedPoints.insert(newCelli, newPointi);
- }
- }
-
- // Copy
- addedPoints_.transfer(newAddedPoints);
-}
-
-
-// ************************************************************************* //
diff --git a/applications/utilities/mesh/advanced/modifyMesh/cellSplitter.H b/applications/utilities/mesh/advanced/modifyMesh/cellSplitter.H
deleted file mode 100644
index 0c141af14a..0000000000
--- a/applications/utilities/mesh/advanced/modifyMesh/cellSplitter.H
+++ /dev/null
@@ -1,153 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration | Website: https://openfoam.org
- \\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see addedPoints_;
-
-
- // Private Member Functions
-
- //- Get patch and zone info for face
- void getFaceInfo
- (
- const label facei,
- label& patchID,
- label& zoneID,
- label& zoneFlip
- ) const;
-
- //- Find the new owner (if any) of the face.
- label newOwner
- (
- const label facei,
- const Map& cellToCells
- ) const;
-
- //- Find the new neighbour (if any) of the face.
- label newNeighbour
- (
- const label facei,
- const Map& cellToCells
- ) const;
-
-
-public:
-
- //- Runtime type information
- ClassName("cellSplitter");
-
- // Constructors
-
- //- Construct from mesh
- cellSplitter(const polyMesh& mesh);
-
- //- Disallow default bitwise copy construction
- cellSplitter(const cellSplitter&) = delete;
-
-
- //- Destructor
- ~cellSplitter();
-
-
- // Member Functions
-
- // Edit
-
- //- Insert mesh changes into meshMod.
- // cellToMidPoint : cell to cut and position of its new midpoint
- void setRefinement
- (
- const Map& cellToMidPoint,
- polyTopoChange& meshMod
- );
-
- //- Force recalculation of locally stored data on topological change
- void topoChange(const polyTopoChangeMap&);
-
-
- // Access
-
- //- Per cell the mid point added.
- const Map& addedPoints() const
- {
- return addedPoints_;
- }
-
-
- // Member Operators
-
- //- Disallow default bitwise assignment
- void operator=(const cellSplitter&) = delete;
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //
diff --git a/applications/utilities/mesh/advanced/modifyMesh/modifyMesh.C b/applications/utilities/mesh/advanced/modifyMesh/modifyMesh.C
deleted file mode 100644
index 37025aabcd..0000000000
--- a/applications/utilities/mesh/advanced/modifyMesh/modifyMesh.C
+++ /dev/null
@@ -1,678 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration | Website: https://openfoam.org
- \\ / A nd | Copyright (C) 2011-2023 OpenFOAM Foundation
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see > pointsToMove(dict.lookup("pointsToMove"));
- List> edgesToSplit(dict.lookup("edgesToSplit"));
- List> facesToTriangulate
- (
- dict.lookup("facesToTriangulate")
- );
-
- bool cutBoundary =
- (
- pointsToMove.size()
- || edgesToSplit.size()
- || facesToTriangulate.size()
- );
-
- List> edgesToCollapse(dict.lookup("edgesToCollapse"));
-
- bool collapseEdge = edgesToCollapse.size();
-
- List> cellsToPyramidise(dict.lookup("cellsToSplit"));
-
- bool cellsToSplit = cellsToPyramidise.size();
-
- // List>
- // cellsToCreate(dict.lookup("cellsToCreate"));
-
- Info<< "Read from " << dict.name() << nl
- << " Boundary cutting module:" << nl
- << " points to move :" << pointsToMove.size() << nl
- << " edges to split :" << edgesToSplit.size() << nl
- << " faces to triangulate:" << facesToTriangulate.size() << nl
- << " Cell splitting module:" << nl
- << " cells to split :" << cellsToPyramidise.size() << nl
- << " Edge collapsing module:" << nl
- << " edges to collapse :" << edgesToCollapse.size() << nl
- //<< " cells to create :" << cellsToCreate.size() << nl
- << endl;
-
- if
- (
- (cutBoundary && collapseEdge)
- || (cutBoundary && cellsToSplit)
- || (collapseEdge && cellsToSplit)
- )
- {
- FatalErrorInFunction
- << "Used more than one mesh modifying module "
- << "(boundary cutting, cell splitting, edge collapsing)" << nl
- << "Please do them in separate passes." << exit(FatalError);
- }
-
-
-
- // Get calculating engine for all of outside
- const SubList outsideFaces
- (
- mesh.faces(),
- mesh.nFaces() - mesh.nInternalFaces(),
- mesh.nInternalFaces()
- );
-
- primitivePatch allBoundary(outsideFaces, mesh.points());
-
-
- // Look up mesh labels and convert to input for boundaryCutter.
-
- bool validInputs = true;
-
-
- Info<< nl << "Looking up points to move ..." << nl << endl;
- Map pointToPos(pointsToMove.size());
- forAll(pointsToMove, i)
- {
- const Pair& pts = pointsToMove[i];
-
- label pointi = findPoint(allBoundary, pts.first());
-
- if (pointi == -1 || !pointToPos.insert(pointi, pts.second()))
- {
- Info<< "Could not insert mesh point " << pointi
- << " for input point " << pts.first() << nl
- << "Perhaps the point is already marked for moving?" << endl;
- validInputs = false;
- }
- }
-
-
- Info<< nl << "Looking up edges to split ..." << nl << endl;
- Map> edgeToCuts(edgesToSplit.size());
- forAll(edgesToSplit, i)
- {
- const Pair& pts = edgesToSplit[i];
-
- label edgeI = findEdge(mesh, allBoundary, pts.first());
-
- if
- (
- edgeI == -1
- || !edgeToCuts.insert(edgeI, List(1, pts.second()))
- )
- {
- Info<< "Could not insert mesh edge " << edgeI
- << " for input point " << pts.first() << nl
- << "Perhaps the edge is already marked for cutting?" << endl;
-
- validInputs = false;
- }
- }
-
-
- Info<< nl << "Looking up faces to triangulate ..." << nl << endl;
- Map faceToDecompose(facesToTriangulate.size());
- forAll(facesToTriangulate, i)
- {
- const Pair& pts = facesToTriangulate[i];
-
- label facei = findFace(mesh, allBoundary, pts.first());
-
- if (facei == -1 || !faceToDecompose.insert(facei, pts.second()))
- {
- Info<< "Could not insert mesh face " << facei
- << " for input point " << pts.first() << nl
- << "Perhaps the face is already marked for splitting?" << endl;
-
- validInputs = false;
- }
- }
-
-
-
- Info<< nl << "Looking up cells to convert to pyramids around"
- << " cell centre ..." << nl << endl;
- Map cellToPyrCentre(cellsToPyramidise.size());
- forAll(cellsToPyramidise, i)
- {
- const Pair& pts = cellsToPyramidise[i];
-
- label celli = findCell(mesh, pts.first());
-
- if (celli == -1 || !cellToPyrCentre.insert(celli, pts.second()))
- {
- Info<< "Could not insert mesh cell " << celli
- << " for input point " << pts.first() << nl
- << "Perhaps the cell is already marked for splitting?" << endl;
-
- validInputs = false;
- }
- }
-
-
- Info<< nl << "Looking up edges to collapse ..." << nl << endl;
- Map edgeToPos(edgesToCollapse.size());
- forAll(edgesToCollapse, i)
- {
- const Pair& pts = edgesToCollapse[i];
-
- label edgeI = findEdge(mesh, allBoundary, pts.first());
-
- if (edgeI == -1 || !edgeToPos.insert(edgeI, pts.second()))
- {
- Info<< "Could not insert mesh edge " << edgeI
- << " for input point " << pts.first() << nl
- << "Perhaps the edge is already marked for collapsing?" << endl;
-
- validInputs = false;
- }
- }
-
-
-
- if (!validInputs)
- {
- Info<< nl << "There was a problem in one of the inputs in the"
- << " dictionary. Not modifying mesh." << endl;
- }
- else if (cellToPyrCentre.size())
- {
- Info<< nl << "All input cells located. Modifying mesh." << endl;
-
- // Mesh change engine
- cellSplitter cutter(mesh);
-
- // Topo change container
- polyTopoChange meshMod(mesh);
-
- // Insert commands into meshMod
- cutter.setRefinement(cellToPyrCentre, meshMod);
-
- // Do changes
- autoPtr map = meshMod.changeMesh(mesh, false);
-
- if (map().hasMotionPoints())
- {
- mesh.setPoints(map().preMotionPoints());
- }
-
- cutter.topoChange(map());
-
- if (!overwrite)
- {
- runTime++;
- }
- else
- {
- mesh.setInstance(oldInstance);
- }
-
- // Write resulting mesh
- Info<< "Writing modified mesh to time " << runTime.name() << endl;
- mesh.write();
- }
- else if (edgeToPos.size())
- {
- Info<< nl << "All input edges located. Modifying mesh." << endl;
-
- // Mesh change engine
- edgeCollapser cutter(mesh);
-
- const edgeList& edges = mesh.edges();
- const pointField& points = mesh.points();
-
- pointField newPoints(points);
-
- PackedBoolList collapseEdge(mesh.nEdges());
- Map collapsePointToLocation(mesh.nPoints());
-
- // Get new positions and construct collapse network
- forAllConstIter(Map, edgeToPos, iter)
- {
- label edgeI = iter.key();
- const edge& e = edges[edgeI];
-
- collapseEdge[edgeI] = true;
- collapsePointToLocation.set(e[1], points[e[0]]);
-
- newPoints[e[0]] = iter();
- }
-
- // Move master point to destination.
- mesh.setPoints(newPoints);
-
- List allPointInfo;
- const globalIndex globalPoints(mesh.nPoints());
- labelList pointPriority(mesh.nPoints(), 0);
-
- cutter.consistentCollapse
- (
- globalPoints,
- pointPriority,
- collapsePointToLocation,
- collapseEdge,
- allPointInfo
- );
-
- // Topo change container
- polyTopoChange meshMod(mesh);
-
- // Insert
- cutter.setRefinement(allPointInfo, meshMod);
-
- // Do changes
- autoPtr map = meshMod.changeMesh(mesh, false);
-
- if (map().hasMotionPoints())
- {
- mesh.setPoints(map().preMotionPoints());
- }
-
- // Not implemented yet:
- // cutter.topoChange(map());
-
-
- if (!overwrite)
- {
- runTime++;
- }
- else
- {
- mesh.setInstance(oldInstance);
- }
-
- // Write resulting mesh
- Info<< "Writing modified mesh to time " << runTime.name() << endl;
- mesh.write();
- }
- else
- {
- Info<< nl << "All input points located. Modifying mesh." << endl;
-
- // Mesh change engine
- boundaryCutter cutter(mesh);
-
- // Topo change container
- polyTopoChange meshMod(mesh);
-
- // Insert commands into meshMod
- cutter.setRefinement
- (
- pointToPos,
- edgeToCuts,
- Map(0), // Faces to split diagonally
- faceToDecompose, // Faces to triangulate
- meshMod
- );
-
- // Do changes
- autoPtr map = meshMod.changeMesh(mesh, false);
-
- if (map().hasMotionPoints())
- {
- mesh.setPoints(map().preMotionPoints());
- }
-
- cutter.topoChange(map());
-
- if (!overwrite)
- {
- runTime++;
- }
- else
- {
- mesh.setInstance(oldInstance);
- }
-
- // Write resulting mesh
- Info<< "Writing modified mesh to time " << runTime.name() << endl;
- mesh.write();
- }
-
-
- Info<< "\nEnd\n" << endl;
- return 0;
-}
-
-
-// ************************************************************************* //
diff --git a/src/dynamicMesh/Make/files b/src/dynamicMesh/Make/files
index 5bd8d18059..267ab7bd27 100644
--- a/src/dynamicMesh/Make/files
+++ b/src/dynamicMesh/Make/files
@@ -4,7 +4,6 @@ polyTopoChange/polyTopoChange.C
polyTopoChange/addPatchCellLayer.C
polyTopoChange/pointEdgeCollapse/pointEdgeCollapse.C
polyTopoChange/edgeCollapser.C
-polyTopoChange/faceCollapser.C
polyTopoChange/removeCells.C
polyTopoChange/removeFaces.C
polyTopoChange/refinementData.C
@@ -12,7 +11,6 @@ polyTopoChange/refinementDistanceData.C
polyTopoChange/removePoints.C
polyTopoChange/combineFaces.C
polyTopoChange/duplicatePoints.C
-polyTopoChange/tetDecomposer.C
hexRef8 = polyTopoChange/hexRef8
$(hexRef8)/hexRef8.C
@@ -31,11 +29,6 @@ polyTopoChanger/attachDetach/attachInterface.C
polyTopoChanger/attachDetach/detachInterface.C
polyTopoChanger/attachDetach/attachDetachPointMatchMap.C
-polyTopoChanger/layerAdditionRemoval/layerAdditionRemoval.C
-polyTopoChanger/layerAdditionRemoval/setLayerPairing.C
-polyTopoChanger/layerAdditionRemoval/addCellLayer.C
-polyTopoChanger/layerAdditionRemoval/removeCellLayer.C
-
enrichedPatch = polyTopoChanger/slidingInterface/enrichedPatch
$(enrichedPatch)/enrichedPatch.C
$(enrichedPatch)/enrichedPatchPointMap.C
diff --git a/src/dynamicMesh/polyTopoChange/faceCollapser.C b/src/dynamicMesh/polyTopoChange/faceCollapser.C
deleted file mode 100644
index 4818665c7c..0000000000
--- a/src/dynamicMesh/polyTopoChange/faceCollapser.C
+++ /dev/null
@@ -1,505 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration | Website: https://openfoam.org
- \\ / A nd | Copyright (C) 2011-2023 OpenFOAM Foundation
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see & splitEdges,
- const label facei,
- polyTopoChange& meshMod
-) const
-{
- const face& f = mesh_.faces()[facei];
- const labelList& fEdges = mesh_.faceEdges()[facei];
-
- // Space for replaced vertices and split edges.
- DynamicList newFace(10 * f.size());
-
- forAll(f, fp)
- {
- label v0 = f[fp];
-
- newFace.append(v0);
-
- // Look ahead one to get edge.
- label fp1 = f.fcIndex(fp);
-
- label v1 = f[fp1];
-
- // Get split on edge if any.
- label edgeI = findEdge(mesh_.edges(), fEdges, v0, v1);
-
- Map::const_iterator edgeFnd =
- splitEdges.find(edgeI);
-
- if (edgeFnd != splitEdges.end())
- {
- // edgeI has been split (by introducing new vertices).
- // Insert new vertices in face in correct order
- // (splitEdges was constructed to be from edge.start() to end())
-
- const labelList& extraVerts = edgeFnd();
-
- if (v0 == mesh_.edges()[edgeI].start())
- {
- forAll(extraVerts, i)
- {
- newFace.append(extraVerts[i]);
- }
- }
- else
- {
- forAllReverse(extraVerts, i)
- {
- newFace.append(extraVerts[i]);
- }
- }
- }
- }
- face newF(newFace.shrink());
-
- // Pout<< "Modifying face:" << facei << " from " << f << " to " << newFace
- // << endl;
-
- if (newF != f)
- {
- label nei = -1;
-
- label patchi = -1;
-
- if (mesh_.isInternalFace(facei))
- {
- nei = mesh_.faceNeighbour()[facei];
- }
- else
- {
- patchi = mesh_.boundaryMesh().whichPatch(facei);
- }
-
- // Get current zone info
- label zoneID = mesh_.faceZones().whichZone(facei);
-
- bool zoneFlip = false;
-
- if (zoneID >= 0)
- {
- const faceZone& fZone = mesh_.faceZones()[zoneID];
-
- zoneFlip = fZone.flipMap()[fZone.whichFace(facei)];
- }
-
- meshMod.setAction
- (
- polyModifyFace
- (
- newF, // modified face
- facei, // label of face being modified
- mesh_.faceOwner()[facei], // owner
- nei, // neighbour
- false, // face flip
- patchi, // patch for face
- false, // remove from zone
- zoneID, // zone for face
- zoneFlip // face flip in zone
- )
- );
- }
-}
-
-
-// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
-
-// Construct from mesh
-Foam::faceCollapser::faceCollapser(const polyMesh& mesh)
-:
- mesh_(mesh)
-{}
-
-
-// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
-
-void Foam::faceCollapser::setRefinement
-(
- const labelList& faceLabels,
- const labelList& fpStart,
- const labelList& fpEnd,
- polyTopoChange& meshMod
-) const
-{
- const pointField& points = mesh_.points();
- const edgeList& edges = mesh_.edges();
- const faceList& faces = mesh_.faces();
- const labelListList& edgeFaces = mesh_.edgeFaces();
-
-
- // From split edge to newly introduced point(s). Can be more than one per
- // edge!
- Map splitEdges(faceLabels.size());
-
- // Mark faces in any way affect by modifying points/edges. Used later on
- // to prevent having to redo all faces.
- labelHashSet affectedFaces(4*faceLabels.size());
-
-
- //
- // Add/remove vertices and construct mapping
- //
-
- forAll(faceLabels, i)
- {
- const label facei = faceLabels[i];
-
- const face& f = faces[facei];
-
- const label fpA = fpStart[i];
- const label fpB = fpEnd[i];
-
- const point& pA = points[f[fpA]];
- const point& pB = points[f[fpB]];
-
- Pout<< "Face:" << f << " collapsed to fp:" << fpA << ' ' << fpB
- << " with points:" << pA << ' ' << pB
- << endl;
-
- // Create line from fpA to fpB
- linePointRef lineAB(pA, pB);
-
- // Get projections of all vertices onto line.
-
- // Distance(squared) to pA for every point on face.
- SortableList dist(f.size());
-
- dist[fpA] = 0;
- dist[fpB] = magSqr(pB - pA);
-
- // Step from fpA to fpB
- // ~~~~~~~~~~~~~~~~~~~~
- // (by incrementing)
-
- label fpMin1 = fpA;
- label fp = f.fcIndex(fpMin1);
-
- while (fp != fpB)
- {
- // See where fp sorts. Make sure it is above fpMin1!
- pointHit near = lineAB.nearestDist(points[f[fp]]);
-
- scalar w = magSqr(near.rawPoint() - pA);
-
- if (w <= dist[fpMin1])
- {
- // Offset.
- w = dist[fpMin1] + 1e-6*(dist[fpB] - dist[fpA]);
-
- point newPoint
- (
- pA + Foam::sqrt(w / (dist[fpB] - dist[fpA]))*(pB - pA)
- );
-
- Pout<< "Adapting position of vertex " << f[fp] << " on face "
- << f << " from " << near.rawPoint() << " to " << newPoint
- << endl;
-
- near.setPoint(newPoint);
- }
-
- // Responsibility of caller to make sure polyModifyPoint is only
- // called once per point. (so max only one collapse face per
- // edge)
- meshMod.setAction
- (
- polyModifyPoint
- (
- f[fp],
- near.rawPoint(),
- false,
- -1,
- true
- )
- );
-
- dist[fp] = w;
-
- // Step to next vertex.
- fpMin1 = fp;
- fp = f.fcIndex(fpMin1);
- }
-
- // Step from fpA to fpB
- // ~~~~~~~~~~~~~~~~~~~~
- // (by decrementing)
-
- fpMin1 = fpA;
- fp = f.rcIndex(fpMin1);
-
- while (fp != fpB)
- {
- // See where fp sorts. Make sure it is below fpMin1!
- pointHit near = lineAB.nearestDist(points[f[fp]]);
-
- scalar w = magSqr(near.rawPoint() - pA);
-
- if (w <= dist[fpMin1])
- {
- // Offset.
- w = dist[fpMin1] + 1e-6*(dist[fpB] - dist[fpA]);
-
- point newPoint
- (
- pA + Foam::sqrt(w / (dist[fpB] - dist[fpA]))*(pB - pA)
- );
-
- Pout<< "Adapting position of vertex " << f[fp] << " on face "
- << f << " from " << near.rawPoint() << " to " << newPoint
- << endl;
-
- near.setPoint(newPoint);
- }
-
- // Responsibility of caller to make sure polyModifyPoint is only
- // called once per point. (so max only one collapse face per
- // edge)
- meshMod.setAction
- (
- polyModifyPoint
- (
- f[fp],
- near.rawPoint(),
- false,
- -1,
- true
- )
- );
-
- dist[fp] = w;
-
- // Step to previous vertex.
- fpMin1 = fp;
- fp = f.rcIndex(fpMin1);
- }
-
- dist.sort();
-
- // Check that fpB sorts latest.
- if (dist.indices()[dist.size()-1] != fpB)
- {
- OFstream str("conflictingFace.obj");
- meshTools::writeOBJ(str, faceList(1, f), points);
-
- FatalErrorInFunction
- << "Trying to collapse face:" << facei << " vertices:" << f
- << " to edges between vertices " << f[fpA] << " and "
- << f[fpB] << " but " << f[fpB] << " does not seem to be the"
- << " vertex furthest away from " << f[fpA] << endl
- << "Dumped conflicting face to obj file conflictingFace.obj"
- << abort(FatalError);
- }
-
-
- // From fp to index in sort:
- Pout<< "Face:" << f << " fpA:" << fpA << " fpB:" << fpB << nl;
-
- labelList sortedFp(f.size());
- forAll(dist.indices(), i)
- {
- label fp = dist.indices()[i];
-
- Pout<< " fp:" << fp << " distance:" << dist[i] << nl;
-
- sortedFp[fp] = i;
- }
-
- const labelList& fEdges = mesh_.faceEdges()[facei];
-
- // Now look up all edges in the face and see if they get extra
- // vertices inserted and build an edge-to-intersected-points table.
-
- // Order of inserted points is in edge order (from e.start to
- // e.end)
-
- forAll(f, fp)
- {
- label fp1 = f.fcIndex(fp);
-
- // Get index in sorted list
- label sorted0 = sortedFp[fp];
- label sorted1 = sortedFp[fp1];
-
- // Get indices in increasing order
- DynamicList edgePoints(f.size());
-
- // Whether edgePoints are from fp to fp1
- bool fpToFp1;
-
- if (sorted0 < sorted1)
- {
- fpToFp1 = true;
-
- for (label j = sorted0+1; j < sorted1; j++)
- {
- edgePoints.append(f[dist.indices()[j]]);
- }
- }
- else
- {
- fpToFp1 = false;
-
- for (label j = sorted1+1; j < sorted0; j++)
- {
- edgePoints.append(f[dist.indices()[j]]);
- }
- }
-
- if (edgePoints.size())
- {
- edgePoints.shrink();
-
- label edgeI = findEdge(edges, fEdges, f[fp], f[fp1]);
-
- const edge& e = edges[edgeI];
-
- if (fpToFp1 == (f[fp] == e.start()))
- {
- splitEdges.insert(edgeI, edgePoints);
- }
- else
- {
- reverse(edgePoints);
- splitEdges.insert(edgeI, edgePoints);
- }
-
- // Mark all faces affected
- insert(edgeFaces[edgeI], facei, affectedFaces);
- }
- }
- }
-
- forAllConstIter(Map, splitEdges, iter)
- {
- Pout<< "Split edge:" << iter.key()
- << " verts:" << mesh_.edges()[iter.key()]
- << " in:" << nl;
- const labelList& edgePoints = iter();
-
- forAll(edgePoints, i)
- {
- Pout<< " " << edgePoints[i] << nl;
- }
- }
-
-
- //
- // Remove faces.
- //
-
- forAll(faceLabels, i)
- {
- const label facei = faceLabels[i];
-
- meshMod.setAction(polyRemoveFace(facei));
-
- // Update list of faces we still have to modify
- affectedFaces.erase(facei);
- }
-
-
- //
- // Modify faces affected (but not removed)
- //
-
- forAllConstIter(labelHashSet, affectedFaces, iter)
- {
- filterFace(splitEdges, iter.key(), meshMod);
- }
-}
-
-
-// ************************************************************************* //
diff --git a/src/dynamicMesh/polyTopoChange/faceCollapser.H b/src/dynamicMesh/polyTopoChange/faceCollapser.H
deleted file mode 100644
index 39f2a695fd..0000000000
--- a/src/dynamicMesh/polyTopoChange/faceCollapser.H
+++ /dev/null
@@ -1,164 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration | Website: https://openfoam.org
- \\ / A nd | Copyright (C) 2011-2023 OpenFOAM Foundation
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see & splitEdges,
- const label facei,
- polyTopoChange& meshMod
- ) const;
-
-
-public:
-
- //- Runtime type information
- ClassName("faceCollapser");
-
-
- // Constructors
-
- //- Construct from mesh.
- faceCollapser(const polyMesh& mesh);
-
- //- Disallow default bitwise copy construction
- faceCollapser(const faceCollapser&) = delete;
-
-
- // Member Functions
-
- // Edit
-
- //- Collapse faces along endpoints. Play commands into
- // polyTopoChange to create mesh.
- void setRefinement
- (
- const labelList& faceLabels,
- const labelList& fpA,
- const labelList& fpB,
- polyTopoChange&
- ) const;
-
- //- Update stored quantities for new mesh labels.
- void topoChange(const polyTopoChangeMap&)
- {}
-
-
- // Member Operators
-
- //- Disallow default bitwise assignment
- void operator=(const faceCollapser&) = delete;
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //
diff --git a/src/dynamicMesh/polyTopoChange/tetDecomposer.C b/src/dynamicMesh/polyTopoChange/tetDecomposer.C
deleted file mode 100644
index 0c9cfff0a0..0000000000
--- a/src/dynamicMesh/polyTopoChange/tetDecomposer.C
+++ /dev/null
@@ -1,713 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration | Website: https://openfoam.org
- \\ / A nd | Copyright (C) 2012-2023 OpenFOAM Foundation
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see ::names[] =
- {
- "faceCentre",
- "faceDiagonal"
- };
-
- const NamedEnum
- tetDecomposer::decompositionTypeNames;
-}
-
-// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
-
-void Foam::tetDecomposer::modifyFace
-(
- polyTopoChange& meshMod,
- const face& f,
- const label facei,
- const label own,
- const label nei,
- const label patchi,
- const label zoneI,
- const bool zoneFlip
-) const
-{
- // First usage of face. Modify.
- if (nei == -1 || own < nei)
- {
- meshMod.modifyFace
- (
- f, // modified face
- facei, // label of face
- own, // owner
- nei, // neighbour
- false, // face flip
- patchi, // patch for face
- zoneI, // zone for face
- zoneFlip // face flip in zone
- );
- }
- else
- {
- meshMod.modifyFace
- (
- f.reverseFace(), // modified face
- facei, // label of face
- nei, // owner
- own, // neighbour
- true, // face flip
- patchi, // patch for face
- zoneI, // zone for face
- !zoneFlip // face flip in zone
- );
- }
-}
-
-
-void Foam::tetDecomposer::addFace
-(
- polyTopoChange& meshMod,
- const face& f,
- const label own,
- const label nei,
- const label masterPointID,
- const label masterEdgeID,
- const label masterFaceID,
- const label patchi,
- const label zoneI,
- const bool zoneFlip
-) const
-{
- // Second or more usage of face. Add.
- if (nei == -1 || own < nei)
- {
- meshMod.addFace
- (
- f, // modified face
- own, // owner
- nei, // neighbour
- masterPointID, // master point
- masterEdgeID, // master edge
- masterFaceID, // master face
- false, // face flip
- patchi, // patch for face
- zoneI, // zone for face
- zoneFlip // face flip in zone
- );
- }
- else
- {
- meshMod.addFace
- (
- f.reverseFace(), // modified face
- nei, // owner
- own, // neighbour
- masterPointID, // master point
- masterEdgeID, // master edge
- masterFaceID, // master face
- true, // face flip
- patchi, // patch for face
- zoneI, // zone for face
- !zoneFlip // face flip in zone
- );
- }
-}
-
-
-// Work out triangle index given the starting vertex in the face
-Foam::label Foam::tetDecomposer::triIndex(const label facei, const label fp)
-const
-{
- const face& f = mesh_.faces()[facei];
- const label fp0 = mesh_.tetBasePtIs()[facei];
-
- // Work out triangle index on this face
- label thisTriI;
- if (fp == fp0)
- {
- thisTriI = 0;
- }
- else if (fp == f.rcIndex(fp0))
- {
- thisTriI = f.size()-3;
- }
- else
- {
- thisTriI = (fp-fp0-1) % (f.size()-2);
- }
- return thisTriI;
-}
-
-
-// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
-
-Foam::tetDecomposer::tetDecomposer(const polyMesh& mesh)
-:
- mesh_(mesh)
-{}
-
-
-// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
-
-void Foam::tetDecomposer::setRefinement
-(
- const decompositionType decomposeType,
- polyTopoChange& meshMod
-)
-{
- cellToPoint_.setSize(mesh_.nCells());
- forAll(mesh_.cellCentres(), celli)
- {
- // Any point on the cell
- label masterPointi = mesh_.faces()[mesh_.cells()[celli][0]][0];
-
- cellToPoint_[celli] = meshMod.addPoint
- (
- mesh_.cellCentres()[celli],
- masterPointi,
- -1,
- true
- );
- }
-
-
- // Add face centre points
- if (decomposeType == FACE_CENTRE_TRIS)
- {
- faceToPoint_.setSize(mesh_.nFaces());
- forAll(mesh_.faceCentres(), facei)
- {
- // Any point on the face
- const label masterPointi = mesh_.faces()[facei][0];
-
- faceToPoint_[facei] = meshMod.addPoint
- (
- mesh_.faceCentres()[facei],
- masterPointi,
- -1,
- true
- );
- }
- }
-
-
- // Per face, per point (faceCentre) or triangle (faceDiag) the added cell
- faceOwnerCells_.setSize(mesh_.nFaces());
- faceNeighbourCells_.setSize(mesh_.nFaces());
-
- if (decomposeType == FACE_CENTRE_TRIS)
- {
- forAll(faceOwnerCells_, facei)
- {
- const face& f = mesh_.faces()[facei];
- faceOwnerCells_[facei].setSize(f.size(), -1);
- faceNeighbourCells_[facei].setSize(f.size(), -1);
- }
- }
- else
- {
- // Force construction of diagonal decomposition
- (void)mesh_.tetBasePtIs();
-
- forAll(faceOwnerCells_, facei)
- {
- const face& f = mesh_.faces()[facei];
- faceOwnerCells_[facei].setSize(f.size()-2, -1);
- faceNeighbourCells_[facei].setSize(f.size()-2, -1);
- }
- }
-
-
- forAll(mesh_.cells(), celli)
- {
- const cell& cFaces = mesh_.cells()[celli];
-
- EdgeMap edgeToFace(8*cFaces.size());
-
- forAll(cFaces, cFacei)
- {
- label facei = cFaces[cFacei];
- const face& f = mesh_.faces()[facei];
-
- // Get reference to either owner or neighbour
- labelList& added =
- (
- (mesh_.faceOwner()[facei] == celli)
- ? faceOwnerCells_[facei]
- : faceNeighbourCells_[facei]
- );
-
- if (decomposeType == FACE_CENTRE_TRIS)
- {
- forAll(f, fp)
- {
- if (cFacei == 0 && fp == 0)
- {
- // Reuse cell itself
- added[fp] = celli;
- }
- else
- {
- added[fp] = meshMod.addCell
- (
- -1, // masterPoint
- -1, // masterEdge
- -1, // masterFace
- celli, // masterCell
- mesh_.cellZones().whichZone(celli)
- );
- }
- }
- }
- else
- {
- for (label triI = 0; triI < f.size()-2; triI++)
- {
- if (cFacei == 0 && triI == 0)
- {
- // Reuse cell itself
- added[triI] = celli;
- }
- else
- {
- added[triI] = meshMod.addCell
- (
- -1, // masterPoint
- -1, // masterEdge
- -1, // masterFace
- celli, // masterCell
- mesh_.cellZones().whichZone(celli)
- );
- }
- }
- }
- }
- }
-
-
-
- // Add triangle faces
- face triangle(3);
-
- forAll(mesh_.faces(), facei)
- {
- label own = mesh_.faceOwner()[facei];
- const labelList& addedOwn = faceOwnerCells_[facei];
- const labelList& addedNei = faceNeighbourCells_[facei];
- const face& f = mesh_.faces()[facei];
-
- label patchi = -1;
- if (facei >= mesh_.nInternalFaces())
- {
- patchi = mesh_.boundaryMesh().whichPatch(facei);
- }
-
- label zoneI = mesh_.faceZones().whichZone(facei);
- bool zoneFlip = false;
- if (zoneI != -1)
- {
- const faceZone& fz = mesh_.faceZones()[zoneI];
- zoneFlip = fz.flipMap()[fz.whichFace(facei)];
- }
-
-
- if (decomposeType == FACE_CENTRE_TRIS)
- {
- forAll(f, fp)
- {
- // 1. Front triangle (decomposition of face itself)
- // (between owner and neighbour cell)
- {
- triangle[0] = f[fp];
- triangle[1] = f[f.fcIndex(fp)];
- triangle[2] = faceToPoint_[facei];
-
- if (fp == 0)
- {
- modifyFace
- (
- meshMod,
- triangle,
- facei,
- addedOwn[fp],
- addedNei[fp],
- patchi,
- zoneI,
- zoneFlip
- );
- }
- else
- {
- addFace
- (
- meshMod,
- triangle,
- addedOwn[fp],
- addedNei[fp],
- -1, // point
- -1, // edge
- facei, // face
- patchi,
- zoneI,
- zoneFlip
- );
- }
- }
-
-
- // 2. Within owner cell - to cell centre
- {
- label newOwn = addedOwn[f.rcIndex(fp)];
- label newNei = addedOwn[fp];
-
- triangle[0] = f[fp];
- triangle[1] = cellToPoint_[own];
- triangle[2] = faceToPoint_[facei];
-
- addFace
- (
- meshMod,
- triangle,
- newOwn,
- newNei,
- f[fp], // point
- -1, // edge
- -1, // face
- -1, // patchi
- zoneI,
- zoneFlip
- );
- }
- // 2b. Within neighbour cell - to cell centre
- if (facei < mesh_.nInternalFaces())
- {
- label newOwn = addedNei[f.rcIndex(fp)];
- label newNei = addedNei[fp];
-
- triangle[0] = f[fp];
- triangle[1] = faceToPoint_[facei];
- triangle[2] = cellToPoint_[mesh_.faceNeighbour()[facei]];
-
- addFace
- (
- meshMod,
- triangle,
- newOwn,
- newNei,
- f[fp], // point
- -1, // edge
- -1, // face
- -1, // patchi
- zoneI,
- zoneFlip
- );
- }
- }
- }
- else
- {
- label fp0 = mesh_.tetBasePtIs()[facei];
- label fp = f.fcIndex(fp0);
-
- for (label triI = 0; triI < f.size()-2; triI++)
- {
- label nextTri = triI+1;
- if (nextTri >= f.size()-2)
- {
- nextTri -= f.size()-2;
- }
- label nextFp = f.fcIndex(fp);
-
-
- // Triangle triI consisting of f[fp0], f[fp], f[nextFp]
-
-
- // 1. Front triangle (decomposition of face itself)
- // (between owner and neighbour cell)
- {
- triangle[0] = f[fp0];
- triangle[1] = f[fp];
- triangle[2] = f[nextFp];
-
- if (triI == 0)
- {
- modifyFace
- (
- meshMod,
- triangle,
- facei,
- addedOwn[triI],
- addedNei[triI],
- patchi,
- zoneI,
- zoneFlip
- );
- }
- else
- {
- addFace
- (
- meshMod,
- triangle,
- addedOwn[triI],
- addedNei[triI],
- -1, // point
- -1, // edge
- facei, // face
- patchi,
- zoneI,
- zoneFlip
- );
- }
- }
-
-
- // 2. Within owner cell - diagonal to cell centre
- if (triI < f.size()-3)
- {
- label newOwn = addedOwn[triI];
- label newNei = addedOwn[nextTri];
-
- triangle[0] = f[fp0];
- triangle[1] = f[nextFp];
- triangle[2] = cellToPoint_[own];
-
- addFace
- (
- meshMod,
- triangle,
- newOwn,
- newNei,
- f[fp], // point
- -1, // edge
- -1, // face
- -1, // patchi
- zoneI,
- zoneFlip
- );
-
- // 2b. Within neighbour cell - to cell centre
- if (facei < mesh_.nInternalFaces())
- {
- label newOwn = addedNei[triI];
- label newNei = addedNei[nextTri];
-
- triangle[0] = f[nextFp];
- triangle[1] = f[fp0];
- triangle[2] =
- cellToPoint_[mesh_.faceNeighbour()[facei]];
-
- addFace
- (
- meshMod,
- triangle,
- newOwn,
- newNei,
- f[fp], // point
- -1, // edge
- -1, // face
- -1, // patchi
- zoneI,
- zoneFlip
- );
- }
- }
-
-
- fp = nextFp;
- }
- }
- }
-
-
-
- // Add triangles for all edges.
- EdgeMap edgeToFace;
-
- forAll(mesh_.cells(), celli)
- {
- const cell& cFaces = mesh_.cells()[celli];
-
- edgeToFace.clear();
-
- forAll(cFaces, cFacei)
- {
- label facei = cFaces[cFacei];
-
- label zoneI = mesh_.faceZones().whichZone(facei);
- bool zoneFlip = false;
- if (zoneI != -1)
- {
- const faceZone& fz = mesh_.faceZones()[zoneI];
- zoneFlip = fz.flipMap()[fz.whichFace(facei)];
- }
-
- const face& f = mesh_.faces()[facei];
- // const labelList& fEdges = mesh_.faceEdges()[facei];
- forAll(f, fp)
- {
- label p0 = f[fp];
- label p1 = f[f.fcIndex(fp)];
- const edge e(p0, p1);
-
- EdgeMap::const_iterator edgeFnd = edgeToFace.find(e);
- if (edgeFnd == edgeToFace.end())
- {
- edgeToFace.insert(e, facei);
- }
- else
- {
- // Found the other face on the edge.
- label otherFacei = edgeFnd();
- const face& otherF = mesh_.faces()[otherFacei];
-
- // Found the other face on the edge. Note that since
- // we are looping in the same order the tets added for
- // otherFacei will be before those of facei
-
- label otherFp = findIndex(otherF, p0);
- if (otherF.nextLabel(otherFp) == p1)
- {
- // ok. otherFp is first vertex of edge.
- }
- else if (otherF.prevLabel(otherFp) == p1)
- {
- otherFp = otherF.rcIndex(otherFp);
- }
- else
- {
- FatalErrorInFunction
- << "problem." << abort(FatalError);
- }
-
-
- // Triangle from edge to cell centre
- if (mesh_.faceOwner()[facei] == celli)
- {
- triangle[0] = p0;
- triangle[1] = p1;
- triangle[2] = cellToPoint_[celli];
- }
- else
- {
- triangle[0] = p1;
- triangle[1] = p0;
- triangle[2] = cellToPoint_[celli];
- }
-
- // Determine tets on either side
- label thisTet, otherTet;
-
- if (decomposeType == FACE_CENTRE_TRIS)
- {
- if (mesh_.faceOwner()[facei] == celli)
- {
- thisTet = faceOwnerCells_[facei][fp];
- }
- else
- {
- thisTet = faceNeighbourCells_[facei][fp];
- }
-
- if (mesh_.faceOwner()[otherFacei] == celli)
- {
- otherTet = faceOwnerCells_[otherFacei][otherFp];
- }
- else
- {
- otherTet =
- faceNeighbourCells_[otherFacei][otherFp];
- }
- }
- else
- {
- label thisTriI = triIndex(facei, fp);
- if (mesh_.faceOwner()[facei] == celli)
- {
- thisTet = faceOwnerCells_[facei][thisTriI];
- }
- else
- {
- thisTet = faceNeighbourCells_[facei][thisTriI];
- }
-
- label otherTriI = triIndex(otherFacei, otherFp);
- if (mesh_.faceOwner()[otherFacei] == celli)
- {
- otherTet = faceOwnerCells_[otherFacei][otherTriI];
- }
- else
- {
- otherTet =
- faceNeighbourCells_[otherFacei][otherTriI];
- }
- }
-
-
- addFace
- (
- meshMod,
- triangle,
- otherTet,
- thisTet,
- -1, // masterPoint
- -1, // fEdges[fp], // masterEdge
- facei, // masterFace
- -1, // patchi
- zoneI,
- zoneFlip
- );
- }
- }
- }
- }
-}
-
-
-void Foam::tetDecomposer::topoChange(const polyTopoChangeMap& map)
-{
- inplaceRenumber(map.reversePointMap(), cellToPoint_);
- inplaceRenumber(map.reversePointMap(), faceToPoint_);
-
- forAll(faceOwnerCells_, facei)
- {
- inplaceRenumber(map.reverseCellMap(), faceOwnerCells_[facei]);
- }
- forAll(faceNeighbourCells_, facei)
- {
- inplaceRenumber(map.reverseCellMap(), faceNeighbourCells_[facei]);
- }
-}
-
-
-// ************************************************************************* //
diff --git a/src/dynamicMesh/polyTopoChange/tetDecomposer.H b/src/dynamicMesh/polyTopoChange/tetDecomposer.H
deleted file mode 100644
index e3f6daf03c..0000000000
--- a/src/dynamicMesh/polyTopoChange/tetDecomposer.H
+++ /dev/null
@@ -1,206 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration | Website: https://openfoam.org
- \\ / A nd | Copyright (C) 2012-2023 OpenFOAM Foundation
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see decompositionTypeNames;
-
-
-private:
-
- // Private Data
-
- const polyMesh& mesh_;
-
- //- From cell to tet point
- labelList cellToPoint_;
-
- //- From face to tet point
- labelList faceToPoint_;
-
-
- // Per face, per point (faceCentre) or triangle (faceDiag)
- // the added tet on the owner side
- labelListList faceOwnerCells_;
-
- // Per face, per point (faceCentre) or triangle (faceDiag)
- // the added tet on the neighbour side
- labelListList faceNeighbourCells_;
-
-
- // Private Member Functions
-
- //- Modify a face
- void modifyFace
- (
- polyTopoChange& meshMod,
- const face& f,
- const label facei,
- const label own,
- const label nei,
- const label patchi,
- const label zoneI,
- const bool zoneFlip
- ) const;
-
- //- Add a face
- void addFace
- (
- polyTopoChange& meshMod,
- const face& f,
- const label own,
- const label nei,
- const label masterPointID,
- const label masterEdgeID,
- const label masterFaceID,
- const label patchi,
- const label zoneI,
- const bool zoneFlip
- ) const;
-
- //- Work out triangle index given the starting vertex in the face
- label triIndex(const label facei, const label fp) const;
-
-
-public:
-
- //- Runtime type information
- ClassName("tetDecomposer");
-
-
- // Constructors
-
- //- Construct from mesh
- tetDecomposer(const polyMesh&);
-
- //- Disallow default bitwise copy construction
- tetDecomposer(const tetDecomposer&) = delete;
-
-
- // Member Functions
-
- // Access
-
- //- From cell to tet point
- const labelList& cellToPoint() const
- {
- return cellToPoint_;
- }
-
- //- From face to tet point
- const labelList& faceToPoint() const
- {
- return faceToPoint_;
- }
-
-
- //- Per face, per point (faceCentre) or triangle (faceDiag)
- // the added tet on the owner side
- const labelListList& faceOwnerCells() const
- {
- return faceOwnerCells_;
- }
-
- //- Per face, per point (faceCentre) or triangle (faceDiag)
- // the added tet on the neighbour side
- const labelListList& faceNeighbourCells() const
- {
- return faceNeighbourCells_;
- }
-
-
- // Edit
-
- //- Insert all changes into meshMod to convert the polyMesh into
- // tets.
- void setRefinement
- (
- const decompositionType decomposeType,
- polyTopoChange& meshMod
- );
-
- //- Force recalculation of locally stored data on topological change
- void topoChange(const polyTopoChangeMap&);
-
-
- // Member Operators
-
- //- Disallow default bitwise assignment
- void operator=(const tetDecomposer&) = delete;
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //
diff --git a/src/dynamicMesh/polyTopoChanger/layerAdditionRemoval/addCellLayer.C b/src/dynamicMesh/polyTopoChanger/layerAdditionRemoval/addCellLayer.C
deleted file mode 100644
index 35e297ed6f..0000000000
--- a/src/dynamicMesh/polyTopoChanger/layerAdditionRemoval/addCellLayer.C
+++ /dev/null
@@ -1,637 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration | Website: https://openfoam.org
- \\ / A nd | Copyright (C) 2011-2023 OpenFOAM Foundation
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see Foam::layerAdditionRemoval::extrusionDir() const
-{
- const polyMesh& mesh = topoChanger().mesh();
- const primitiveFacePatch& masterFaceLayer =
- mesh.faceZones()[faceZoneID_.index()]();
-
- const pointField& points = mesh.points();
- const labelList& mp = masterFaceLayer.meshPoints();
-
- tmp textrusionDir(new vectorField(mp.size()));
- vectorField& extrusionDir = textrusionDir.ref();
-
- if (setLayerPairing())
- {
- if (debug)
- {
- Pout<< "void layerAdditionRemoval::extrusionDir() const "
- << " for object " << name() << " : "
- << "Using edges for point insertion" << endl;
- }
-
- // Detected a valid layer. Grab the point and face collapse mapping
- const labelList& ptc = pointsPairing();
-
- forAll(extrusionDir, mpI)
- {
- extrusionDir[mpI] = points[ptc[mpI]] - points[mp[mpI]];
- }
- }
- else
- {
- if (debug)
- {
- Pout<< "void layerAdditionRemoval::extrusionDir() const "
- << " for object " << name() << " : "
- << "A valid layer could not be found in front of "
- << "the addition face layer. Using face-based "
- << "point normals for point addition"
- << endl;
- }
-
- extrusionDir = minLayerThickness_*masterFaceLayer.pointNormals();
- }
-
- return textrusionDir;
-}
-
-
-void Foam::layerAdditionRemoval::addCellLayer
-(
- polyTopoChange& ref
-) const
-{
- // Insert the layer addition instructions into the topological change
-
- // Algorithm:
- // 1. For every point in the master zone add a new point
- // 2. For every face in the master zone add a cell
- // 3. Go through all the edges of the master zone. For all internal edges,
- // add a face with the correct orientation and make it internal.
- // For all boundary edges, find the patch it belongs to and add the face
- // to this patch
- // 4. Create a set of new faces on the patch image and assign them to be
- // between the old master cells and the newly created cells
- // 5. Modify all the faces in the patch such that they are located
- // between the old slave cells and newly created cells
-
- if (debug)
- {
- Pout<< "void layerAdditionRemoval::addCellLayer("
- << "polyTopoChange& ref) const for object " << name() << " : "
- << "Adding cell layer" << endl;
- }
-
- // Create the points
-
- const polyMesh& mesh = topoChanger().mesh();
- const primitiveFacePatch& masterFaceLayer =
- mesh.faceZones()[faceZoneID_.index()]();
-
- const pointField& points = mesh.points();
- const labelList& mp = masterFaceLayer.meshPoints();
-
- // Get the extrusion direction for the added points
-
- const vectorField pointOffsets(extrusionDir());
-
- // Add the new points
- labelList addedPoints(mp.size());
-
- forAll(mp, pointi)
- {
- // Add the point nominal thickness away from the master zone point
- // and grab the label
- addedPoints[pointi] =
- ref.setAction
- (
- polyAddPoint
- (
- points[mp[pointi]] // point
- + addDelta_*pointOffsets[pointi],
- mp[pointi], // master point
- -1, // zone for point
- true // supports a cell
- )
- );
- }
-
- if (debug > 1)
- {
- Pout<< "mp: " << mp << " addedPoints: " << addedPoints << endl;
- }
-
- // Create the cells
-
- const labelList& mc =
- mesh.faceZones()[faceZoneID_.index()].masterCells();
- const labelList& sc =
- mesh.faceZones()[faceZoneID_.index()].slaveCells();
-
- const labelList& mf = mesh.faceZones()[faceZoneID_.index()];
- const boolList& mfFlip = mesh.faceZones()[faceZoneID_.index()].flipMap();
-
- const labelList& own = mesh.faceOwner();
- const labelList& nei = mesh.faceNeighbour();
-
- labelList addedCells(mf.size());
-
- forAll(mf, facei)
- {
- label celli = mc[facei];
- label zoneI = mesh.cellZones().whichZone(celli);
-
- addedCells[facei] =
- ref.setAction
- (
- polyAddCell
- (
- -1, // master point
- -1, // master edge
- mf[facei], // master face
- -1, // master cell
- zoneI // zone for cell
- )
- );
- }
-
- // Create the new faces
-
- // Grab the local faces from the master zone
- const faceList& zoneFaces = masterFaceLayer.localFaces();
-
- // Create the new faces by copying the master zone. All the added
- // faces need to point into the newly created cells, which means
- // all the faces from the master layer are flipped. The flux flip
- // is determined from the original master layer face and the face
- // owner: if the master cell is equal to the face owner the flux
- // remains the same; otherwise it is flipped
-
- forAll(zoneFaces, facei)
- {
- const face oldFace = zoneFaces[facei].reverseFace();
-
- face newFace(oldFace.size());
-
- forAll(oldFace, pointi)
- {
- newFace[pointi] = addedPoints[oldFace[pointi]];
- }
-
- bool flipFaceFlux = false;
-
- // Flip the face as necessary
- if
- (
- !mesh.isInternalFace(mf[facei])
- || mc[facei] == nei[mf[facei]]
- )
- {
- flipFaceFlux = true;
- }
-
- // Add the face
- ref.setAction
- (
- polyAddFace
- (
- newFace, // face
- mc[facei], // owner
- addedCells[facei], // neighbour
- -1, // master point
- -1, // master edge
- mf[facei], // master face for addition
- flipFaceFlux, // flux flip
- -1, // patch for face
- -1, // zone for face
- false // face zone flip
- )
- );
-
- if (debug > 1)
- {
- Pout<< "adding face: " << newFace
- << " own: " << mc[facei]
- << " nei: " << addedCells[facei]
- << endl;
- }
- }
-
- // Modify the faces from the master zone for the new neighbour
-
- const faceList& faces = mesh.faces();
-
- // Pout<< "mfFlip: " << mfFlip << endl;
-
- forAll(mf, facei)
- {
- const label curfaceID = mf[facei];
-
- // If the face is internal, modify its owner to be the newly
- // created cell. No flip is necessary
- if (!mesh.isInternalFace(curfaceID))
- {
- ref.setAction
- (
- polyModifyFace
- (
- faces[curfaceID], // modified face
- curfaceID, // label of face being modified
- addedCells[facei], // owner
- -1, // neighbour
- false, // face flip
- mesh.boundaryMesh().whichPatch(curfaceID),// patch for face
- false, // remove from zone
- faceZoneID_.index(), // zone for face
- mfFlip[facei] // face flip in zone
- )
- );
-
- if (debug > 1)
- {
- Pout<< "Modifying a boundary face. Face: " << curfaceID
- << " flip: " << mfFlip[facei]
- << endl;
- }
- }
-
- // If slave cell is owner, the face remains the same (but with
- // a new neighbour - the newly created cell). Otherwise, the
- // face is flipped.
- else if (sc[facei] == own[curfaceID])
- {
- // Orientation is good, just change neighbour
- ref.setAction
- (
- polyModifyFace
- (
- faces[curfaceID], // modified face
- curfaceID, // label of face being modified
- own[curfaceID], // owner
- addedCells[facei], // neighbour
- false, // face flip
- mesh.boundaryMesh().whichPatch(curfaceID),// patch for face
- false, // remove from zone
- faceZoneID_.index(), // zone for face
- mfFlip[facei] // face flip in zone
- )
- );
-
- if (debug > 1)
- {
- Pout<< "modify face, no flip " << curfaceID
- << " own: " << own[curfaceID]
- << " nei: " << addedCells[facei]
- << endl;
- }
- }
- else
- {
- // Change in orientation; flip face
- ref.setAction
- (
- polyModifyFace
- (
- faces[curfaceID].reverseFace(), // modified face
- curfaceID, // label of face being modified
- nei[curfaceID], // owner
- addedCells[facei], // neighbour
- true, // face flip
- mesh.boundaryMesh().whichPatch(curfaceID), // patch for face
- false, // remove from zone
- faceZoneID_.index(), // zone for face
- !mfFlip[facei] // face flip in zone
- )
- );
-
- if (debug > 1)
- {
- Pout<< "modify face, with flip " << curfaceID
- << " own: " << own[curfaceID]
- << " nei: " << addedCells[facei]
- << endl;
- }
- }
- }
-
- // Create new faces from edges
-
- const edgeList& zoneLocalEdges = masterFaceLayer.edges();
-
- const labelListList& edgeFaces = masterFaceLayer.edgeFaces();
-
- label nInternalEdges = masterFaceLayer.nInternalEdges();
-
- const labelList& meshEdges =
- mesh.faceZones()[faceZoneID_.index()].meshEdges();
-
- // Do all internal edges
- for (label curEdgeID = 0; curEdgeID < nInternalEdges; curEdgeID++)
- {
- face newFace(4);
-
- newFace[0] = mp[zoneLocalEdges[curEdgeID].start()];
- newFace[1] = mp[zoneLocalEdges[curEdgeID].end()];
- newFace[2] = addedPoints[zoneLocalEdges[curEdgeID].end()];
- newFace[3] = addedPoints[zoneLocalEdges[curEdgeID].start()];
-
- ref.setAction
- (
- polyAddFace
- (
- newFace, // face
- addedCells[edgeFaces[curEdgeID][0]], // owner
- addedCells[edgeFaces[curEdgeID][1]], // neighbour
- -1, // master point
- meshEdges[curEdgeID], // master edge
- -1, // master face
- false, // flip flux
- -1, // patch for face
- -1, // zone for face
- false // face zone flip
- )
- );
-
- if (debug > 1)
- {
- Pout<< "Add internal face off edge: " << newFace
- << " own: " << addedCells[edgeFaces[curEdgeID][0]]
- << " nei: " << addedCells[edgeFaces[curEdgeID][1]]
- << endl;
- }
- }
-
- // Prepare creation of faces from boundary edges.
- // Note:
- // A tricky part of the algorithm is finding the patch into which the
- // newly created face will be added. For this purpose, take the edge
- // and grab all the faces coming from it. From the set of faces
- // eliminate the internal faces and find the boundary face from the rest.
- // If there is more than one boundary face (excluding the ones in
- // the master zone), the patch with the lower label is selected.
-
- const labelListList& meshEdgeFaces = mesh.edgeFaces();
-
- const meshFaceZones& meshZones = mesh.faceZones();
-
- // Do all boundary edges
-
- for
- (
- label curEdgeID = nInternalEdges;
- curEdgeID < zoneLocalEdges.size();
- curEdgeID++
- )
- {
- face newFace(4);
- newFace[0] = mp[zoneLocalEdges[curEdgeID].start()];
- newFace[1] = mp[zoneLocalEdges[curEdgeID].end()];
- newFace[2] = addedPoints[zoneLocalEdges[curEdgeID].end()];
- newFace[3] = addedPoints[zoneLocalEdges[curEdgeID].start()];
-
- // Determine the patch for insertion
- const labelList& curFaces = meshEdgeFaces[meshEdges[curEdgeID]];
-
- label patchID = -1;
- label zoneID = -1;
-
- forAll(curFaces, facei)
- {
- const label cf = curFaces[facei];
-
- if (!mesh.isInternalFace(cf))
- {
- // Face not internal. Check if it is in the zone
- if (meshZones.whichZone(cf) != faceZoneID_.index())
- {
- // Found the face in a boundary patch which is not in zone
- patchID = mesh.boundaryMesh().whichPatch(cf);
- zoneID = mesh.faceZones().whichZone(cf);
-
- break;
- }
- }
- }
-
- if (patchID < 0)
- {
- FatalErrorInFunction
- << "Cannot find patch for edge " << meshEdges[curEdgeID]
- << ". Edge: " << mesh.edges()[meshEdges[curEdgeID]]
- << abort(FatalError);
- }
-
- ref.setAction
- (
- polyAddFace
- (
- newFace, // face
- addedCells[edgeFaces[curEdgeID][0]], // owner
- -1, // neighbour
- -1, // master point
- meshEdges[curEdgeID], // master edge
- -1, // master face
- false, // flip flux
- patchID, // patch for face
- zoneID, // zone
- false // zone face flip
- )
- );
-
- if (debug > 1)
- {
- Pout<< "add boundary face: " << newFace
- << " into patch " << patchID
- << " own: " << addedCells[edgeFaces[curEdgeID][0]]
- << endl;
- }
- }
-
- // Modify the remaining faces of the master cells to reconnect to the new
- // layer of faces.
- // Algorithm: Go through all the cells of the master zone and make
- // a map of faces to avoid duplicates. Do not insert the faces in
- // the master patch (as they have already been dealt with). Make
- // a master layer point renumbering map, which for every point in
- // the master layer gives its new label. Loop through all faces in
- // the map and attempt to renumber them using the master layer
- // point renumbering map. Once the face is renumbered, compare it
- // with the original face; if they are the same, the face has not
- // changed; if not, modify the face but keep all of its old
- // attributes (apart from the vertex numbers).
-
- // Create the map of faces in the master cell layer
- labelHashSet masterCellFaceMap(primitiveMesh::facesPerCell_*mc.size());
-
- const cellList& cells = mesh.cells();
-
- forAll(mc, celli)
- {
- const labelList& curFaces = cells[mc[celli]];
-
- forAll(curFaces, facei)
- {
- // Check if the face belongs to the master zone; if not add it
- if (meshZones.whichZone(curFaces[facei]) != faceZoneID_.index())
- {
- masterCellFaceMap.insert(curFaces[facei]);
- }
- }
- }
-
- // Create the master layer point map
- Map masterLayerPointMap(2*mp.size());
-
- forAll(mp, pointi)
- {
- masterLayerPointMap.insert
- (
- mp[pointi],
- addedPoints[pointi]
- );
- }
-
- // Grab the list of faces of the master layer
- const labelList masterCellFaces = masterCellFaceMap.toc();
-
- forAll(masterCellFaces, facei)
- {
- // Attempt to renumber the face using the masterLayerPointMap.
- // Missing point remain the same
-
- const label curFaceID = masterCellFaces[facei];
-
- const face& oldFace = faces[curFaceID];
-
- face newFace(oldFace.size());
-
- bool changed = false;
-
- forAll(oldFace, pointi)
- {
- if (masterLayerPointMap.found(oldFace[pointi]))
- {
- changed = true;
-
- newFace[pointi] = masterLayerPointMap.find(oldFace[pointi])();
- }
- else
- {
- newFace[pointi] = oldFace[pointi];
- }
- }
-
- // If the face has changed, create a modification entry
- if (changed)
- {
- // Get face zone and its flip
- label modifiedFaceZone = mesh.faceZones().whichZone(curFaceID);
- bool modifiedFaceZoneFlip = false;
-
- if (modifiedFaceZone >= 0)
- {
- modifiedFaceZoneFlip =
- mesh.faceZones()[modifiedFaceZone].flipMap()
- [
- mesh.faceZones()[modifiedFaceZone].whichFace(curFaceID)
- ];
- }
-
- if (mesh.isInternalFace(curFaceID))
- {
- ref.setAction
- (
- polyModifyFace
- (
- newFace, // modified face
- curFaceID, // label of face being modified
- own[curFaceID], // owner
- nei[curFaceID], // neighbour
- false, // face flip
- -1, // patch for face
- false, // remove from zone
- modifiedFaceZone, // zone for face
- modifiedFaceZoneFlip // face flip in zone
- )
- );
-
- if (debug > 1)
- {
- Pout<< "modifying stick-out face. Internal Old face: "
- << oldFace
- << " new face: " << newFace
- << " own: " << own[curFaceID]
- << " nei: " << nei[curFaceID]
- << endl;
- }
- }
- else
- {
- ref.setAction
- (
- polyModifyFace
- (
- newFace, // modified face
- curFaceID, // label of face being modified
- own[curFaceID], // owner
- -1, // neighbour
- false, // face flip
- mesh.boundaryMesh().whichPatch(curFaceID),
- // patch for face
- false, // remove from zone
- modifiedFaceZone, // zone for face
- modifiedFaceZoneFlip // face flip in zone
- )
- );
-
- if (debug > 1)
- {
- Pout<< "modifying stick-out face. Boundary Old face: "
- << oldFace
- << " new face: " << newFace
- << " own: " << own[curFaceID]
- << " patch: "
- << mesh.boundaryMesh().whichPatch(curFaceID)
- << endl;
- }
- }
- }
- }
-
- if (debug)
- {
- Pout<< "void layerAdditionRemoval::addCellLayer(polyTopoChange&) const "
- << " for object " << name() << ": "
- << "Finished adding cell layer" << endl;
- }
-}
-
-
-// ************************************************************************* //
diff --git a/src/dynamicMesh/polyTopoChanger/layerAdditionRemoval/layerAdditionRemoval.C b/src/dynamicMesh/polyTopoChanger/layerAdditionRemoval/layerAdditionRemoval.C
deleted file mode 100644
index 9860ea64f9..0000000000
--- a/src/dynamicMesh/polyTopoChanger/layerAdditionRemoval/layerAdditionRemoval.C
+++ /dev/null
@@ -1,518 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration | Website: https://openfoam.org
- \\ / A nd | Copyright (C) 2011-2023 OpenFOAM Foundation
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see ());
-
- if (nFaces == 0)
- {
- FatalErrorInFunction
- << "Face extrusion zone contains no faces. "
- << "Please check your mesh definition."
- << abort(FatalError);
- }
-
- if (debug)
- {
- Pout<< "Cell layer addition/removal object " << name() << " :" << nl
- << " faceZoneID: " << faceZoneID_ << endl;
- }
-}
-
-
-Foam::scalar Foam::layerAdditionRemoval::readOldThickness
-(
- const dictionary& dict
-)
-{
- return dict.lookupOrDefault("oldLayerThickness", -1.0);
-}
-
-
-void Foam::layerAdditionRemoval::clearAddressing() const
-{
- if (pointsPairingPtr_)
- {
- if (debug)
- {
- Pout<< "layerAdditionRemoval::clearAddressing()" << nl
- << " clearing pointsPairingPtr_" << endl;
- }
-
- deleteDemandDrivenData(pointsPairingPtr_);
- }
-
- if (facesPairingPtr_)
- {
- if (debug)
- {
- Pout<< "layerAdditionRemoval::clearAddressing()" << nl
- << " clearing facesPairingPtr_" << endl;
- }
-
- deleteDemandDrivenData(facesPairingPtr_);
- }
-}
-
-
-// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
-
-Foam::layerAdditionRemoval::layerAdditionRemoval
-(
- const word& name,
- const label index,
- const polyTopoChanger& ptc,
- const word& zoneName,
- const scalar minThickness,
- const scalar maxThickness,
- const Switch thicknessFromVolume
-)
-:
- polyMeshModifier(name, index, ptc, true),
- faceZoneID_(zoneName, ptc.mesh().faceZones()),
- minLayerThickness_(minThickness),
- maxLayerThickness_(maxThickness),
- thicknessFromVolume_(thicknessFromVolume),
- oldLayerThickness_(-1.0),
- pointsPairingPtr_(nullptr),
- facesPairingPtr_(nullptr),
- triggerRemoval_(-1),
- triggerAddition_(-1)
-{
- checkDefinition();
-}
-
-
-Foam::layerAdditionRemoval::layerAdditionRemoval
-(
- const word& name,
- const dictionary& dict,
- const label index,
- const polyTopoChanger& ptc
-)
-:
- polyMeshModifier(name, index, ptc, Switch(dict.lookup("active"))),
- faceZoneID_(dict.lookup("faceZoneName"), ptc.mesh().faceZones()),
- minLayerThickness_(dict.lookup("minLayerThickness")),
- maxLayerThickness_(dict.lookup("maxLayerThickness")),
- thicknessFromVolume_
- (
- dict.lookupOrDefault("thicknessFromVolume", true)
- ),
- oldLayerThickness_(readOldThickness(dict)),
- pointsPairingPtr_(nullptr),
- facesPairingPtr_(nullptr),
- triggerRemoval_(-1),
- triggerAddition_(-1)
-{
- checkDefinition();
-}
-
-
-// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
-
-Foam::layerAdditionRemoval::~layerAdditionRemoval()
-{
- clearAddressing();
-}
-
-
-// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
-
-bool Foam::layerAdditionRemoval::changeTopology() const
-{
- // Protect from multiple calculation in the same time-step
- if (triggerRemoval_ > -1 || triggerAddition_ > -1)
- {
- return true;
- }
-
- // Go through all the cells in the master layer and calculate
- // approximate layer thickness as the ratio of the cell volume and
- // face area in the face zone.
- // Layer addition:
- // When the max thickness exceeds the threshold, trigger refinement.
- // Layer removal:
- // When the min thickness falls below the threshold, trigger removal.
-
- const polyMesh& mesh = topoChanger().mesh();
-
- const faceZone& fz = mesh.faceZones()[faceZoneID_.index()];
- const labelList& mc = fz.masterCells();
-
- const scalarField& V = mesh.cellVolumes();
- const vectorField& S = mesh.faceAreas();
-
- if (min(V) < -vSmall)
- {
- FatalErrorInFunction
- << "negative cell volume. Error in mesh motion before "
- << "topological change.\n V: " << V
- << abort(FatalError);
- }
-
- scalar avgDelta = 0;
- scalar minDelta = great;
- scalar maxDelta = 0;
- label nDelta = 0;
-
- if (thicknessFromVolume_)
- {
- // Thickness calculated from cell volume/face area
- forAll(fz, facei)
- {
- scalar curDelta = V[mc[facei]]/mag(S[fz[facei]]);
- avgDelta += curDelta;
- minDelta = min(minDelta, curDelta);
- maxDelta = max(maxDelta, curDelta);
- }
-
- nDelta = fz.size();
- }
- else
- {
- // Thickness calculated from edges on layer
- const Map& meshZonePointMap = fz().meshPointMap();
-
- // Edges with only one point on zone
- forAll(mc, facei)
- {
- const cell& cFaces = mesh.cells()[mc[facei]];
- const edgeList cellEdges(cFaces.edges(mesh.faces()));
-
- forAll(cellEdges, i)
- {
- const edge& e = cellEdges[i];
-
- if (meshZonePointMap.found(e[0]))
- {
- if (!meshZonePointMap.found(e[1]))
- {
- scalar curDelta = e.mag(mesh.points());
- avgDelta += curDelta;
- nDelta++;
- minDelta = min(minDelta, curDelta);
- maxDelta = max(maxDelta, curDelta);
- }
- }
- else
- {
- if (meshZonePointMap.found(e[1]))
- {
- scalar curDelta = e.mag(mesh.points());
- avgDelta += curDelta;
- nDelta++;
- minDelta = min(minDelta, curDelta);
- maxDelta = max(maxDelta, curDelta);
- }
- }
- }
- }
- }
-
- reduce(minDelta, minOp());
- reduce(maxDelta, maxOp());
- reduce(avgDelta, sumOp());
- reduce(nDelta, sumOp());
-
- avgDelta /= nDelta;
-
- if (debug)
- {
- Pout<< "bool layerAdditionRemoval::changeTopology() const "
- << " for object " << name() << " : " << nl
- << "Layer thickness: min: " << minDelta
- << " max: " << maxDelta << " avg: " << avgDelta
- << " old thickness: " << oldLayerThickness_ << nl
- << "Removal threshold: " << minLayerThickness_
- << " addition threshold: " << maxLayerThickness_ << endl;
- }
-
- bool topologicalChange = false;
-
- // If the thickness is decreasing and crosses the min thickness,
- // trigger removal
- if (oldLayerThickness_ < 0)
- {
- if (debug)
- {
- Pout<< "First step. No addition/removal" << endl;
- }
-
- // No topological changes allowed before first mesh motion
- oldLayerThickness_ = avgDelta;
-
- topologicalChange = false;
- }
- else if (avgDelta < oldLayerThickness_)
- {
- // Layers moving towards removal
- if (minDelta < minLayerThickness_)
- {
- // Check layer pairing
- if (setLayerPairing())
- {
- // A mesh layer detected. Check that collapse is valid
- if (validCollapse())
- {
- // At this point, info about moving the old mesh
- // in a way to collapse the cells in the removed
- // layer is available. Not sure what to do with
- // it.
-
- if (debug)
- {
- Pout<< "bool layerAdditionRemoval::changeTopology() "
- << " const for object " << name() << " : "
- << "Triggering layer removal" << endl;
- }
-
- triggerRemoval_ = mesh.time().timeIndex();
-
- // Old thickness looses meaning.
- // Set it up to indicate layer removal
- oldLayerThickness_ = great;
-
- topologicalChange = true;
- }
- else
- {
- // No removal, clear addressing
- clearAddressing();
- }
- }
- }
- else
- {
- oldLayerThickness_ = avgDelta;
- }
- }
- else
- {
- // Layers moving towards addition
- if (maxDelta > maxLayerThickness_)
- {
- if (debug)
- {
- Pout<< "bool layerAdditionRemoval::changeTopology() const "
- << " for object " << name() << " : "
- << "Triggering layer addition" << endl;
- }
-
- triggerAddition_ = mesh.time().timeIndex();
-
- // Old thickness looses meaning.
- // Set it up to indicate layer removal
- oldLayerThickness_ = 0;
-
- topologicalChange = true;
- }
- else
- {
- oldLayerThickness_ = avgDelta;
- }
- }
-
- return topologicalChange;
-}
-
-
-void Foam::layerAdditionRemoval::setRefinement(polyTopoChange& ref) const
-{
- // Insert the layer addition/removal instructions
- // into the topological change
-
- if (triggerRemoval_ == topoChanger().mesh().time().timeIndex())
- {
- removeCellLayer(ref);
-
- // Clear addressing. This also resets the addition/removal data
- if (debug)
- {
- Pout<< "layerAdditionRemoval::setRefinement(polyTopoChange&) "
- << "for object " << name() << " : "
- << "Clearing addressing after layer removal" << endl;
- }
-
- triggerRemoval_ = -1;
- clearAddressing();
- }
-
- if (triggerAddition_ == topoChanger().mesh().time().timeIndex())
- {
- addCellLayer(ref);
-
- // Clear addressing. This also resets the addition/removal data
- if (debug)
- {
- Pout<< "layerAdditionRemoval::setRefinement(polyTopoChange&) "
- << "for object " << name() << " : "
- << "Clearing addressing after layer addition" << endl;
- }
-
- triggerAddition_ = -1;
- clearAddressing();
- }
-}
-
-
-void Foam::layerAdditionRemoval::topoChange(const polyTopoChangeMap&)
-{
- if (debug)
- {
- Pout<< "layerAdditionRemoval::topoChange(const polyTopoChangeMap&) "
- << "for object " << name() << " : "
- << "Clearing addressing on external request";
-
- if (pointsPairingPtr_ || facesPairingPtr_)
- {
- Pout<< "Pointers set." << endl;
- }
- else
- {
- Pout<< "Pointers not set." << endl;
- }
- }
-
- // Mesh has changed topologically. Update local topological data
- faceZoneID_.update(topoChanger().mesh().faceZones());
-
- clearAddressing();
-}
-
-
-void Foam::layerAdditionRemoval::setMinLayerThickness(const scalar t) const
-{
- if (t < vSmall || maxLayerThickness_ < t)
- {
- FatalErrorInFunction
- << "Incorrect layer thickness definition."
- << abort(FatalError);
- }
-
- minLayerThickness_ = t;
-}
-
-
-void Foam::layerAdditionRemoval::setMaxLayerThickness(const scalar t) const
-{
- if (t < minLayerThickness_)
- {
- FatalErrorInFunction
- << "Incorrect layer thickness definition."
- << abort(FatalError);
- }
-
- maxLayerThickness_ = t;
-}
-
-
-void Foam::layerAdditionRemoval::write(Ostream& os) const
-{
- os << nl << type() << nl
- << name()<< nl
- << faceZoneID_ << nl
- << minLayerThickness_ << nl
- << oldLayerThickness_ << nl
- << maxLayerThickness_ << nl
- << thicknessFromVolume_ << endl;
-}
-
-
-void Foam::layerAdditionRemoval::writeDict(Ostream& os) const
-{
- os << nl << name() << nl << token::BEGIN_BLOCK << nl
- << " type " << type()
- << token::END_STATEMENT << nl
- << " faceZoneName " << faceZoneID_.name()
- << token::END_STATEMENT << nl
- << " minLayerThickness " << minLayerThickness_
- << token::END_STATEMENT << nl
- << " maxLayerThickness " << maxLayerThickness_
- << token::END_STATEMENT << nl
- << " thicknessFromVolume " << thicknessFromVolume_
- << token::END_STATEMENT << nl
- << " oldLayerThickness " << oldLayerThickness_
- << token::END_STATEMENT << nl
- << " active " << active()
- << token::END_STATEMENT << nl
- << token::END_BLOCK << endl;
-}
-
-
-// ************************************************************************* //
diff --git a/src/dynamicMesh/polyTopoChanger/layerAdditionRemoval/layerAdditionRemoval.H b/src/dynamicMesh/polyTopoChanger/layerAdditionRemoval/layerAdditionRemoval.H
deleted file mode 100644
index 95bc3f5dcf..0000000000
--- a/src/dynamicMesh/polyTopoChanger/layerAdditionRemoval/layerAdditionRemoval.H
+++ /dev/null
@@ -1,234 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration | Website: https://openfoam.org
- \\ / A nd | Copyright (C) 2011-2023 OpenFOAM Foundation
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see extrusionDir() const;
-
- //- Add a layer of cells
- void addCellLayer(polyTopoChange&) const;
-
- //- Remove a layer of cells
- void removeCellLayer(polyTopoChange&) const;
-
- //- Clear addressing
- void clearAddressing() const;
-
-
- // Helpers
-
- //- Optionally read old thickness
- static scalar readOldThickness(const dictionary&);
-
-
- // Static Data Members
-
- //- Thickness insertion fraction for the pre-motion
- static const scalar addDelta_;
-
- //- Thickness removal fraction for the cell collapse
- // Note: the cell will be collapsed to this relative
- // thickness before the layer is removed.
- static const scalar removeDelta_;
-
-public:
-
- //- Runtime type information
- TypeName("layerAdditionRemoval");
-
-
- // Constructors
-
- //- Construct from components
- layerAdditionRemoval
- (
- const word& name,
- const label index,
- const polyTopoChanger& ptc,
- const word& zoneName,
- const scalar minThickness,
- const scalar maxThickness,
- const Switch thicknessFromVolume = true
- );
-
- //- Construct from dictionary
- layerAdditionRemoval
- (
- const word& name,
- const dictionary& dict,
- const label index,
- const polyTopoChanger& ptc
- );
-
- //- Disallow default bitwise copy construction
- layerAdditionRemoval(const layerAdditionRemoval&) = delete;
-
-
- //- Destructor
- virtual ~layerAdditionRemoval();
-
-
- // Member Functions
-
- //- Check for topology change
- virtual bool changeTopology() const;
-
- //- Insert the layer addition/removal instructions
- // into the topological change
- virtual void setRefinement(polyTopoChange&) const;
-
- //- Modify motion points to comply with the topological change
- virtual void modifyMotionPoints(pointField& motionPoints) const;
-
- //- Force recalculation of locally stored data on topological change
- virtual void topoChange(const polyTopoChangeMap&);
-
-
- // Edit
-
- //- Return min layer thickness which triggers removal
- scalar minLayerThickness() const
- {
- return minLayerThickness_;
- }
-
- //- Set min layer thickness which triggers removal
- void setMinLayerThickness(const scalar t) const;
-
- //- Return max layer thickness which triggers removal
- scalar maxLayerThickness() const
- {
- return maxLayerThickness_;
- }
-
- //- Set max layer thickness which triggers removal
- void setMaxLayerThickness(const scalar t) const;
-
-
- //- Write
- virtual void write(Ostream&) const;
-
- //- Write dictionary
- virtual void writeDict(Ostream&) const;
-
-
- // Member Operators
-
- //- Disallow default bitwise assignment
- void operator=(const layerAdditionRemoval&) = delete;
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //
diff --git a/src/dynamicMesh/polyTopoChanger/layerAdditionRemoval/removeCellLayer.C b/src/dynamicMesh/polyTopoChanger/layerAdditionRemoval/removeCellLayer.C
deleted file mode 100644
index 4970dabd0c..0000000000
--- a/src/dynamicMesh/polyTopoChanger/layerAdditionRemoval/removeCellLayer.C
+++ /dev/null
@@ -1,420 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration | Website: https://openfoam.org
- \\ / A nd | Copyright (C) 2011-2023 OpenFOAM Foundation
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see ()) > 0)
- {
- return false;
- }
- else
- {
- return true;
- }
-}
-
-
-void Foam::layerAdditionRemoval::removeCellLayer
-(
- polyTopoChange& ref
-) const
-{
- // Algorithm for layer removal. Second phase: topological change
- // 2) Go through all the faces of the master cell layer and remove
- // the ones that are not in the master face zone.
- //
- // 3) Grab all the faces coming out of points that are collapsed
- // and select the ones that are not marked for removal. Go
- // through the remaining faces and replace the point to remove by
- // the equivalent point in the master face zone.
- if (debug)
- {
- Pout<< "Removing the cell layer for object " << name() << endl;
- }
-
- const polyMesh& mesh = topoChanger().mesh();
-
- const labelList& own = mesh.faceOwner();
- const labelList& nei = mesh.faceNeighbour();
-
- const labelList& ptc = pointsPairing();
- const labelList& ftc = facesPairing();
-
- // Remove all the cells from the master layer
- const labelList& mc =
- mesh.faceZones()[faceZoneID_.index()].masterCells();
-
- forAll(mc, facei)
- {
- label slaveSideCell = own[ftc[facei]];
-
- if (mesh.isInternalFace(ftc[facei]) && slaveSideCell == mc[facei])
- {
- // Owner cell of the face is being removed.
- // Grab the neighbour instead
- slaveSideCell = nei[ftc[facei]];
- }
-
- ref.setAction(polyRemoveCell(mc[facei], slaveSideCell));
- }
-
- // Remove all the faces from the master layer cells which are not in
- // the master face layer
- labelHashSet facesToRemoveMap(mc.size()*primitiveMesh::facesPerCell_);
-
- const cellList& cells = mesh.cells();
-
- forAll(mc, celli)
- {
- const cell& curCell = cells[mc[celli]];
-
- forAll(curCell, facei)
- {
- // Check if the face is in the master zone. If not, remove it
- if
- (
- mesh.faceZones().whichZone(curCell[facei])
- != faceZoneID_.index()
- )
- {
- facesToRemoveMap.insert(curCell[facei]);
- }
- }
- }
-
- forAllConstIter(labelHashSet, facesToRemoveMap, iter)
- {
- ref.setAction(polyRemoveFace(iter.key()));
- }
-
- // Remove all points that will be collapsed
- forAll(ptc, pointi)
- {
- ref.setAction(polyRemovePoint(ptc[pointi]));
- }
-
- // Grab all faces coming off points to be deleted. If the face
- // has not been deleted, replace the removed point with the
- // equivalent from the master layer.
-
- // Make a map of all point to be removed, giving the master point label
- // for each of them
-
- Map removedPointMap(2*ptc.size());
-
- const labelList& meshPoints =
- mesh.faceZones()[faceZoneID_.index()]().meshPoints();
-
- forAll(ptc, pointi)
- {
- removedPointMap.insert(ptc[pointi], meshPoints[pointi]);
- }
-
- const labelListList& pf = mesh.pointFaces();
-
- const faceList& faces = mesh.faces();
-
- // Make a list of faces to be modified using the map to avoid duplicates
- labelHashSet facesToModify(ptc.size()*primitiveMesh::facesPerPoint_);
-
- forAll(ptc, pointi)
- {
- const labelList& curFaces = pf[ptc[pointi]];
-
- forAll(curFaces, facei)
- {
- if (!facesToRemoveMap.found(curFaces[facei]))
- {
- facesToModify.insert(curFaces[facei]);
- }
- }
- }
-
- labelList ftm = facesToModify.toc();
-
- if (debug > 1)
- {
- Pout<< "faces to modify: " << ftm << endl;
- }
-
- forAll(ftm, facei)
- {
- // For every face to modify, copy the face and re-map the vertices.
- // It is known all the faces will be changed since they hang off
- // re-mapped vertices
- label curFaceID = ftm[facei];
-
- face newFace(faces[curFaceID]);
-
- forAll(newFace, pointi)
- {
- Map::iterator rpmIter =
- removedPointMap.find(newFace[pointi]);
-
- if (rpmIter != removedPointMap.end())
- {
- // Point mapped. Replace it
- newFace[pointi] = rpmIter();
- }
- }
-
- if (debug > 1)
- {
- Pout<< "face label: " << curFaceID
- << " old face: " << faces[curFaceID]
- << " new face: " << newFace << endl;
- }
-
- // Get face zone and its flip
- label modifiedFaceZone = mesh.faceZones().whichZone(curFaceID);
- bool modifiedFaceZoneFlip = false;
-
- if (modifiedFaceZone >= 0)
- {
- const faceZone& fz = mesh.faceZones()[modifiedFaceZone];
- modifiedFaceZoneFlip = fz.flipMap()[fz.whichFace(curFaceID)];
- }
-
- label newNeighbour = -1;
-
- if (curFaceID < mesh.nInternalFaces())
- {
- newNeighbour = nei[curFaceID];
- }
-
- // Modify the face
- ref.setAction
- (
- polyModifyFace
- (
- newFace, // modified face
- curFaceID, // label of face being modified
- own[curFaceID], // owner
- newNeighbour, // neighbour
- false, // face flip
- mesh.boundaryMesh().whichPatch(curFaceID),// patch for face
- false, // remove from zone
- modifiedFaceZone, // zone for face
- modifiedFaceZoneFlip // face flip in zone
- )
- );
- }
-
- // Modify the faces in the master layer to point past the removed cells
-
- const labelList& mf = mesh.faceZones()[faceZoneID_.index()];
- const boolList& mfFlip = mesh.faceZones()[faceZoneID_.index()].flipMap();
-
- forAll(mf, facei)
- {
- // Grab the owner and neighbour of the faces to be collapsed and get rid
- // of the cell to be removed
- label masterSideCell = own[mf[facei]];
-
- if (masterSideCell == mc[facei])
- {
- if (mesh.isInternalFace(mf[facei]))
- {
- // Owner cell of the face is being removed.
- // Grab the neighbour instead
- masterSideCell = nei[mf[facei]];
- }
- else
- {
- masterSideCell = -1;
- }
- }
-
- label slaveSideCell = own[ftc[facei]];
-
- if (slaveSideCell == mc[facei])
- {
- if (mesh.isInternalFace(ftc[facei]))
- {
- // Owner cell of the face is being removed.
- // Grab the neighbour instead
- slaveSideCell = nei[ftc[facei]];
- }
- else
- {
- slaveSideCell = -1;
- }
- }
-
- // Find out if the face needs to be flipped
- label newOwner = -1;
- label newNeighbour = -1;
- bool flipFace = false;
- label newPatchID = -1;
- label newZoneID = -1;
-
- // Is any of the faces a boundary face? If so, grab the patch
- // A boundary-to-boundary collapse is checked for in validCollapse()
- // and cannot happen here.
-
- if (!mesh.isInternalFace(mf[facei]))
- {
- // Master is the boundary face: it gets a new owner but no flip
- newOwner = slaveSideCell;
- newNeighbour = -1;
- flipFace = false;
- newPatchID = mesh.boundaryMesh().whichPatch(mf[facei]);
- newZoneID = mesh.faceZones().whichZone(mf[facei]);
- }
- else if (!mesh.isInternalFace(ftc[facei]))
- {
- // Slave is the boundary face: grab its patch
- newOwner = slaveSideCell;
- newNeighbour = -1;
-
- // Find out if the face flip is necessary
- if (own[mf[facei]] == slaveSideCell)
- {
- flipFace = false;
- }
- else
- {
- flipFace = true;
- }
-
- newPatchID = mesh.boundaryMesh().whichPatch(ftc[facei]);
-
- // The zone of the master face is preserved
- newZoneID = mesh.faceZones().whichZone(mf[facei]);
- }
- else
- {
- // Both faces are internal: flip is decided based on which of the
- // new cells around it has a lower label
- newOwner = min(masterSideCell, slaveSideCell);
- newNeighbour = max(masterSideCell, slaveSideCell);
-
- if (newOwner == own[mf[facei]] || newNeighbour == nei[mf[facei]])
- {
- flipFace = false;
- }
- else
- {
- flipFace = true;
- }
-
- newPatchID = -1;
-
- // The zone of the master face is preserved
- newZoneID = mesh.faceZones().whichZone(mf[facei]);
- }
-
- // Modify the face and flip if necessary
- face newFace = faces[mf[facei]];
- bool zoneFlip = mfFlip[facei];
-
- if (flipFace)
- {
- newFace.flip();
- zoneFlip = !zoneFlip;
- }
-
- if (debug > 1)
- {
- Pout<< "Modifying face " << mf[facei]
- << " newFace: " << newFace << nl
- << " newOwner: " << newOwner
- << " newNeighbour: " << newNeighbour
- << " flipFace: " << flipFace
- << " newPatchID: " << newPatchID
- << " newZoneID: " << newZoneID << nl
- << " oldOwn: " << own[mf[facei]]
- << " oldNei: " << nei[mf[facei]] << endl;
- }
-
- ref.setAction
- (
- polyModifyFace
- (
- newFace, // modified face
- mf[facei], // label of face being modified
- newOwner, // owner
- newNeighbour, // neighbour
- flipFace, // flip
- newPatchID, // patch for face
- false, // remove from zone
- newZoneID, // new zone
- zoneFlip // face flip in zone
- )
- );
- }
-}
-
-
-// ************************************************************************* //
diff --git a/src/dynamicMesh/polyTopoChanger/layerAdditionRemoval/setLayerPairing.C b/src/dynamicMesh/polyTopoChanger/layerAdditionRemoval/setLayerPairing.C
deleted file mode 100644
index ced1db7dcc..0000000000
--- a/src/dynamicMesh/polyTopoChanger/layerAdditionRemoval/setLayerPairing.C
+++ /dev/null
@@ -1,232 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration | Website: https://openfoam.org
- \\ / A nd | Copyright (C) 2011-2023 OpenFOAM Foundation
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see ());
- reduce(nFaceErrors, sumOp());
-
- if (nPointErrors > 0 || nFaceErrors > 0)
- {
- clearAddressing();
-
- return false;
- }
- else
- {
- // Valid layer
- return true;
- }
-}
-
-
-const Foam::labelList& Foam::layerAdditionRemoval::pointsPairing() const
-{
- if (!pointsPairingPtr_)
- {
- FatalErrorInFunction
- << "Problem with layer pairing data for object " << name()
- << abort(FatalError);
- }
-
- return *pointsPairingPtr_;
-}
-
-const Foam::labelList& Foam::layerAdditionRemoval::facesPairing() const
-{
- if (!facesPairingPtr_)
- {
- FatalErrorInFunction
- << "Problem with layer pairing data for object " << name()
- << abort(FatalError);
- }
-
- return *facesPairingPtr_;
-}
-
-
-// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
-
-void Foam::layerAdditionRemoval::modifyMotionPoints
-(
- pointField& motionPoints
-) const
-{
- if (debug)
- {
- Pout<< "void layerAdditionRemoval::modifyMotionPoints("
- << "pointField& motionPoints) const for object "
- << name() << " : ";
- }
-
- if (debug)
- {
- Pout<< "No motion point adjustment" << endl;
- }
-}
-
-
-// ************************************************************************* //