From 950e6672597e63148499bc79b877925e4e133a2d Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Mon, 13 Jul 2020 15:29:12 +0100
Subject: [PATCH] ENH: snappyHexMesh: optionally remove 'small' regions. Fixes
#1772.
This is for a very specific use case where the faceZones are
imprinted after meshing the normal geometry. This sometimes
splits off badly connected bits of the mesh. One way to remove
these is to use e.g. subsetMesh. This embeds the
same functionality inside snappyHexMesh.
---
.../meshRefinement/meshRefinement.H | 28 ++--
.../meshRefinement/meshRefinementBaffles.C | 113 +++++++------
.../refinementParameters.C | 1 +
.../refinementParameters.H | 11 +-
.../snappyHexMeshDriver/snappyRefineDriver.C | 151 +++++++++++++++++-
.../snappyHexMeshDriver/snappyRefineDriver.H | 6 +-
6 files changed, 252 insertions(+), 58 deletions(-)
diff --git a/src/mesh/snappyHexMesh/meshRefinement/meshRefinement.H b/src/mesh/snappyHexMesh/meshRefinement/meshRefinement.H
index 8b0fbaa8567..0fdf183b78b 100644
--- a/src/mesh/snappyHexMesh/meshRefinement/meshRefinement.H
+++ b/src/mesh/snappyHexMesh/meshRefinement/meshRefinement.H
@@ -244,15 +244,6 @@ private:
labelList& minLevel
) const;
- //- Remove cells. Put exposedFaces into exposedPatchIDs.
- autoPtr<mapPolyMesh> doRemoveCells
- (
- const labelList& cellsToRemove,
- const labelList& exposedFaces,
- const labelList& exposedPatchIDs,
- removeCells& cellRemover
- );
-
//- Get cells which are inside any closed surface. Note that
// all closed surfaces
// will have already been oriented to have keepPoint outside.
@@ -1336,6 +1327,16 @@ public:
const List<labelPair>& baffles
);
+ //- Get per-face information (faceZone, master/slave patch)
+ void getZoneFaces
+ (
+ const labelList& zoneIDs,
+ labelList& faceZoneID,
+ labelList& ownPatch,
+ labelList& neiPatch,
+ labelList& nBaffles
+ ) const;
+
//- Create baffles for faces on faceZones. Return created baffles
// (= pairs of faces) and corresponding faceZone
autoPtr<mapPolyMesh> createZoneBaffles
@@ -1458,6 +1459,15 @@ public:
const writer<scalar>& leakPathFormatter
);
+ //- Remove cells. Put exposedFaces into exposedPatchIDs.
+ autoPtr<mapPolyMesh> doRemoveCells
+ (
+ const labelList& cellsToRemove,
+ const labelList& exposedFaces,
+ const labelList& exposedPatchIDs,
+ removeCells& cellRemover
+ );
+
//- Split faces into two
void doSplitFaces
(
diff --git a/src/mesh/snappyHexMesh/meshRefinement/meshRefinementBaffles.C b/src/mesh/snappyHexMesh/meshRefinement/meshRefinementBaffles.C
index eaa2e901691..f6db6804bfd 100644
--- a/src/mesh/snappyHexMesh/meshRefinement/meshRefinementBaffles.C
+++ b/src/mesh/snappyHexMesh/meshRefinement/meshRefinementBaffles.C
@@ -715,6 +715,67 @@ Foam::List<Foam::labelPair> Foam::meshRefinement::subsetBaffles
}
+void Foam::meshRefinement::getZoneFaces
+(
+ const labelList& zoneIDs,
+ labelList& faceZoneID,
+ labelList& ownPatch,
+ labelList& neiPatch,
+ labelList& nBaffles
+) const
+{
+ const faceZoneMesh& faceZones = mesh_.faceZones();
+
+ // Per (internal) face the patch related to the faceZone
+ ownPatch.setSize(mesh_.nFaces());
+ ownPatch= -1;
+ neiPatch.setSize(mesh_.nFaces());
+ neiPatch = -1;
+ faceZoneID.setSize(mesh_.nFaces());
+ faceZoneID = -1;
+ nBaffles.setSize(zoneIDs.size());
+ nBaffles = Zero;
+
+ forAll(zoneIDs, j)
+ {
+ label zoneI = zoneIDs[j];
+ const faceZone& fz = faceZones[zoneI];
+ const word& masterName = faceZoneToMasterPatch_[fz.name()];
+ label masterPatchI = mesh_.boundaryMesh().findPatchID(masterName);
+ const word& slaveName = faceZoneToSlavePatch_[fz.name()];
+ label slavePatchI = mesh_.boundaryMesh().findPatchID(slaveName);
+
+ if (masterPatchI == -1 || slavePatchI == -1)
+ {
+ FatalErrorInFunction
+ << "Problem: masterPatchI:" << masterPatchI
+ << " slavePatchI:" << slavePatchI << exit(FatalError);
+ }
+
+ forAll(fz, i)
+ {
+ label faceI = fz[i];
+ if (mesh_.isInternalFace(faceI))
+ {
+ if (fz.flipMap()[i])
+ {
+ ownPatch[faceI] = slavePatchI;
+ neiPatch[faceI] = masterPatchI;
+ }
+ else
+ {
+ ownPatch[faceI] = masterPatchI;
+ neiPatch[faceI] = slavePatchI;
+ }
+ faceZoneID[faceI] = zoneI;
+
+ nBaffles[j]++;
+ }
+ }
+ }
+}
+
+
Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::createZoneBaffles
(
const labelList& zoneIDs,
@@ -731,52 +792,12 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::createZoneBaffles
// Split internal faces on interface surfaces
Info<< "Converting zoned faces into baffles ..." << endl;
- // Per (internal) face the patch it should go into
- labelList ownPatch(mesh_.nFaces(), -1);
- labelList neiPatch(mesh_.nFaces(), -1);
- labelList faceZoneID(mesh_.nFaces(), -1);
-
- labelList nBaffles(zoneIDs.size(), Zero);
-
- forAll(zoneIDs, j)
- {
- label zoneI = zoneIDs[j];
-
- const faceZone& fz = faceZones[zoneI];
-
- const word& masterName = faceZoneToMasterPatch_[fz.name()];
- label masterPatchI = mesh_.boundaryMesh().findPatchID(masterName);
- const word& slaveName = faceZoneToSlavePatch_[fz.name()];
- label slavePatchI = mesh_.boundaryMesh().findPatchID(slaveName);
-
- if (masterPatchI == -1 || slavePatchI == -1)
- {
- FatalErrorInFunction
- << "Problem: masterPatchI:" << masterPatchI
- << " slavePatchI:" << slavePatchI << exit(FatalError);
- }
-
- forAll(fz, i)
- {
- label faceI = fz[i];
- if (mesh_.isInternalFace(faceI))
- {
- if (fz.flipMap()[i])
- {
- ownPatch[faceI] = slavePatchI;
- neiPatch[faceI] = masterPatchI;
- }
- else
- {
- ownPatch[faceI] = masterPatchI;
- neiPatch[faceI] = slavePatchI;
- }
- faceZoneID[faceI] = zoneI;
-
- nBaffles[j]++;
- }
- }
- }
+ // Get faceZone and patch(es) per face (or -1 if face not on faceZone)
+ labelList faceZoneID;
+ labelList ownPatch;
+ labelList neiPatch;
+ labelList nBaffles;
+ getZoneFaces(zoneIDs, faceZoneID, ownPatch, neiPatch, nBaffles);
label nLocalBaffles = sum(nBaffles);
diff --git a/src/mesh/snappyHexMesh/snappyHexMeshDriver/refinementParameters/refinementParameters.C b/src/mesh/snappyHexMesh/snappyHexMeshDriver/refinementParameters/refinementParameters.C
index 47c2bae13ad..f774db14ce4 100644
--- a/src/mesh/snappyHexMesh/snappyHexMeshDriver/refinementParameters/refinementParameters.C
+++ b/src/mesh/snappyHexMesh/snappyHexMeshDriver/refinementParameters/refinementParameters.C
@@ -94,6 +94,7 @@ Foam::refinementParameters::refinementParameters
),
nErodeCellZone_(dict.getOrDefault<label>("nCellZoneErodeIter", 0)),
nFilterIter_(dict.getOrDefault<label>("nFilterIter", 2)),
+ minCellFraction_(dict.getOrDefault<scalar>("minCellFraction", 0)),
dryRun_(dryRun)
{
point locationInMesh;
diff --git a/src/mesh/snappyHexMesh/snappyHexMeshDriver/refinementParameters/refinementParameters.H b/src/mesh/snappyHexMesh/snappyHexMeshDriver/refinementParameters/refinementParameters.H
index 4e301de0e4d..d05bce28e47 100644
--- a/src/mesh/snappyHexMesh/snappyHexMeshDriver/refinementParameters/refinementParameters.H
+++ b/src/mesh/snappyHexMesh/snappyHexMeshDriver/refinementParameters/refinementParameters.H
@@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2014 OpenFOAM Foundation
- Copyright (C) 2015-2019 OpenCFD Ltd.
+ Copyright (C) 2015-2020 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@@ -113,6 +113,8 @@ class refinementParameters
const label nFilterIter_;
+ const scalar minCellFraction_;
+
const bool dryRun_;
@@ -235,6 +237,13 @@ public:
return nFilterIter_;
}
+ //- When are disconnected regions small. Fraction of overall size
+ // of a zone or background. Default 0.
+ scalar minCellFraction() const
+ {
+ return minCellFraction_;
+ }
+
// Other
diff --git a/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappyRefineDriver.C b/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappyRefineDriver.C
index 5dd6392ac40..374827f252a 100644
--- a/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappyRefineDriver.C
+++ b/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappyRefineDriver.C
@@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2015 OpenFOAM Foundation
- Copyright (C) 2015-2019 OpenCFD Ltd.
+ Copyright (C) 2015-2020 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@@ -47,6 +47,8 @@ License
#include "fvMeshSubset.H"
#include "interpolationTable.H"
#include "snappyVoxelMeshDriver.H"
+#include "regionSplit.H"
+#include "removeCells.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@@ -2869,6 +2871,142 @@ void Foam::snappyRefineDriver::mergePatchFaces
}
+void Foam::snappyRefineDriver::deleteSmallRegions
+(
+ const refinementParameters& refineParams
+)
+{
+ const fvMesh& mesh = meshRefiner_.mesh();
+ const cellZoneMesh& czm = mesh.cellZones();
+
+ //const labelList zoneIDs
+ //(
+ // meshRefiner_.getZones
+ // (
+ // List<surfaceZonesInfo::faceZoneType> fzTypes
+ // ({
+ // surfaceZonesInfo::BAFFLE,
+ // surfaceZonesInfo::BOUNDARY,
+ // });
+ // )
+ //);
+ const labelList zoneIDs(identity(mesh.faceZones().size()));
+
+ // Get faceZone and patch(es) per face (or -1 if face not on faceZone)
+ labelList faceZoneID;
+ labelList ownPatch;
+ labelList neiPatch;
+ labelList nB; // local count of faces per faceZone
+ meshRefiner_.getZoneFaces(zoneIDs, faceZoneID, ownPatch, neiPatch, nB);
+
+
+ // Mark all faces on outside of zones. Note: assumes that faceZones
+ // are consistent with the outside of cellZones ...
+
+ boolList isBlockedFace(mesh.nFaces(), false);
+ meshRefiner_.selectSeparatedCoupledFaces(isBlockedFace);
+
+ forAll(ownPatch, facei)
+ {
+ if (ownPatch[facei] != -1)
+ {
+ isBlockedFace[facei] = true;
+ }
+ }
+
+ // Map from cell to zone. Not that background cells are not -1 but
+ // cellZones.size()
+ labelList cellToZone(mesh.nCells(), czm.size());
+ for (const auto& cz : czm)
+ {
+ UIndirectList<label>(cellToZone, cz) = cz.index();
+ }
+
+ // Walk to split into regions
+ const regionSplit cellRegion(mesh, isBlockedFace);
+
+ // Count number of cells per zone and per region
+ labelList nCellsPerRegion(cellRegion.nRegions(), 0);
+ labelList regionToZone(cellRegion.nRegions(), -2);
+ labelList nCellsPerZone(czm.size()+1, 0);
+ forAll(cellRegion, celli)
+ {
+ const label regioni = cellRegion[celli];
+ const label zonei = cellToZone[celli];
+
+ // Zone for this region
+ regionToZone[regioni] = zonei;
+
+ nCellsPerRegion[regioni]++;
+ nCellsPerZone[zonei]++;
+ }
+ Pstream::listCombineGather(nCellsPerRegion, plusEqOp<label>());
+ Pstream::listCombineGather(regionToZone, maxEqOp<label>());
+ Pstream::listCombineGather(nCellsPerZone, plusEqOp<label>());
+
+
+ // Mark small regions
+ forAll(nCellsPerRegion, regioni)
+ {
+ const label zonei = regionToZone[regioni];
+
+ if
+ (
+ nCellsPerRegion[regioni]
+ < refineParams.minCellFraction()*nCellsPerZone[zonei]
+ )
+ {
+ Info<< "Deleting region " << regioni
+ << " (size " << nCellsPerRegion[regioni]
+ << ") of zone size " << nCellsPerZone[zonei]
+ << endl;
+
+ // Mark region to be deleted. 0 size (= global) should never
+ // occur.
+ nCellsPerRegion[regioni] = 0;
+ }
+ }
+
+ DynamicList<label> cellsToRemove(mesh.nCells()/128);
+ forAll(cellRegion, celli)
+ {
+ if (nCellsPerRegion[cellRegion[celli]] == 0)
+ {
+ cellsToRemove.append(celli);
+ }
+ }
+ const label nTotCellsToRemove = returnReduce
+ (
+ cellsToRemove.size(),
+ sumOp<label>()
+ );
+ if (nTotCellsToRemove > 0)
+ {
+ Info<< "Deleting " << nTotCellsToRemove
+ << " cells in small regions" << endl;
+
+ removeCells cellRemover(mesh);
+
+ cellsToRemove.shrink();
+ const labelList exposedFaces
+ (
+ cellRemover.getExposedFaces(cellsToRemove)
+ );
+ const labelList exposedPatch
+ (
+ UIndirectList<label>(ownPatch, exposedFaces)
+ );
+ (void)meshRefiner_.doRemoveCells
+ (
+ cellsToRemove,
+ exposedFaces,
+ exposedPatch,
+ cellRemover
+ );
+ }
+}
+
+
void Foam::snappyRefineDriver::doRefine
(
const dictionary& refineDict,
@@ -3099,6 +3237,17 @@ void Foam::snappyRefineDriver::doRefine
}
+ if (refineParams.minCellFraction() > 0)
+ {
+ // Some small disconnected bits of mesh might remain since at
+ // this point faceZones have not been converted into e.g. baffles.
+ // We don't know whether e.g. the baffles are reset to be cyclicAMI
+ // thus reconnecting. For now check if there are any particularly
+ // small regions.
+ deleteSmallRegions(refineParams);
+ }
+
+
if (!dryRun_ && Pstream::parRun())
{
Info<< nl
diff --git a/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappyRefineDriver.H b/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappyRefineDriver.H
index d3667968b7a..787190069e9 100644
--- a/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappyRefineDriver.H
+++ b/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappyRefineDriver.H
@@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2014 OpenFOAM Foundation
- Copyright (C) 2015-2018 OpenCFD Ltd.
+ Copyright (C) 2015-2020 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@@ -232,6 +232,10 @@ class snappyRefineDriver
const dictionary& motionDict
);
+ //- Optionally delete some small regions
+ void deleteSmallRegions(const refinementParameters&);
+
+
//- No copy construct
snappyRefineDriver(const snappyRefineDriver&) = delete;
--
GitLab