diff --git a/src/meshTools/regionSplit/regionSplit.C b/src/meshTools/regionSplit/regionSplit.C
index 2583c2ddfe8540b730b718c0ea17448301e9b98d..1a8de0fa152526c37959375d28df6f9fcb638c6e 100644
--- a/src/meshTools/regionSplit/regionSplit.C
+++ b/src/meshTools/regionSplit/regionSplit.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -28,6 +28,8 @@ License
#include "processorPolyPatch.H"
#include "globalIndex.H"
#include "syncTools.H"
+#include "FaceCellWave.H"
+#include "minData.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@@ -40,524 +42,243 @@ defineTypeNameAndDebug(regionSplit, 0);
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
-// Handle (non-processor) coupled faces.
-void Foam::regionSplit::transferCoupledFaceRegion
+void Foam::regionSplit::calcNonCompactRegionSplit
(
- const label faceI,
- const label otherFaceI,
+ const globalIndex& globalFaces,
+ const boolList& blockedFace,
+ const List<labelPair>& explicitConnections,
- labelList& faceRegion,
- DynamicList<label>& newChangedFaces
+ labelList& cellRegion
) const
{
- if (faceRegion[faceI] >= 0)
- {
- if (faceRegion[otherFaceI] == -1)
- {
- faceRegion[otherFaceI] = faceRegion[faceI];
- newChangedFaces.append(otherFaceI);
- }
- else if (faceRegion[otherFaceI] == -2)
- {
- // otherFaceI blocked but faceI is not. Is illegal for coupled
- // faces, not for explicit connections.
- }
- else if (faceRegion[otherFaceI] != faceRegion[faceI])
- {
- FatalErrorIn
- (
- "regionSplit::transferCoupledFaceRegion"
- "(const label, const label, labelList&, labelList&) const"
- ) << "Problem : coupled face " << faceI
- << " on patch " << mesh().boundaryMesh().whichPatch(faceI)
- << " has region " << faceRegion[faceI]
- << " but coupled face " << otherFaceI
- << " has region " << faceRegion[otherFaceI]
- << endl
- << "Is your blocked faces specification"
- << " synchronized across coupled boundaries?"
- << abort(FatalError);
- }
- }
- else if (faceRegion[faceI] == -1)
+ // Field on cells and faces.
+ List<minData> cellData(mesh().nCells());
+ List<minData> faceData(mesh().nFaces());
+
+ // Take over blockedFaces by seeding a negative number
+ // (so is always less than the decomposition)
+ label nUnblocked = 0;
+ forAll(faceData, faceI)
{
- if (faceRegion[otherFaceI] >= 0)
+ if (blockedFace.size() && blockedFace[faceI])
{
- faceRegion[faceI] = faceRegion[otherFaceI];
- newChangedFaces.append(faceI);
+ faceData[faceI] = minData(-2);
}
- else if (faceRegion[otherFaceI] == -2)
+ else
{
- // otherFaceI blocked but faceI is not. Is illegal for coupled
- // faces, not for explicit connections.
+ nUnblocked++;
}
}
-}
-
-
-void Foam::regionSplit::fillSeedMask
-(
- const List<labelPair>& explicitConnections,
- labelList& cellRegion,
- labelList& faceRegion,
- const label seedCellID,
- const label markValue
-) const
-{
- // Do seed cell
- cellRegion[seedCellID] = markValue;
-
- // Collect faces on seed cell
- const cell& cFaces = mesh().cells()[seedCellID];
+ // Seed unblocked faces
+ labelList seedFaces(nUnblocked);
+ List<minData> seedData(nUnblocked);
+ nUnblocked = 0;
- label nFaces = 0;
- labelList changedFaces(cFaces.size());
-
- forAll(cFaces, i)
+ forAll(faceData, faceI)
{
- label faceI = cFaces[i];
-
- if (faceRegion[faceI] == -1)
+ if (blockedFace.empty() || !blockedFace[faceI])
{
- faceRegion[faceI] = markValue;
- changedFaces[nFaces++] = faceI;
+ seedFaces[nUnblocked] = faceI;
+ // Seed face with globally unique number
+ seedData[nUnblocked] = minData(globalFaces.toGlobal(faceI));
+ nUnblocked++;
}
}
- changedFaces.setSize(nFaces);
-
- // Loop over changed faces. FaceCellWave in small.
- while (changedFaces.size())
- {
- //if (debug)
- //{
- // Pout<< "regionSplit::fillSeedMask : changedFaces:"
- // << changedFaces.size() << endl;
- //}
+ // Propagate information inwards
+ FaceCellWave<minData> deltaCalc
+ (
+ mesh(),
+ explicitConnections,
+ false, // disable walking through cyclicAMI for backwards compatibility
+ seedFaces,
+ seedData,
+ faceData,
+ cellData,
+ mesh().globalData().nTotalCells()+1
+ );
- DynamicList<label> changedCells(changedFaces.size());
- forAll(changedFaces, i)
+ // And extract
+ cellRegion.setSize(mesh().nCells());
+ forAll(cellRegion, cellI)
+ {
+ if (cellData[cellI].valid(deltaCalc.data()))
{
- label faceI = changedFaces[i];
-
- label own = mesh().faceOwner()[faceI];
-
- if (cellRegion[own] == -1)
- {
- cellRegion[own] = markValue;
- changedCells.append(own);
- }
-
- if (mesh().isInternalFace(faceI))
- {
- label nei = mesh().faceNeighbour()[faceI];
-
- if (cellRegion[nei] == -1)
- {
- cellRegion[nei] = markValue;
- changedCells.append(nei);
- }
- }
+ cellRegion[cellI] = cellData[cellI].data();
}
-
-
- //if (debug)
- //{
- // Pout<< "regionSplit::fillSeedMask : changedCells:"
- // << changedCells.size() << endl;
- //}
-
- // Loop over changedCells and collect faces
- DynamicList<label> newChangedFaces(changedCells.size());
-
- forAll(changedCells, i)
+ else
{
- label cellI = changedCells[i];
-
+ // Unvisited cell -> only possible if surrounded by blocked faces.
+ // If so make up region from any of the faces
const cell& cFaces = mesh().cells()[cellI];
+ label faceI = cFaces[0];
- forAll(cFaces, cFaceI)
+ if (blockedFace.size() && !blockedFace[faceI])
{
- label faceI = cFaces[cFaceI];
-
- if (faceRegion[faceI] == -1)
- {
- faceRegion[faceI] = markValue;
- newChangedFaces.append(faceI);
- }
- }
- }
-
-
- //if (debug)
- //{
- // Pout<< "regionSplit::fillSeedMask : changedFaces before sync:"
- // << changedFaces.size() << endl;
- //}
-
-
- // Check for changes to any locally coupled face.
- // Global connections are done later.
-
- const polyBoundaryMesh& patches = mesh().boundaryMesh();
-
- forAll(patches, patchI)
- {
- const polyPatch& pp = patches[patchI];
-
- if
- (
- isA<cyclicPolyPatch>(pp)
- && refCast<const cyclicPolyPatch>(pp).owner()
- )
- {
- // Transfer from neighbourPatch to here or vice versa.
-
- const cyclicPolyPatch& cycPatch =
- refCast<const cyclicPolyPatch>(pp);
-
- label faceI = cycPatch.start();
-
- forAll(cycPatch, i)
- {
- label otherFaceI = cycPatch.transformGlobalFace(faceI);
-
- transferCoupledFaceRegion
- (
- faceI,
- otherFaceI,
- faceRegion,
- newChangedFaces
- );
-
- faceI++;
- }
+ FatalErrorIn("regionSplit::calcNonCompactRegionSplit(..)")
+ << "Problem: unblocked face " << faceI
+ << " at " << mesh().faceCentres()[faceI]
+ << " on unassigned cell " << cellI
+ << mesh().cellCentres()[faceI]
+ << exit(FatalError);
}
+ cellRegion[cellI] = globalFaces.toGlobal(faceI);
}
- forAll(explicitConnections, i)
- {
- transferCoupledFaceRegion
- (
- explicitConnections[i][0],
- explicitConnections[i][1],
- faceRegion,
- newChangedFaces
- );
- }
-
- //if (debug)
- //{
- // Pout<< "regionSplit::fillSeedMask : changedFaces after sync:"
- // << newChangedFaces.size() << endl;
- //}
-
- changedFaces.transfer(newChangedFaces);
}
}
-Foam::label Foam::regionSplit::calcLocalRegionSplit
+Foam::autoPtr<Foam::globalIndex> Foam::regionSplit::calcRegionSplit
(
+ const bool doGlobalRegions,
const boolList& blockedFace,
const List<labelPair>& explicitConnections,
labelList& cellRegion
) const
{
- if (debug)
- {
- if (blockedFace.size())
- {
- // Check that blockedFace is synced.
- boolList syncBlockedFace(blockedFace);
- syncTools::swapFaceList(mesh(), syncBlockedFace);
-
- forAll(syncBlockedFace, faceI)
- {
- if (syncBlockedFace[faceI] != blockedFace[faceI])
- {
- FatalErrorIn
- (
- "regionSplit::calcLocalRegionSplit(..)"
- ) << "Face " << faceI << " not synchronised. My value:"
- << blockedFace[faceI] << " coupled value:"
- << syncBlockedFace[faceI]
- << abort(FatalError);
- }
- }
- }
- }
-
- // Region per face.
- // -1 unassigned
- // -2 blocked
- labelList faceRegion(mesh().nFaces(), -1);
-
- if (blockedFace.size())
- {
- forAll(blockedFace, faceI)
- {
- if (blockedFace[faceI])
- {
- faceRegion[faceI] = -2;
- }
- }
- }
-
-
- // Assign local regions
- // ~~~~~~~~~~~~~~~~~~~~
-
- // Start with region 0
- label nLocalRegions = 0;
+ // See header in regionSplit.H
- label unsetCellI = 0;
- do
+ if (!doGlobalRegions)
{
- // Find first unset cell
+ // Block all parallel faces to avoid comms across
+ boolList coupledOrBlockedFace(blockedFace);
+ const polyBoundaryMesh& pbm = mesh().boundaryMesh();
- for (; unsetCellI < mesh().nCells(); unsetCellI++)
+ if (coupledOrBlockedFace.size())
{
- if (cellRegion[unsetCellI] == -1)
+ forAll(pbm, patchI)
{
- break;
+ const polyPatch& pp = pbm[patchI];
+ if (isA<processorPolyPatch>(pp))
+ {
+ label faceI = pp.start();
+ forAll(pp, i)
+ {
+ coupledOrBlockedFace[faceI++] = true;
+ }
+ }
}
}
- if (unsetCellI >= mesh().nCells())
+ // Create dummy (local only) globalIndex
+ labelList offsets(Pstream::nProcs()+1, 0);
+ for (label i = Pstream::myProcNo()+1; i < offsets.size(); i++)
{
- break;
+ offsets[i] = mesh().nFaces();
}
+ const globalIndex globalRegions(offsets.xfer());
- fillSeedMask
+ // Minimise regions across connected cells
+ // Note: still uses global decisions so all processors are running
+ // in lock-step, i.e. slowest determines overall time.
+ // To avoid this we could switch off Pstream::parRun.
+ calcNonCompactRegionSplit
(
+ globalRegions,
+ coupledOrBlockedFace,
explicitConnections,
- cellRegion,
- faceRegion,
- unsetCellI,
- nLocalRegions
+ cellRegion
);
- // Current unsetCell has now been handled. Go to next region.
- nLocalRegions++;
- unsetCellI++;
- }
- while (true);
-
-
- if (debug)
- {
+ // Compact
+ Map<label> globalToCompact(mesh().nCells()/8);
forAll(cellRegion, cellI)
{
- if (cellRegion[cellI] < 0)
+ label region = cellRegion[cellI];
+
+ label globalRegion;
+
+ Map<label>::const_iterator fnd = globalToCompact.find(region);
+ if (fnd == globalToCompact.end())
{
- FatalErrorIn("regionSplit::calcLocalRegionSplit(..)")
- << "cell:" << cellI << " region:" << cellRegion[cellI]
- << abort(FatalError);
+ globalRegion = globalRegions.toGlobal(globalToCompact.size());
+ globalToCompact.insert(region, globalRegion);
}
+ else
+ {
+ globalRegion = fnd();
+ }
+ cellRegion[cellI] = globalRegion;
}
- forAll(faceRegion, faceI)
+
+ // Return globalIndex with size = localSize and all regions local
+ labelList compactOffsets(Pstream::nProcs()+1, 0);
+ for (label i = Pstream::myProcNo()+1; i < compactOffsets.size(); i++)
{
- if (faceRegion[faceI] == -1)
- {
- FatalErrorIn("regionSplit::calcLocalRegionSplit(..)")
- << "face:" << faceI << " region:" << faceRegion[faceI]
- << abort(FatalError);
- }
+ compactOffsets[i] = globalToCompact.size();
}
+
+ return autoPtr<globalIndex>(new globalIndex(compactOffsets.xfer()));
}
- return nLocalRegions;
-}
-Foam::autoPtr<Foam::globalIndex> Foam::regionSplit::calcRegionSplit
-(
- const bool doGlobalRegions,
- const boolList& blockedFace,
- const List<labelPair>& explicitConnections,
-
- labelList& cellRegion
-) const
-{
- // See header in regionSplit.H
+ // Initial global region numbers
+ const globalIndex globalRegions(mesh().nFaces());
- // 1. Do local analysis
- label nLocalRegions = calcLocalRegionSplit
+ // Minimise regions across connected cells (including parallel)
+ calcNonCompactRegionSplit
(
+ globalRegions,
blockedFace,
explicitConnections,
cellRegion
);
- if (!doGlobalRegions)
- {
- return autoPtr<globalIndex>(new globalIndex(nLocalRegions));
- }
-
- // 2. Assign global regions
- // ~~~~~~~~~~~~~~~~~~~~~~~~
- // Offset local regions to create unique global regions.
+ // Now our cellRegion will have
+ // - non-local regions (i.e. originating from other processors)
+ // - non-compact locally originating regions
+ // so we'll need to compact
- globalIndex globalRegions(nLocalRegions);
-
-
- // Convert regions to global ones
- forAll(cellRegion, cellI)
+ // 4a: count per originating processor the number of regions
+ labelList nOriginating(Pstream::nProcs(), 0);
{
- cellRegion[cellI] = globalRegions.toGlobal(cellRegion[cellI]);
- }
-
-
- // 3. Merge global regions
- // ~~~~~~~~~~~~~~~~~~~~~~~
- // Regions across non-blocked proc patches get merged.
- // This will set merged global regions to be the min of both.
- // (this will create gaps in the global region list so they will get
- // merged later on)
-
- while (true)
- {
- if (debug)
- {
- Pout<< nl << "-- Starting Iteration --" << endl;
- }
-
+ labelHashSet haveRegion(mesh().nCells()/8);
- const polyBoundaryMesh& patches = mesh().boundaryMesh();
-
- labelList nbrRegion(mesh().nFaces()-mesh().nInternalFaces(), -1);
- forAll(patches, patchI)
+ forAll(cellRegion, cellI)
{
- const polyPatch& pp = patches[patchI];
+ label region = cellRegion[cellI];
- if (pp.coupled())
+ // Count originating processor. Use isLocal as efficiency since
+ // most cells are locally originating.
+ if (globalRegions.isLocal(region))
{
- const labelList& patchCells = pp.faceCells();
- SubList<label> patchNbrRegion
- (
- nbrRegion,
- pp.size(),
- pp.start()-mesh().nInternalFaces()
- );
-
- forAll(patchCells, i)
+ if (haveRegion.insert(region))
{
- label faceI = pp.start()+i;
- if (!blockedFace.size() || !blockedFace[faceI])
- {
- patchNbrRegion[i] = cellRegion[patchCells[i]];
- }
+ nOriginating[Pstream::myProcNo()]++;
}
}
- }
- syncTools::swapBoundaryFaceList(mesh(), nbrRegion);
-
- Map<label> globalToMerged(mesh().nFaces()-mesh().nInternalFaces());
-
- forAll(patches, patchI)
- {
- const polyPatch& pp = patches[patchI];
-
- if (pp.coupled())
+ else
{
- const labelList& patchCells = pp.faceCells();
- SubList<label> patchNbrRegion
- (
- nbrRegion,
- pp.size(),
- pp.start()-mesh().nInternalFaces()
- );
-
- forAll(patchCells, i)
+ label procI = globalRegions.whichProcID(region);
+ if (haveRegion.insert(region))
{
- label faceI = pp.start()+i;
-
- if (!blockedFace.size() || !blockedFace[faceI])
- {
- if (patchNbrRegion[i] < cellRegion[patchCells[i]])
- {
- //Pout<< "on patch:" << pp.name()
- // << " cell:" << patchCells[i]
- // << " at:"
- // << mesh().cellCentres()[patchCells[i]]
- // << " was:" << cellRegion[patchCells[i]]
- // << " nbr:" << patchNbrRegion[i]
- // << endl;
-
- globalToMerged.insert
- (
- cellRegion[patchCells[i]],
- patchNbrRegion[i]
- );
- }
- }
+ nOriginating[procI]++;
}
}
}
-
-
- label nMerged = returnReduce(globalToMerged.size(), sumOp<label>());
-
- if (debug)
- {
- Pout<< "nMerged:" << nMerged << endl;
- }
-
- if (nMerged == 0)
- {
- break;
- }
-
- // Renumber the regions according to the globalToMerged
- forAll(cellRegion, cellI)
- {
- label regionI = cellRegion[cellI];
- Map<label>::const_iterator iter = globalToMerged.find(regionI);
- if (iter != globalToMerged.end())
- {
- cellRegion[cellI] = iter();
- }
- }
- }
-
-
- // Now our cellRegion will have non-local elements in it. So compact
- // it.
-
- // 4a: count. Use a labelHashSet to count regions only once.
- label nCompact = 0;
- {
- labelHashSet localRegion(mesh().nFaces()-mesh().nInternalFaces());
- forAll(cellRegion, cellI)
- {
- if
- (
- globalRegions.isLocal(cellRegion[cellI])
- && localRegion.insert(cellRegion[cellI])
- )
- {
- nCompact++;
- }
- }
}
if (debug)
{
- Pout<< "Compacted from " << nLocalRegions
- << " down to " << nCompact << " local regions." << endl;
+ Pout<< "Counted " << nOriginating[Pstream::myProcNo()]
+ << " local regions." << endl;
}
// Global numbering for compacted local regions
- autoPtr<globalIndex> globalCompactPtr(new globalIndex(nCompact));
+ autoPtr<globalIndex> globalCompactPtr
+ (
+ new globalIndex(nOriginating[Pstream::myProcNo()])
+ );
const globalIndex& globalCompact = globalCompactPtr();
@@ -568,12 +289,15 @@ Foam::autoPtr<Foam::globalIndex> Foam::regionSplit::calcRegionSplit
// labelList.
// Local compaction map
- Map<label> globalToCompact(2*nCompact);
+ Map<label> globalToCompact(2*nOriginating[Pstream::myProcNo()]);
// Remote regions we want the compact number for
List<labelHashSet> nonLocal(Pstream::nProcs());
forAll(nonLocal, procI)
{
- nonLocal[procI].resize((nLocalRegions-nCompact)/Pstream::nProcs());
+ if (procI != Pstream::myProcNo())
+ {
+ nonLocal[procI].resize(2*nOriginating[procI]);
+ }
}
forAll(cellRegion, cellI)
@@ -581,15 +305,12 @@ Foam::autoPtr<Foam::globalIndex> Foam::regionSplit::calcRegionSplit
label region = cellRegion[cellI];
if (globalRegions.isLocal(region))
{
- Map<label>::const_iterator iter = globalToCompact.find(region);
- if (iter == globalToCompact.end())
- {
- label compactRegion = globalCompact.toGlobal
- (
- globalToCompact.size()
- );
- globalToCompact.insert(region, compactRegion);
- }
+ // Insert new compact region (if not yet present)
+ globalToCompact.insert
+ (
+ region,
+ globalCompact.toGlobal(globalToCompact.size())
+ );
}
else
{
@@ -603,22 +324,17 @@ Foam::autoPtr<Foam::globalIndex> Foam::regionSplit::calcRegionSplit
// Convert the nonLocal (labelHashSets) to labelLists.
labelListList sendNonLocal(Pstream::nProcs());
- labelList nNonLocal(Pstream::nProcs(), 0);
forAll(sendNonLocal, procI)
{
- sendNonLocal[procI].setSize(nonLocal[procI].size());
- forAllConstIter(labelHashSet, nonLocal[procI], iter)
- {
- sendNonLocal[procI][nNonLocal[procI]++] = iter.key();
- }
+ sendNonLocal[procI] = nonLocal[procI].toc();
}
if (debug)
{
- forAll(nNonLocal, procI)
+ forAll(sendNonLocal, procI)
{
Pout<< " from processor " << procI
- << " want " << nNonLocal[procI]
+ << " want " << sendNonLocal[procI].size()
<< " region numbers."
<< endl;
}
@@ -627,7 +343,7 @@ Foam::autoPtr<Foam::globalIndex> Foam::regionSplit::calcRegionSplit
// Get the wanted region labels into recvNonLocal
- labelListList recvNonLocal;
+ labelListList recvNonLocal(Pstream::nProcs());
labelListList sizes;
Pstream::exchange<labelList, label>
(
@@ -655,6 +371,7 @@ Foam::autoPtr<Foam::globalIndex> Foam::regionSplit::calcRegionSplit
// Send back (into recvNonLocal)
recvNonLocal.clear();
+ recvNonLocal.setSize(sendWantedLocal.size());
sizes.clear();
Pstream::exchange<labelList, label>
(
diff --git a/src/meshTools/regionSplit/regionSplit.H b/src/meshTools/regionSplit/regionSplit.H
index 83302529df24349890b4c4425d21fe4a0c5dfdae..543259a10bfeeb6345e4b84ab9138e84676605d0 100644
--- a/src/meshTools/regionSplit/regionSplit.H
+++ b/src/meshTools/regionSplit/regionSplit.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -89,6 +89,9 @@ Description
Can optionally keep all regions local to the processor.
+ Note: does not walk across cyclicAMI/cyclicACMI - since not 'coupled()'
+ at the patch level.
+
SourceFiles
regionSplit.C
@@ -126,30 +129,10 @@ class regionSplit
// Private Member Functions
- //- Transfer faceRegion data from one face to the other (or vice versa)
- void transferCoupledFaceRegion
- (
- const label faceI,
- const label otherFaceI,
-
- labelList& faceRegion,
- DynamicList<label>& newChangedFaces
- ) const;
-
- //- Given a seed cell label, fill cellRegion/faceRegion with markValue
- // for contiguous region around it
- void fillSeedMask
- (
- const List<labelPair>& explicitConnections,
- labelList& cellRegion,
- labelList& faceRegion,
- const label seedCellID,
- const label markValue
- ) const;
-
- //- Calculate local region split. Return number of regions.
- label calcLocalRegionSplit
+ //- Calculate region split in non-compact (global) numbering.
+ void calcNonCompactRegionSplit
(
+ const globalIndex& globalFaces,
const boolList& blockedFace,
const List<labelPair>& explicitConnections,
labelList& cellRegion