From 145cad56f1da9b1131f7c6907b4a454697a0329e Mon Sep 17 00:00:00 2001
From: Mark Olesen <Mark.Olesen@esi-group.com>
Date: Fri, 1 Oct 2021 18:26:49 +0200
Subject: [PATCH] BUG: parallel construct finiteArea fails with arbitrary
connections (#2152)
- the case of 'fan-like' processor was previously assumed to be
rare (see merge-request !490 and issue #2084).
However, Vaggelis Papoutsis noticed that even fairly normal geometries
can trigger problems.
- replaced the old patch/patch matching style with a more general
edge-based synchronisation and matching that appears to handle
the corner cases inherently. The internal communication overhead
is essentially unchanged, and the logic is simpler.
ENH: additional framework for managing patch connectivity
---
src/finiteArea/Make/files | 1 +
src/finiteArea/faMesh/faMesh.C | 19 +-
src/finiteArea/faMesh/faMesh.H | 150 ++-
src/finiteArea/faMesh/faMeshI.H | 11 +
src/finiteArea/faMesh/faMeshPatches.C | 953 +++++-------------
src/finiteArea/faMesh/faMeshTopology.C | 803 +++++++++++++++
.../faMesh/faPatches/faPatch/faPatch.C | 76 +-
.../faMesh/faPatches/faPatch/faPatch.H | 23 +-
.../faMesh/faPatches/faPatch/faPatchData.C | 19 +
.../faMesh/faPatches/faPatch/faPatchData.H | 3 +
10 files changed, 1337 insertions(+), 721 deletions(-)
create mode 100644 src/finiteArea/faMesh/faMeshTopology.C
diff --git a/src/finiteArea/Make/files b/src/finiteArea/Make/files
index 14d538ace1f..85d56a5e87c 100644
--- a/src/finiteArea/Make/files
+++ b/src/finiteArea/Make/files
@@ -2,6 +2,7 @@ faMesh/faGlobalMeshData/faGlobalMeshData.C
faMesh/faMesh.C
faMesh/faMeshDemandDrivenData.C
faMesh/faMeshPatches.C
+faMesh/faMeshTopology.C
faMesh/faMeshUpdate.C
faMesh/faBoundaryMesh/faBoundaryMesh.C
diff --git a/src/finiteArea/faMesh/faMesh.C b/src/finiteArea/faMesh/faMesh.C
index ef64e0541ec..021b83a3270 100644
--- a/src/finiteArea/faMesh/faMesh.C
+++ b/src/finiteArea/faMesh/faMesh.C
@@ -111,15 +111,15 @@ static labelList selectPatchFaces
void Foam::faMesh::initPatch() const
{
- if (patchPtr_)
- {
- delete patchPtr_;
- }
- patchPtr_ = new uindirectPrimitivePatch
+ patchPtr_.reset
(
- UIndirectList<face>(mesh().faces(), faceLabels_),
- mesh().points()
+ new uindirectPrimitivePatch
+ (
+ UIndirectList<face>(mesh().faces(), faceLabels_),
+ mesh().points()
+ )
);
+ bndConnectPtr_.reset(nullptr);
}
@@ -172,8 +172,9 @@ void Foam::faMesh::clearGeomNotAreas() const
{
DebugInFunction << "Clearing geometry" << endl;
+ patchPtr_.reset(nullptr);
+ bndConnectPtr_.reset(nullptr);
deleteDemandDrivenData(SPtr_);
- deleteDemandDrivenData(patchPtr_);
deleteDemandDrivenData(patchStartsPtr_);
deleteDemandDrivenData(LePtr_);
deleteDemandDrivenData(magLePtr_);
@@ -256,6 +257,7 @@ Foam::faMesh::faMesh(const polyMesh& pMesh)
),
comm_(Pstream::worldComm),
patchPtr_(nullptr),
+ bndConnectPtr_(nullptr),
lduPtr_(nullptr),
curTimeIndex_(time().timeIndex()),
SPtr_(nullptr),
@@ -349,6 +351,7 @@ Foam::faMesh::faMesh
),
comm_(Pstream::worldComm),
patchPtr_(nullptr),
+ bndConnectPtr_(nullptr),
lduPtr_(nullptr),
curTimeIndex_(time().timeIndex()),
SPtr_(nullptr),
diff --git a/src/finiteArea/faMesh/faMesh.H b/src/finiteArea/faMesh/faMesh.H
index 8a9b1650118..313db0517aa 100644
--- a/src/finiteArea/faMesh/faMesh.H
+++ b/src/finiteArea/faMesh/faMesh.H
@@ -33,6 +33,9 @@ Description
SourceFiles
faMesh.C
faMeshDemandDrivenData.C
+ faMeshPatches.C
+ faMeshTopology.C
+ faMeshUpdate.C
Author
Zeljko Tukovic, FMENA
@@ -71,7 +74,6 @@ namespace Foam
class faMeshLduAddressing;
class faMeshMapper;
class faPatchData;
-template<class T> class LabelledItem;
/*---------------------------------------------------------------------------*\
Class faMesh Declaration
@@ -86,6 +88,107 @@ class faMesh
public faSolution,
public data
{
+ // Private (internal) classes/structures
+
+ //- A (proc, patchi, patchEdgei) tuple used internally for managing
+ //- patch/patch bookkeeping during construction.
+ // Finite-area patches are stored with negated indices, which makes
+ // them readily identifiable and always sort before normal patches.
+ // Note
+ struct patchTuple
+ :
+ public FixedList<label, 4>
+ {
+
+ // Constructors
+
+ // Inherit constructors
+ using FixedList<label, 4>::FixedList;
+
+ //- Default construct as 'invalid'
+ patchTuple()
+ {
+ clear();
+ }
+
+
+ // Static Member Functions
+
+ // Globally consistent ordering
+ //
+ // 1. sort left/right as lower/higher processor connection
+ // 2. sort by proc/patch/patch index
+ static void sort(UList<Pair<patchTuple>>& list)
+ {
+ for (auto& tuples : list)
+ {
+ tuples.sort();
+ }
+ Foam::stableSort(list);
+ }
+
+
+ // Member Functions
+
+ //- Reset to 'invalid'
+ void clear()
+ {
+ procNo(-1);
+ patchi(labelMax);
+ patchEdgei(-1);
+ meshFacei(-1);
+ }
+
+ //- Valid if proc and edge are non-negative
+ bool valid() const noexcept
+ {
+ return (procNo() >= 0 && patchEdgei() >= 0);
+ }
+
+ // Processor is the first sort index
+ label procNo() const { return (*this)[0]; }
+ void procNo(label val) { (*this)[0] = val; }
+
+ // PatchId (-ve for finiteArea patches) is the second sort index
+ label patchi() const { return (*this)[1]; }
+ void patchi(label val) { (*this)[1] = val; }
+
+ // The patch edge index (on the finiteArea patch)
+ // is the third sort index
+ label patchEdgei() const { return (*this)[2]; }
+ void patchEdgei(label val) { (*this)[2] = val; }
+
+ // The processor-local mesh face is the fourth sort index
+ label meshFacei() const { return (*this)[3]; }
+ void meshFacei(label val) { (*this)[3] = val; }
+
+ //- Return the real patchId
+ label realPatchi() const
+ {
+ const label id = patchi();
+ return (id < 0 ? -(id + 1) : id);
+ }
+
+ //- Set patchId as finiteArea
+ void faPatchi(label val)
+ {
+ patchi(-(val + 1));
+ }
+
+ //- Considered to be finiteArea if patchi < 0
+ bool is_finiteArea() const noexcept
+ {
+ return (patchi() < 0);
+ }
+
+ //- Considered to be processor local
+ bool is_localProc() const noexcept
+ {
+ return (procNo() == Pstream::myProcNo());
+ }
+ };
+
+
// Private Data
//- Face labels
@@ -131,7 +234,10 @@ class faMesh
// Demand-driven data
//- Primitive patch
- mutable uindirectPrimitivePatch* patchPtr_;
+ mutable std::unique_ptr<uindirectPrimitivePatch> patchPtr_;
+
+ //- List of proc/mesh-face for boundary edge neighbours
+ mutable std::unique_ptr<List<labelPair>> bndConnectPtr_;
//- Ldu addressing data
mutable faMeshLduAddressing* lduPtr_;
@@ -211,6 +317,20 @@ class faMesh
//- Set primitive mesh data
void setPrimitiveMeshData();
+ //- Get list of (proc/patchi/patchEdgei/meshFacei) tuple pairs in an
+ //- globally consistent ordering
+ List<Pair<patchTuple>> getBoundaryEdgeConnections() const;
+
+ //- Determine the boundary edge neighbour connections
+ void calcBoundaryConnections() const;
+
+ //- Define boundary edge neighbours (proc/face) based on
+ //- gathered topology information
+ void setBoundaryConnections
+ (
+ const List<Pair<patchTuple>>& bndEdgeConnections
+ ) const;
+
// Private member functions to calculate demand driven data
@@ -268,9 +388,6 @@ class faMesh
// Helpers
- //- Get the polyPatch pairs for the boundary edges (natural order)
- List<LabelledItem<edge>> getBoundaryEdgePatchPairs() const;
-
//- Create a single patch
PtrList<faPatch> createOnePatch
(
@@ -286,14 +403,6 @@ class faMesh
const dictionary* defaultPatchDefinition = nullptr
) const;
- //- Reorder processor edges using order of the
- //- neighbour processorPolyPatch
- void reorderProcEdges
- (
- faPatchData& patchDef,
- const List<LabelledItem<edge>>& bndEdgePatchPairs
- ) const;
-
public:
@@ -451,7 +560,7 @@ public:
//- Return constant reference to boundary mesh
inline const faBoundaryMesh& boundary() const noexcept;
- //- Return faMesh face labels
+ //- Return the underlying polyMesh face labels
inline const labelList& faceLabels() const noexcept;
@@ -486,6 +595,18 @@ public:
return (edgeIndex < nInternalEdges_);
}
+ //- List of proc/face for the boundary edge neighbours
+ //- using primitive patch edge numbering.
+ inline const List<labelPair>& boundaryConnections() const;
+
+ //- Boundary edge neighbour processors
+ //- (does not include own proc)
+ labelList boundaryProcs() const;
+
+ //- List of proc/size for the boundary edge neighbour processors
+ //- (does not include own proc)
+ List<labelPair> boundaryProcSizes() const;
+
// Mesh motion and morphing
@@ -590,6 +711,7 @@ public:
} // End namespace Foam
+
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "faMeshI.H"
diff --git a/src/finiteArea/faMesh/faMeshI.H b/src/finiteArea/faMesh/faMeshI.H
index a3d79328606..dc1b238c3d4 100644
--- a/src/finiteArea/faMesh/faMeshI.H
+++ b/src/finiteArea/faMesh/faMeshI.H
@@ -139,4 +139,15 @@ inline Foam::uindirectPrimitivePatch& Foam::faMesh::patch()
}
+inline const Foam::List<Foam::labelPair>&
+Foam::faMesh::boundaryConnections() const
+{
+ if (!bndConnectPtr_)
+ {
+ calcBoundaryConnections();
+ }
+ return *bndConnectPtr_;
+}
+
+
// ************************************************************************* //
diff --git a/src/finiteArea/faMesh/faMeshPatches.C b/src/finiteArea/faMesh/faMeshPatches.C
index 436f3f5678e..42d8d99d6aa 100644
--- a/src/finiteArea/faMesh/faMeshPatches.C
+++ b/src/finiteArea/faMesh/faMeshPatches.C
@@ -27,311 +27,10 @@ License
\*---------------------------------------------------------------------------*/
#include "faMesh.H"
-#include "IndirectList.H"
#include "faPatchData.H"
#include "processorPolyPatch.H"
#include "processorFaPatch.H"
-#include "globalMeshData.H"
-#include "indirectPrimitivePatch.H"
-#include "edgeHashes.H"
-#include "LabelledItem.H"
-
-// * * * * * * * * * * * * * * * Local Functions * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
-// Manage patch pairs with a 'labelled' edge.
-// The edge first/second correspond to the owner/neighbour patches.
-// The index is a face index on the neighbour patch (FUTURE).
-
-// Local typedefs
-
-typedef LabelledItem<edge> patchPairInfo;
-typedef List<patchPairInfo> patchPairInfoList;
-typedef UIndirectList<patchPairInfo> patchPairInfoUIndList;
-
-
-// Handling of dangling coupled edges.
-// Tag values to "push" with special -(patchId+2)
-struct combineDanglingEdge
-{
- const label upperLimit;
-
- // Set dangling patchId from real patchId
- static void setDangling(patchPairInfo& pairing, const label patchId)
- {
- pairing.first() = pairing.second() = -(patchId + 2);
- pairing.setIndex(-1); // Invalidate
- }
-
- // Convert dangling patchId to real patchId
- static void correct(patchPairInfo& pairing)
- {
- if (pairing.first() < -1)
- {
- pairing.first() = -(pairing.first() + 2);
- }
- if (pairing.second() < -1)
- {
- pairing.second() = -(pairing.second() + 2);
- }
- }
-
- //- Construct with upper limit (the number of non-processor patches)
- explicit combineDanglingEdge(const label nNonProcessor)
- :
- upperLimit(nNonProcessor)
- {}
-
-
- // Combine operation: overwrite unused or processor patches with
- // 'dangling' patch information only
- void operator()(patchPairInfo& x, const patchPairInfo& y) const
- {
- if (y.first() < -1 && edge::compare(x, y) == 0)
- {
- if (x.first() == -1 || x.first() >= upperLimit)
- {
- x.first() = y.first();
- }
- if (x.second() == -1 || x.second() >= upperLimit)
- {
- x.second() = y.first();
- x.index() = y.index();
- }
- }
- }
-};
-
-
-// Populate patch pairings according to the boundary edges
-void findEdgePatchPairing
-(
- const polyBoundaryMesh& pbm,
- const EdgeMap<label>& edgeToIndex,
- patchPairInfoList& patchPairs,
- label nMissing = -1
-)
-{
- // Count how many slots (both sides) to be filled
- if (nMissing < 0)
- {
- nMissing = 0;
- for (const patchPairInfo& pairing : patchPairs)
- {
- if (pairing.first() == -1) ++nMissing;
- if (pairing.second() == -1) ++nMissing;
- }
- }
-
- forAll(pbm, patchID)
- {
- if (!nMissing) break; // Everything filled
-
- const polyPatch& pp = pbm[patchID];
-
- const bool isProcPatch = isA<processorPolyPatch>(pp);
-
- // Examine neighbour boundary edges
- for (label edgei = pp.nInternalEdges(); edgei < pp.nEdges(); ++edgei)
- {
- // Lookup global edge of this neighbour,
- // find matching owner boundary edge
-
- const label edgeIndex =
- edgeToIndex.lookup(pp.meshEdge(edgei), -1);
-
- if (edgeIndex != -1)
- {
- // Add patchId.
- // - hash-like so will only insert once
- // - also add in the attached patch face (there is only one),
- // saved in mesh face numbering for processor patches
-
- patchPairInfo& pairing = patchPairs[edgeIndex];
-
- if (pairing.insert(patchID))
- {
- if (isProcPatch)
- {
- // Save the mesh face
-
- const label patchFacei = pp.edgeFaces()[edgei][0];
-
- pairing.index() = (pp.start() + patchFacei);
- }
-
- --nMissing;
-
- if (!nMissing) break; // Early exit
- }
- }
- }
- }
-}
-
-} // End namespace Foam
-
-
-// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
-
-void Foam::faMesh::reorderProcEdges
-(
- faPatchData& patchDef,
- const List<LabelledItem<edge>>& bndEdgePatchPairs
-) const
-{
- if (!patchDef.coupled() || patchDef.edgeLabels_.empty())
- {
- return;
- }
-
- const polyBoundaryMesh& pbm = mesh().boundaryMesh();
-
- const label procPatchID = patchDef.neighPolyPatchId_;
-
- const auto* procPatch = isA<processorPolyPatch>(pbm[procPatchID]);
-
- if (!procPatch)
- {
- FatalErrorInFunction
- << "Internal addressing error. Patch " << procPatchID
- << " is not a processor patch" << nl
- << abort(FatalError);
- }
-
- // Reorder processor edges using order of neighbour processorPolyPatch
-
- const label nProcEdges = patchDef.edgeLabels_.size();
-
- labelList procFaces(nProcEdges, -1);
-
- forAll(procFaces, edgei)
- {
- const label bndEdgei =
- (patchDef.edgeLabels_[edgei] - patch().nInternalEdges());
-
- procFaces[edgei] = bndEdgePatchPairs[bndEdgei].index();
- }
-
- // Ascending proc-face numbering
- const labelList sortIndices(Foam::sortedOrder(procFaces));
-
- if (procFaces.size() && procFaces[sortIndices[0]] < 0)
- {
- FatalErrorInFunction
- << "Internal addressing error. Patch " << procPatchID
- << " with negative face" << nl
- << abort(FatalError);
- }
-
-
- const labelList& oldEdgeLabels = patchDef.edgeLabels_;
- labelList newEdgeLabels(oldEdgeLabels.size());
-
- // Most of the time, an individual proc-face will only be singly
- // attached to the finite-area patch. In rarer case, there could
- // multiple connections. For these cases, need to walk the face
- // edges - the direction depends on owner vs neighbour side.
-
- EdgeMap<label> multihit;
-
- for (label edgei = 0; edgei < nProcEdges; /*nil*/)
- {
- const label meshFacei = procFaces[sortIndices[edgei]];
-
- // Find all identical faces
- label endEdgei = edgei + 1; // one beyond
- while
- (
- (endEdgei < nProcEdges)
- && (meshFacei == procFaces[sortIndices[endEdgei]])
- )
- {
- ++endEdgei;
- }
-
- if (edgei + 1 == endEdgei)
- {
- // Simplest case - a single connection
-
- newEdgeLabels[edgei] = oldEdgeLabels[sortIndices[edgei]];
- }
- else
- {
- multihit.clear();
-
- // Map from global edge to patch local edgeId
- for (label i = edgei; i < endEdgei; ++i)
- {
- const label patchEdgei = oldEdgeLabels[sortIndices[i]];
-
- // The edge in mesh numbering
- multihit.insert
- (
- patch().meshEdge(patchEdgei),
- patchEdgei
- );
- }
-
- if (multihit.size() != (endEdgei - edgei))
- {
- FatalErrorInFunction
- << "Could only hash " << multihit.size()
- << " edges from " << (endEdgei - edgei)
- << " ... indicates a non-manifold connection" << nl
- << multihit << nl
- << abort(FatalError);
- }
-
- const face& f = mesh().faces()[meshFacei];
-
- forAll(f, fedgei) // Note size() == nEdges()
- {
- edge e =
- (
- patchDef.owner()
- ? f.edge(fedgei) // Forward walk
- : f.rcEdge(fedgei) // Reverse walk
- );
-
- auto iter = multihit.find(e);
- if (iter.found())
- {
- newEdgeLabels[edgei++] = iter.val();
- multihit.erase(iter);
- if (multihit.empty())
- {
- break;
- }
- }
- }
-
- if (edgei != endEdgei)
- {
- FatalErrorInFunction
- << "Missed " << (edgei < endEdgei)
- << " edges for face: " << meshFacei
- << " ... indicates serious geometry issue" << nl
- << multihit << nl
- << abort(FatalError);
- }
- if (!multihit.empty())
- {
- FatalErrorInFunction
- << "Missed edges for face: " << meshFacei
- << " ... indicates serious geometry issue" << nl
- << multihit << nl
- << abort(FatalError);
- }
- }
- edgei = endEdgei;
- }
-
- patchDef.edgeLabels_.transfer(newEdgeLabels);
-}
-
+#include "foamVtkLineWriter.H"
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
@@ -399,276 +98,6 @@ Foam::PtrList<Foam::faPatch> Foam::faMesh::createOnePatch
}
-Foam::List<Foam::LabelledItem<Foam::edge>>
-Foam::faMesh::getBoundaryEdgePatchPairs() const
-{
- const polyBoundaryMesh& pbm = mesh().boundaryMesh();
-
- const label nInternalEdges = patch().nInternalEdges();
- const label nBoundaryEdges = patch().nBoundaryEdges();
-
- // Map edges (mesh numbering) back to a boundary index
- EdgeMap<label> edgeToBoundaryIndex(2*nBoundaryEdges);
-
- // Use labelled 'edge' for accounting of patch pairs
- patchPairInfoList bndEdgePatchPairs(nBoundaryEdges);
-
-
- // Pass 1:
- // - setup lookup (edge -> bnd index)
- // - add owner patch for each boundary edge
- {
- const SubList<labelList> bndEdgeToPatchFace
- (
- patch().edgeFaces(),
- patch().nBoundaryEdges(),
- patch().nInternalEdges()
- );
-
- for (label bndEdgei = 0; bndEdgei < nBoundaryEdges; ++bndEdgei)
- {
- edgeToBoundaryIndex.insert
- (
- patch().meshEdge(bndEdgei + nInternalEdges),
- edgeToBoundaryIndex.size()
- );
-
- // The attached patch face (there is only one):
- const label patchFacei = bndEdgeToPatchFace[bndEdgei][0];
-
- const label meshFacei = faceLabels_[patchFacei];
-
- const label patchId = pbm.whichPatch(meshFacei);
-
- bndEdgePatchPairs[bndEdgei].insert(patchId);
- }
- }
-
- // Pass 2:
- // - Add in first neighbour patch for the boundary edges
- // - examine all possible connecting neighbours
-
- findEdgePatchPairing
- (
- pbm,
- edgeToBoundaryIndex,
- bndEdgePatchPairs,
- edgeToBoundaryIndex.size() // Number of places to still fill
- );
-
-
- // Nothing dangling if running in serial - can return already
- if (!Pstream::parRun())
- {
- return bndEdgePatchPairs;
- }
-
- // In parallel need to check for "dangling" edges, which are finiteArea
- // boundary edges that only exist on one side of a proc boundary.
- // Eg, proc boundary coincides with the outer extent of the finiteArea
-
- const globalMeshData& globalData = mesh().globalData();
- const indirectPrimitivePatch& cpp = globalData.coupledPatch();
- const mapDistribute& map = globalData.globalEdgeSlavesMap();
-
- // Construct coupled edge usage with all data
- List<unsigned char> coupledEdgesUsed(map.constructSize(), 0u);
-
- forAll(cpp.edges(), coupledEdgei)
- {
- const auto iter =
- edgeToBoundaryIndex.cfind(cpp.meshEdge(coupledEdgei));
-
- // Used from this side or other other side
- coupledEdgesUsed[coupledEdgei] = (iter.found() ? 1 : 2);
- }
-
- // Save the original (pre-sync) coupling state
- const List<unsigned char> coupledEdgesOrig(coupledEdgesUsed);
-
- globalData.syncData
- (
- coupledEdgesUsed,
- globalData.globalEdgeSlaves(),
- globalData.globalEdgeTransformedSlaves(), // probably not used
- map,
- bitOrEqOp<unsigned char>()
- );
-
-
- // Check for one-sided edge coupling (coupled value == 3)
- // original == 1:
- // - coupled to a real finiteArea edge.
- // - receive a patch-pair value
- //
- // original == 2:
- // - a "dangled" edge. Information required for other procs.
- // - push a patch-pair value
-
- // Map edges (mesh numbering) back to a coupled index.
- // These are the edges to 'push' information for.
-
- EdgeMap<label> edgeToCoupledIndex;
-
- label nEdgesPull = 0;
-
- forAll(coupledEdgesUsed, coupledEdgei)
- {
- if (coupledEdgesUsed[coupledEdgei] == 3)
- {
- if (coupledEdgesOrig[coupledEdgei] == 1)
- {
- // Coupled side with finiteArea
- ++nEdgesPull;
- }
- else if (coupledEdgesOrig[coupledEdgei] == 2)
- {
- // Coupled side without finiteArea
- edgeToCoupledIndex.insert
- (
- cpp.meshEdge(coupledEdgei),
- coupledEdgei
- );
- }
- }
- }
-
- // Nothing to do - can return already
- if (returnReduce(edgeToCoupledIndex.empty(), andOp<bool>()))
- {
- return bndEdgePatchPairs;
- }
-
- // Data locations to pull
- labelList patchEdgeLabels(nEdgesPull);
- labelList coupledEdgeLabels(nEdgesPull);
-
- // Populate the locations
- {
- nEdgesPull = 0;
-
- forAll(cpp.edges(), coupledEdgei)
- {
- if
- (
- coupledEdgesUsed[coupledEdgei] == 3
- && coupledEdgesOrig[coupledEdgei] == 1
- )
- {
- // Pull this edge
- const auto iter =
- edgeToBoundaryIndex.cfind(cpp.meshEdge(coupledEdgei));
-
- if (iter.found())
- {
- patchEdgeLabels[nEdgesPull] = iter.val();
- coupledEdgeLabels[nEdgesPull] = coupledEdgei;
- ++nEdgesPull;
- }
- else
- {
- // Should be impossible to fail here
- FatalErrorInFunction
- << "Failed on second lookup of "
- << cpp.meshEdge(coupledEdgei) << nl
- << abort(FatalError);
- }
- }
- }
-
- if (nEdgesPull != coupledEdgeLabels.size())
- {
- FatalErrorInFunction
- << "Failed lookup of some coupled edges" << nl
- << abort(FatalError);
- }
- }
-
- //- Construct edge sync with all data
- patchPairInfoList cppEdgeData(map.constructSize());
-
- // Fill in for 'push' locations. Only really interested in the owner
- // (corresponds to non-proc connection), but grab everything
- findEdgePatchPairing
- (
- pbm,
- edgeToCoupledIndex,
- cppEdgeData,
- 2*edgeToCoupledIndex.size() // Accept both sides?
- );
-
-
- const label nNonProcessor = pbm.nNonProcessor();
-
- // Adjust patch information to reflect dangling patch neighbour
- // Tag with -(value+2)
- forAllConstIters(edgeToCoupledIndex, iter)
- {
- const edge& e = iter.key();
- const label coupledEdgei = iter.val();
-
- patchPairInfo& pairing = cppEdgeData[coupledEdgei];
- const label ownerPatchId = pairing.first();
-
- // Some sanity checks
- if (ownerPatchId < 0)
- {
- FatalErrorInFunction
- << "Error finding dangling edge at "
- << cpp.points()[e.first()] << ' '
- << cpp.points()[e.second()] << nl
- << abort(FatalError);
- }
- else if (ownerPatchId >= nNonProcessor)
- {
- FatalErrorInFunction
- << "Cannot handle edge on processor-processor connection at "
- << cpp.points()[e.first()] << ' '
- << cpp.points()[e.second()] << nl
- << abort(FatalError);
- }
-
- combineDanglingEdge::setDangling(pairing, ownerPatchId);
-
- // TBD:
- // may wish to remember the corresponding proc number,
- // if we wish to bridge across 'fan-like' connections.
- //
- // pairing.setIndex(-(Pstream::myProcNo() + 2));
- }
-
-
- // Synchronize edge information
-
- const combineDanglingEdge edgeCombineOp(nNonProcessor);
-
- globalMeshData::syncData
- (
- cppEdgeData,
- globalData.globalEdgeSlaves(),
- globalData.globalEdgeTransformedSlaves(), // probably not used
- map,
- edgeCombineOp
- );
-
-
- // Combine back from pushed cpp-edge data
-
- forAll(patchEdgeLabels, i)
- {
- patchPairInfo& pairing = bndEdgePatchPairs[patchEdgeLabels[i]];
- const patchPairInfo& other = cppEdgeData[coupledEdgeLabels[i]];
-
- edgeCombineOp(pairing, other);
-
- // Resolve special tagging
- combineDanglingEdge::correct(pairing);
- }
-
- return bndEdgePatchPairs;
-}
-
-
Foam::PtrList<Foam::faPatch> Foam::faMesh::createPatchList
(
const dictionary& bndDict,
@@ -679,7 +108,7 @@ Foam::PtrList<Foam::faPatch> Foam::faMesh::createPatchList
const polyBoundaryMesh& pbm = mesh().boundaryMesh();
// Transcribe into patch definitions
- DynamicList<faPatchData> faPatchDefs(bndDict.size() + 4);
+ DynamicList<faPatchData> faPatchDefs(bndDict.size() + 8);
for (const entry& dEntry : bndDict)
{
if (!dEntry.isDict())
@@ -734,8 +163,10 @@ Foam::PtrList<Foam::faPatch> Foam::faMesh::createPatchList
patchDef.type_ = "empty";
}
+ label nWarnUndefinedPatch(5);
+
// Placeholder for any undefined edges
- const label undefPatchId = faPatchDefs.size();
+ const label undefPatchIndex = faPatchDefs.size();
{
faPatchDefs.append(faPatchData());
@@ -745,198 +176,321 @@ Foam::PtrList<Foam::faPatch> Foam::faMesh::createPatchList
if (defaultPatchDefinition)
{
- (*defaultPatchDefinition).readIfPresent("name", patchDef.name_);
+ if
+ (
+ (*defaultPatchDefinition).readIfPresent("name", patchDef.name_)
+ )
+ {
+ // Suppress warnings if defaultPatch name was specified
+ // - probably means we want to use this mechanism
+ nWarnUndefinedPatch = 0;
+ }
(*defaultPatchDefinition).readIfPresent("type", patchDef.type_);
}
}
// ----------------------------------------------------------------------
- const label nInternalEdges = patch().nInternalEdges();
- const label nBoundaryEdges = patch().nBoundaryEdges();
+ // Get edge connections (in globally consistent ordering)
+ List<Pair<patchTuple>> bndEdgeConnections
+ (
+ this->getBoundaryEdgeConnections()
+ );
+
+ // Update accordingly
+ this->setBoundaryConnections(bndEdgeConnections);
+
- patchPairInfoList bndEdgePatchPairs(this->getBoundaryEdgePatchPairs());
+ // Lookup of patchDef for each connection. Initially all -1 (unvisited)
+ labelList patchDefLookup(bndEdgeConnections.size(), -1);
- labelList bndEdgeFaPatchIDs(nBoundaryEdges, -1);
+ Map<labelHashSet> procConnections;
+ labelHashSet patchDefsUsed;
- for (label bndEdgei = 0; bndEdgei < nBoundaryEdges; ++bndEdgei)
+ label nBadEdges(0);
+ labelHashSet badEdges(2*bndEdgeConnections.size());
+
+ forAll(bndEdgeConnections, connecti)
{
- const patchPairInfo& patchPair = bndEdgePatchPairs[bndEdgei];
+ const Pair<patchTuple>& connection = bndEdgeConnections[connecti];
+ const auto& a = connection.first();
+ const auto& b = connection.second();
+
+ edge patchPair;
+
+ if (a.is_finiteArea())
+ {
+ if (b.is_finiteArea())
+ {
+ // A processor-processor connection
+ if (a.procNo() == b.procNo())
+ {
+ FatalErrorInFunction
+ << "Processor-processor addressing error:" << nl
+ << "Both connections have the same processor: "
+ << a.procNo() << nl
+ << abort(FatalError);
+ }
+ else if (a.is_localProc())
+ {
+ procConnections(b.procNo()).insert(connecti);
+ }
+ else
+ {
+ procConnections(a.procNo()).insert(connecti);
+ }
+
+ continue;
+ }
+ else if (a.is_localProc())
+ {
+ patchPair.first() = a.realPatchi();
+ patchPair.second() = b.realPatchi();
+ }
+ }
+ else if (b.is_finiteArea() && b.is_localProc())
+ {
+ patchPair.first() = b.realPatchi();
+ patchPair.second() = a.realPatchi();
+ }
+
// Find first definition with a matching neighbour and
// possibly with a matching owner.
- label bestPatchi = -1;
+ label bestPatchDefi = -1;
+
+ const label nPatchDefs = (patchPair.valid() ? faPatchDefs.size() : 0);
- for (label patchi = 0; patchi < faPatchDefs.size(); ++patchi)
+ for (label patchDefi = 0; patchDefi < nPatchDefs; ++patchDefi)
{
- const int match = faPatchDefs[patchi].matchPatchPair(patchPair);
+ const int match = faPatchDefs[patchDefi].matchPatchPair(patchPair);
if (match == 3)
{
- // Match (owner/neighbour) - done!
- bestPatchi = patchi;
+ // Exact match (owner/neighbour) - done!
+ bestPatchDefi = patchDefi;
break;
}
- else if (match == 2 && bestPatchi < 0)
+ else if (match == 2 && bestPatchDefi < 0)
{
// Match (neighbour) - keep looking for exact match
- bestPatchi = patchi;
+ bestPatchDefi = patchDefi;
}
}
- bndEdgeFaPatchIDs[bndEdgei] = bestPatchi;
+ if (bestPatchDefi < 0)
+ {
+ bestPatchDefi = undefPatchIndex; // Missed?
+ }
+
+ patchDefLookup[connecti] = bestPatchDefi;
+ patchDefsUsed.insert(bestPatchDefi);
}
+ // Remove undefPatchIndex if not actually needed
+ if (!returnReduce(patchDefsUsed.found(undefPatchIndex), orOp<bool>()))
+ {
+ faPatchDefs.remove(undefPatchIndex);
+ }
+ else
+ {
+ patchDefsUsed.insert(undefPatchIndex); // Parallel consistency
- // Extract which edges map to which patch
- // and set edgeLabels for each faPatch
+ // Report locations of undefined edges
- DynamicList<label> selectEdges(bndEdgeFaPatchIDs.size());
+ badEdges.clear();
+ forAll(patchDefLookup, connecti)
+ {
+ if (patchDefLookup[connecti] == undefPatchIndex)
+ {
+ const auto& connection = bndEdgeConnections[connecti];
- for (label patchi = 0; patchi < faPatchDefs.size(); ++patchi)
- {
- auto& patchDef = faPatchDefs[patchi];
+ const auto& a = connection.first();
+ const auto& b = connection.second();
- selectEdges.clear();
+ if (a.is_localProc() && a.is_finiteArea())
+ {
+ badEdges.insert(a.patchEdgei());
+ }
+ else if (b.is_localProc() && b.is_finiteArea())
+ {
+ badEdges.insert(b.patchEdgei());
+ }
+
+ if (badEdges.size() <= nWarnUndefinedPatch)
+ {
+ Pout<< "Undefined connection: "
+ << "(patch:" << a.realPatchi()
+ << " face:" << a.meshFacei()
+ << ") and (patch:" << b.realPatchi()
+ << " face:" << b.meshFacei() << ") connects: "
+ << pbm[a.realPatchi()].name() << " to "
+ << pbm[b.realPatchi()].name() << nl;
+ }
+ }
+ }
- forAll(bndEdgeFaPatchIDs, bndEdgei)
+ if ((nBadEdges = returnReduce(badEdges.size(), sumOp<label>())) > 0)
{
- if (bndEdgeFaPatchIDs[bndEdgei] == patchi)
+ // Report directly as Info, not InfoInFunction
+ // since it can also be an expected result when
+ // nWarnUndefinedPatch == 0
+ Info<< "Had " << nBadEdges << '/'
+ << returnReduce(patch().nBoundaryEdges(), sumOp<label>())
+ << " undefined edge connections, added to defaultPatch: "
+ << faPatchDefs[undefPatchIndex].name_ << nl;
+
+ if (nWarnUndefinedPatch)
{
- selectEdges.append(bndEdgei + nInternalEdges);
+ edgeList dumpEdges(patch().edges(), badEdges.sortedToc());
+
+ vtk::lineWriter writer
+ (
+ patch().localPoints(),
+ dumpEdges,
+ fileName
+ (
+ mesh().time().globalPath()
+ / ("faMesh-construct.undefEdges")
+ )
+ );
+
+ writer.writeGeometry();
+
+ // CellData
+ writer.beginCellData();
+ writer.writeProcIDs();
+
+ InfoInFunction
+ << "(debug) wrote " << writer.output().name() << nl;
}
}
-
- patchDef.edgeLabels_ = selectEdges;
}
- // Check for undefined edges
- selectEdges.clear();
-
- forAll(bndEdgeFaPatchIDs, bndEdgei)
+ // Create processor-processor definitions
+ Map<label> procToDefLookup(2*procConnections.size());
{
- if (bndEdgeFaPatchIDs[bndEdgei] == -1)
+ faPatchDefs.reserve(faPatchDefs.size() + procConnections.size());
+
+ for (const label otherProci : procConnections.sortedToc())
{
- selectEdges.append(bndEdgei + nInternalEdges);
- }
- }
+ const label patchDefi = faPatchDefs.size();
+ procToDefLookup.insert(otherProci, patchDefi);
- // Save the information
- faPatchDefs[undefPatchId].edgeLabels_ = selectEdges;
+ // Add entry
+ faPatchDefs.append(faPatchData());
+ auto& patchDef = faPatchDefs.last();
- bool hasUndefined = returnReduce(!selectEdges.empty(), orOp<bool>());
+ patchDef.assign_coupled(Pstream::myProcNo(), otherProci);
+ }
+ }
- if (hasUndefined)
- {
- // The initial edges to consider
- const labelList& undefinedEdges =
- faPatchDefs[undefPatchId].edgeLabels_;
- // Check for edges that butt against a processor (or other) patch
+ // Extract which edges map to which patch
+ // and set edgeLabels for each faPatch
- labelList edgeNbrPolyPatch(undefinedEdges.size(), -1);
- forAll(edgeNbrPolyPatch, edgei)
- {
- const label patchEdgei = undefinedEdges[edgei];
- const label bndEdgei = (patchEdgei - nInternalEdges);
+ DynamicList<label> selectEdges(bndEdgeConnections.size());
+ label nOffProcessorEdges = 0;
- edgeNbrPolyPatch[edgei] = bndEdgePatchPairs[bndEdgei].second();
- }
+ for (const label patchDefi : patchDefsUsed.sortedToc())
+ {
+ auto& patchDef = faPatchDefs[patchDefi];
+ selectEdges.clear();
- // Categorize as processor/non-processor associations
- labelHashSet procPatchIDs;
- labelHashSet nonProcPatchIDs;
+ // Find the corresponding entries
+ // and extract the patchEdgeId
- for (const label polyPatchID : edgeNbrPolyPatch)
+ forAll(patchDefLookup, connecti)
{
- if (polyPatchID == -1)
- {
- nonProcPatchIDs.insert(polyPatchID);
- }
- else if
- (
- !nonProcPatchIDs.found(polyPatchID)
- && !procPatchIDs.found(polyPatchID)
- )
+ if (patchDefLookup[connecti] == patchDefi)
{
- if (isA<processorPolyPatch>(pbm[polyPatchID]))
+ const auto& a = bndEdgeConnections[connecti].first();
+ const auto& b = bndEdgeConnections[connecti].second();
+
+ if (a.is_localProc() && a.is_finiteArea())
{
- procPatchIDs.insert(polyPatchID);
+ selectEdges.append(a.patchEdgei());
+ }
+ else if (b.is_localProc() && b.is_finiteArea())
+ {
+ selectEdges.append(b.patchEdgei());
}
else
{
- nonProcPatchIDs.insert(polyPatchID);
+ FatalErrorInFunction
+ << "Error in programming logic" << nl
+ << abort(FatalError);
}
- }
- }
- // Select by processor association
- for (const label polyPatchID : procPatchIDs.sortedToc())
- {
- selectEdges.clear();
-
- forAll(edgeNbrPolyPatch, edgei)
- {
- if (edgeNbrPolyPatch[edgei] == polyPatchID)
+ if (a.is_localProc() != b.is_localProc())
{
- selectEdges.append(undefinedEdges[edgei]);
+ ++nOffProcessorEdges;
}
+
+ patchDefLookup[connecti] = -2; // Mark as handled
}
+ }
- faPatchDefs.append(faPatchData());
+ // Remove any cosmetic sorting artifacts from off-processor
+ // termination by doing using a regular sort here.
- auto& patchDef = faPatchDefs.last();
- patchDef.name_ = pbm[polyPatchID].name();
- patchDef.type_ = processorFaPatch::typeName;
- patchDef.neighPolyPatchId_ = polyPatchID; // Needed for reorder
+ Foam::sort(selectEdges);
+ patchDef.edgeLabels_ = selectEdges;
+ }
- const auto* ppp = isA<processorPolyPatch>(pbm[polyPatchID]);
- if (ppp)
- {
- patchDef.ownerProcId_ = ppp->myProcNo();
- patchDef.neighProcId_ = ppp->neighbProcNo();
- }
+ if (debug)
+ {
+ Pout<< "Had " << nOffProcessorEdges
+ << " patch edges connected off-processor" << endl;
+
+ InfoInFunction
+ << "Total "
+ << returnReduce(nOffProcessorEdges, sumOp<label>())
+ << " patch edges connected off-processor" << endl;
+ }
- patchDef.edgeLabels_ = selectEdges;
- }
+ // Processor patches
+ for (const label otherProci : procToDefLookup.sortedToc())
+ {
+ const label patchDefi = procToDefLookup[otherProci];
- // Check for any remaining undefined edges
+ auto& patchDef = faPatchDefs[patchDefi];
selectEdges.clear();
- // Simply grab any/all (don't worry about which patch)
- if (!nonProcPatchIDs.empty())
+ // Find the corresponding entries
+ // and extract the patchEdgeId
+
+ for (const label connecti : procConnections(otherProci).sortedToc())
{
- forAll(edgeNbrPolyPatch, edgei)
+ const auto& connection = bndEdgeConnections[connecti];
+ const auto& a = connection.first();
+ const auto& b = connection.second();
+
+ if (a.is_localProc())
{
- const label polyPatchID = edgeNbrPolyPatch[edgei];
- if (nonProcPatchIDs.found(polyPatchID))
- {
- selectEdges.append(undefinedEdges[edgei]);
- }
+ selectEdges.append(a.patchEdgei());
+ }
+ else if (b.is_localProc())
+ {
+ selectEdges.append(b.patchEdgei());
+ }
+ else
+ {
+ FatalErrorInFunction
+ << "Error in programming logic" << nl
+ << abort(FatalError);
}
- }
- // Complete the information
- faPatchDefs[undefPatchId].edgeLabels_ = selectEdges;
+ patchDefLookup[connecti] = -2; // Mark as handled
+ }
- hasUndefined = returnReduce(!selectEdges.empty(), orOp<bool>());
- }
+ // The edge order is guaranteed to be consistent from the original
+ // getBoundaryEdgeConnections() - sorted by proc/patch/edge
- // Remove unnecessary entry
- if (!hasUndefined)
- {
- faPatchDefs.remove(undefPatchId);
+ patchDef.edgeLabels_ = selectEdges;
}
- for (auto& patchDef : faPatchDefs)
- {
- if (patchDef.coupled())
- {
- reorderProcEdges(patchDef, bndEdgePatchPairs);
- patchDef.neighPolyPatchId_ = -1; // No lookup of neighbour faces
- }
- }
// Now convert list of definitions to list of patches
@@ -962,6 +516,23 @@ Foam::PtrList<Foam::faPatch> Foam::faMesh::createPatchList
++nPatches;
}
+
+ if (debug > 1)
+ {
+ label nTotal = 0;
+ Pout<< "Created new finiteArea patches:" << nl;
+ for (const faPatch& p : newPatches)
+ {
+ Pout<< " size: " << p.size()
+ << " name:" << p.name()
+ << " type:" << p.type() << nl;
+ nTotal += p.size();
+ }
+
+ Pout<< "addressed: " << nTotal
+ << '/' << patch().nBoundaryEdges() << " edges" << endl;
+ }
+
return newPatches;
}
diff --git a/src/finiteArea/faMesh/faMeshTopology.C b/src/finiteArea/faMesh/faMeshTopology.C
new file mode 100644
index 00000000000..179054ac0b6
--- /dev/null
+++ b/src/finiteArea/faMesh/faMeshTopology.C
@@ -0,0 +1,803 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | www.openfoam.com
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+ Copyright (C) 2021 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "faMesh.H"
+#include "globalMeshData.H"
+#include "indirectPrimitivePatch.H"
+#include "edgeHashes.H"
+#include "foamVtkLineWriter.H"
+#include "foamVtkIndPatchWriter.H"
+#include "foamVtkUIndPatchWriter.H"
+
+// * * * * * * * * * * * * * * * Local Functions * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// Print out edges as point pairs
+template<class PatchType>
+static void printPatchEdges
+(
+ Ostream& os,
+ const PatchType& p,
+ const labelList& edgeIds,
+ label maxOutput = 10
+)
+{
+ label nOutput = 0;
+
+ for (const label patchEdgei : edgeIds)
+ {
+ const edge e(p.meshEdge(patchEdgei));
+
+ os << " "
+ << p.points()[e.first()] << ' '
+ << p.points()[e.second()] << nl;
+
+ ++nOutput;
+ if (maxOutput > 0 && nOutput >= maxOutput)
+ {
+ os << " ... suppressing further output" << nl;
+ break;
+ }
+ }
+}
+
+} // End namespace Foam
+
+
+// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
+
+Foam::List<Foam::Pair<Foam::faMesh::patchTuple>>
+Foam::faMesh::getBoundaryEdgeConnections() const
+{
+ const polyBoundaryMesh& pbm = mesh().boundaryMesh();
+
+ const label nNonProcessor = pbm.nNonProcessor();
+
+ const label nInternalEdges = patch().nInternalEdges();
+ const label nBoundaryEdges = patch().nBoundaryEdges();
+
+ // The output result:
+ List<Pair<patchTuple>> bndEdgeConnections(nBoundaryEdges);
+
+ // Map edges (mesh numbering) back to a boundary index
+ EdgeMap<label> edgeToBoundaryIndex(2*nBoundaryEdges);
+
+ label nBadEdges(0);
+ labelHashSet badEdges(2*nBoundaryEdges);
+
+ {
+ // Local collection structure for accounting of patch pairs.
+ // Based on 'edge' which has some hash-like insertion properties
+ // that are useful.
+ struct patchPairingType : public Foam::edge
+ {
+ label patchEdgei_ = -1;
+ label meshFacei_ = -1;
+
+ void clear()
+ {
+ Foam::edge::clear();
+ patchEdgei_ = -1;
+ meshFacei_ = -1;
+ }
+ };
+
+ List<patchPairingType> patchPairings(nBoundaryEdges);
+
+ DebugInFunction
+ << "Determining required boundary edge connections, "
+ << "resolving locally attached boundary edges." << endl;
+
+
+ // Pass 1:
+ // - setup lookup (edge -> bnd index)
+ // - add owner patch for each boundary edge
+ for (label bndEdgei = 0; bndEdgei < nBoundaryEdges; ++bndEdgei)
+ {
+ const label patchEdgei = (bndEdgei + nInternalEdges);
+
+ edgeToBoundaryIndex.insert
+ (
+ patch().meshEdge(patchEdgei),
+ bndEdgei
+ );
+
+ // The attached patch face. Should only be one!
+ const labelList& edgeFaces = patch().edgeFaces()[patchEdgei];
+
+ if (edgeFaces.size() != 1)
+ {
+ badEdges.insert(patchEdgei);
+ continue;
+ }
+
+ const label patchFacei = edgeFaces[0];
+ const label meshFacei = faceLabels_[patchFacei];
+ const label bndFacei = (meshFacei - mesh().nInternalFaces());
+
+ /// const label patchId = pbm.whichPatch(meshFacei);
+ const label patchId = pbm.patchID()[bndFacei];
+
+ // Primary bookkeeping
+ {
+ auto& tuple = bndEdgeConnections[bndEdgei].first();
+
+ tuple.procNo(Pstream::myProcNo());
+ tuple.faPatchi(patchId); // Tag as finiteArea patch
+ tuple.patchEdgei(patchEdgei);
+ tuple.meshFacei(meshFacei);
+ }
+
+ // Temporary local bookkeeping
+ {
+ auto& pairing = patchPairings[bndEdgei];
+
+ pairing.clear(); // invalidate
+ pairing.insert(patchId); // 'hash' into first location
+ pairing.patchEdgei_ = patchEdgei;
+ pairing.meshFacei_ = meshFacei;
+ }
+ }
+
+ if ((nBadEdges = returnReduce(badEdges.size(), sumOp<label>())) > 0)
+ {
+ edgeList dumpEdges(patch().edges(), badEdges.sortedToc());
+
+ vtk::lineWriter writer
+ (
+ patch().localPoints(),
+ dumpEdges,
+ fileName
+ (
+ mesh().time().globalPath()
+ / ("faMesh-construct.nonManifoldEdges")
+ )
+ );
+
+ writer.writeGeometry();
+
+ // CellData
+ writer.beginCellData();
+ writer.writeProcIDs();
+
+ FatalErrorInFunction
+ << "Boundary edges not singly connected: "
+ << nBadEdges << '/' << nBoundaryEdges << nl;
+
+ printPatchEdges
+ (
+ FatalError,
+ patch(),
+ badEdges.sortedToc()
+ );
+
+ InfoInFunction
+ << "(debug) wrote " << writer.output().name() << nl;
+
+ FatalError << abort(FatalError);
+ }
+ badEdges.clear();
+
+ // Pass 2:
+ // Add in first connecting neighbour patch for the boundary edges.
+ // Need to examine all possibly connecting (non-processor) neighbours,
+ // but only need to check their boundary edges.
+
+ label nMissing = patchPairings.size();
+
+ for (label patchi = 0; patchi < nNonProcessor; ++patchi)
+ {
+ if (!nMissing) break; // Early exit
+
+ const polyPatch& pp = pbm[patchi];
+
+ // Check boundary edges
+ for
+ (
+ label patchEdgei = pp.nInternalEdges();
+ patchEdgei < pp.nEdges();
+ ++patchEdgei
+ )
+ {
+ const label bndEdgei =
+ edgeToBoundaryIndex.lookup(pp.meshEdge(patchEdgei), -1);
+
+ if (bndEdgei != -1)
+ {
+ // Has a matching owner boundary edge
+
+ auto& pairing = patchPairings[bndEdgei];
+
+ // Add neighbour (patchId, patchEdgei, meshFacei)
+ // 'hash' into the second location
+ // which does not insert the same value twice
+ if (pairing.insert(patchi))
+ {
+ // The attached patch face. Should only be one!
+ const labelList& edgeFaces = pp.edgeFaces()[patchEdgei];
+
+ if (edgeFaces.size() != 1)
+ {
+ pairing.erase(patchi);
+ badEdges.insert(badEdges.size());
+ continue;
+ }
+
+ const label patchFacei = edgeFaces[0];
+ const label meshFacei = patchFacei + pp.start();
+
+ // The neighbour information
+ pairing.patchEdgei_ = patchEdgei;
+ pairing.meshFacei_ = meshFacei;
+
+ --nMissing;
+ if (!nMissing) break; // Early exit
+ }
+ }
+ }
+ }
+
+ if ((nBadEdges = returnReduce(badEdges.size(), sumOp<label>())) > 0)
+ {
+ FatalErrorInFunction
+ << "Had " << nBadEdges
+ << " boundary edges with missing or multiple edge connections"
+ << abort(FatalError);
+ }
+
+ // Combine local bookkeeping into final list
+ badEdges.clear();
+ for (label bndEdgei = 0; bndEdgei < nBoundaryEdges; ++bndEdgei)
+ {
+ const auto& pairing = patchPairings[bndEdgei];
+ const label nbrPatchi = pairing.second();
+ const label nbrPatchEdgei = pairing.patchEdgei_;
+ const label nbrMeshFacei = pairing.meshFacei_;
+
+ if (nbrMeshFacei >= 0)
+ {
+ // Add into primary bookkeeping
+ auto& tuple = bndEdgeConnections[bndEdgei].second();
+
+ tuple.procNo(Pstream::myProcNo());
+ tuple.patchi(nbrPatchi);
+ tuple.patchEdgei(nbrPatchEdgei);
+ tuple.meshFacei(nbrMeshFacei);
+ }
+ else if (!Pstream::parRun())
+ {
+ badEdges.insert(nInternalEdges + bndEdgei);
+ }
+ }
+ }
+
+
+ // ~~~~~~
+ // Serial - can return already
+ // ~~~~~~
+ if (!Pstream::parRun())
+ {
+ // Verbose report of missing edges - in serial
+
+ nBadEdges = badEdges.size();
+ if (nBadEdges)
+ {
+ edgeList dumpEdges(patch().edges(), badEdges.sortedToc());
+
+ vtk::lineWriter writer
+ (
+ patch().localPoints(),
+ dumpEdges,
+ fileName
+ (
+ mesh().time().globalPath()
+ / ("faMesh-construct.invalidEdges")
+ )
+ );
+
+ writer.writeGeometry();
+
+ // CellData
+ writer.beginCellData();
+ writer.writeProcIDs();
+
+ FatalErrorInFunction
+ << "Boundary edges with missing/invalid neighbours: "
+ << nBadEdges << '/' << nBoundaryEdges << nl;
+
+ printPatchEdges
+ (
+ FatalError,
+ patch(),
+ badEdges.sortedToc()
+ );
+
+ InfoInFunction
+ << "(debug) wrote " << writer.output().name() << nl;
+
+ FatalError << abort(FatalError);
+ }
+
+ // Globally consistent ordering
+ patchTuple::sort(bndEdgeConnections);
+ return bndEdgeConnections;
+ }
+
+
+ // ~~~~~~~~
+ // Parallel
+ // ~~~~~~~~
+
+ DebugInFunction
+ << "Creating global coupling data" << endl;
+
+ const globalMeshData& globalData = mesh().globalData();
+ const indirectPrimitivePatch& cpp = globalData.coupledPatch();
+ const mapDistribute& map = globalData.globalEdgeSlavesMap();
+ const label nCoupledEdges = cpp.nEdges();
+
+ // Construct coupled edge usage with all data
+ List<bool> coupledEdgesUsed(map.constructSize(), false);
+
+ // Markup finiteArea boundary edges that are coupled across processors
+ for (label cppEdgei = 0; cppEdgei < nCoupledEdges; ++cppEdgei)
+ {
+ coupledEdgesUsed[cppEdgei] =
+ edgeToBoundaryIndex.found(cpp.meshEdge(cppEdgei));
+ }
+
+ DebugInFunction
+ << "Starting sync of boundary edge topology" << endl;
+
+ globalMeshData::syncData
+ (
+ coupledEdgesUsed,
+ globalData.globalEdgeSlaves(),
+ globalData.globalEdgeTransformedSlaves(), // probably not used
+ map,
+ orEqOp<bool>()
+ );
+
+ if (debug)
+ {
+ label nAttached = 0;
+ for (label cppEdgei = 0; cppEdgei < nCoupledEdges; ++cppEdgei)
+ {
+ if (coupledEdgesUsed[cppEdgei])
+ {
+ ++nAttached;
+ }
+ }
+
+ InfoInFunction
+ << "Approx "
+ << returnReduce(nAttached, sumOp<label>())
+ << " connected boundary edges (total, some duplicates)" << endl;
+ }
+
+ // Combine information
+
+ // Normally 0 or 2 connections
+ List<DynamicList<patchTuple, 2>> gatheredConnections(map.constructSize());
+
+ // Map edges (mesh numbering) back to a coupled index in use
+ EdgeMap<label> edgeToCoupledIndex(2*nCoupledEdges);
+
+ // Pass 1
+ // Look for attached boundary edges
+ // - boundary edges from this side go into gathered connections
+ // - boundary edges connected from the other side are noted for later
+
+ for (label cppEdgei = 0; cppEdgei < nCoupledEdges; ++cppEdgei)
+ {
+ if (coupledEdgesUsed[cppEdgei])
+ {
+ const edge meshEdge(cpp.meshEdge(cppEdgei));
+
+ const label bndEdgei =
+ edgeToBoundaryIndex.lookup(meshEdge, -1);
+
+ if (bndEdgei != -1)
+ {
+ // A boundary finiteEdge edge (known from this side)
+
+ auto& gathered = gatheredConnections[cppEdgei];
+ gathered.setCapacity(2);
+ gathered.resize(1);
+ auto& tuple = gathered.last();
+
+ tuple = bndEdgeConnections[bndEdgei].first();
+ }
+ else
+ {
+ // Boundary edge connected from the other side
+ // - mark for it to be added later
+
+ edgeToCoupledIndex.insert(meshEdge, cppEdgei);
+ }
+ }
+ }
+
+ // Pass 2
+ // - add previously noted boundary edges (connected from other side)
+ // into gathered connections
+
+ badEdges.clear();
+ for (label patchi = 0; patchi < nNonProcessor; ++patchi)
+ {
+ if (edgeToCoupledIndex.empty()) break; // Early exit
+
+ const polyPatch& pp = pbm[patchi];
+
+ // Check boundary edges
+ for
+ (
+ label patchEdgei = pp.nInternalEdges();
+ patchEdgei < pp.nEdges();
+ ++patchEdgei
+ )
+ {
+ const edge meshEdge(pp.meshEdge(patchEdgei));
+
+ const label cppEdgei =
+ edgeToCoupledIndex.lookup(meshEdge, -1);
+
+ if (cppEdgei != -1)
+ {
+ // A known connection
+
+ // The attached patch face. Should only be one!
+ const labelList& edgeFaces = pp.edgeFaces()[patchEdgei];
+
+ if (edgeFaces.size() != 1)
+ {
+ badEdges.insert(cppEdgei);
+ continue;
+ }
+
+ const label patchFacei = edgeFaces[0];
+ const label meshFacei = patchFacei + pp.start();
+
+ auto& gathered = gatheredConnections[cppEdgei];
+ gathered.setCapacity(2);
+ gathered.resize(1);
+ auto& tuple = gathered.last();
+
+ tuple.procNo(Pstream::myProcNo());
+ tuple.patchi(patchi);
+ tuple.patchEdgei(patchEdgei);
+ tuple.meshFacei(meshFacei);
+
+ // Do not consider again
+ edgeToCoupledIndex.erase(meshEdge);
+
+ if (edgeToCoupledIndex.empty()) break; // Early exit
+ }
+ }
+ }
+
+ if ((nBadEdges = returnReduce(badEdges.size(), sumOp<label>())) > 0)
+ {
+ FatalErrorInFunction
+ << "Had " << nBadEdges << " coupled boundary edges"
+ << " with missing or multiple edge connections"
+ << abort(FatalError);
+ }
+
+ DebugInFunction
+ << "Starting sync of boundary edge information" << endl;
+
+ globalMeshData::syncData
+ (
+ gatheredConnections,
+ globalData.globalEdgeSlaves(),
+ globalData.globalEdgeTransformedSlaves(), // probably not used
+ map,
+ ListOps::appendEqOp<patchTuple>()
+ );
+
+
+ DebugInFunction
+ << "Collating sync information" << endl;
+
+ // Pick out gathered connections and add into primary bookkeeping
+ for (label cppEdgei = 0; cppEdgei < nCoupledEdges; ++cppEdgei)
+ {
+ const auto& gathered = gatheredConnections[cppEdgei];
+
+ const label bndEdgei =
+ edgeToBoundaryIndex.lookup(cpp.meshEdge(cppEdgei), -1);
+
+ if
+ (
+ // A boundary finiteEdge edge (known from this side)
+ bndEdgei != -1
+
+ // Gathered a one-to-one connection
+ && gathered.size() == 2
+ )
+ {
+ const auto& a = gathered[0];
+ const auto& b = gathered[1];
+
+ // Copy second side of connection
+ auto& connection = bndEdgeConnections[bndEdgei];
+
+ connection.second() = (connection.first() == b) ? a : b;
+ }
+ }
+
+ // Check missing/invalid
+ badEdges.clear();
+ for (label bndEdgei = 0; bndEdgei < nBoundaryEdges; ++bndEdgei)
+ {
+ const auto& connection = bndEdgeConnections[bndEdgei];
+
+ if (!connection.second().valid())
+ {
+ badEdges.insert(nInternalEdges + bndEdgei);
+ }
+ }
+
+
+ if (debug & 8)
+ {
+ // Boundary edges
+ {
+ vtk::lineWriter writer
+ (
+ patch().localPoints(),
+ patch().boundaryEdges(),
+ fileName
+ (
+ mesh().time().globalPath()
+ / ("faMesh-construct.boundaryEdges")
+ )
+ );
+
+ writer.writeGeometry();
+
+ // CellData
+ writer.beginCellData();
+ writer.writeProcIDs();
+
+ // For each boundary edge - the associate neighbour patch
+ labelList neighProc(nBoundaryEdges);
+ labelList neighPatch(nBoundaryEdges);
+ for (label bndEdgei = 0; bndEdgei < nBoundaryEdges; ++bndEdgei)
+ {
+ const auto& connection = bndEdgeConnections[bndEdgei];
+
+ neighProc[bndEdgei] = connection.second().procNo();
+ neighPatch[bndEdgei] = connection.second().patchi();
+ }
+
+ writer.write("neighProc", neighProc);
+ writer.write("neighPatch", neighPatch);
+ }
+
+ // finiteArea
+ {
+ vtk::uindirectPatchWriter writer
+ (
+ patch(),
+ fileName
+ (
+ mesh().time().globalPath()
+ / ("faMesh-construct.faPatch")
+ )
+ );
+
+ writer.writeGeometry();
+
+ // CellData
+ writer.beginCellData();
+ writer.writeProcIDs();
+ }
+ }
+
+ // Verbose report of missing edges
+ if ((nBadEdges = returnReduce(badEdges.size(), sumOp<label>())) > 0)
+ {
+ edgeList dumpEdges(patch().edges(), badEdges.sortedToc());
+
+ vtk::lineWriter writer
+ (
+ patch().localPoints(),
+ dumpEdges,
+ fileName
+ (
+ mesh().time().globalPath()
+ / ("faMesh-construct.invalidEdges")
+ )
+ );
+
+ writer.writeGeometry();
+
+ // CellData
+ writer.beginCellData();
+ writer.writeProcIDs();
+
+ FatalErrorInFunction
+ << "Boundary edges with missing/invalid neighbours: "
+ << nBadEdges << '/' << nBoundaryEdges << nl;
+
+ printPatchEdges
+ (
+ FatalError,
+ patch(),
+ badEdges.sortedToc()
+ );
+
+ InfoInFunction
+ << "(debug) wrote " << writer.output().name() << nl;
+
+ FatalError << abort(FatalError);
+ }
+
+
+ // Globally consistent ordering
+ patchTuple::sort(bndEdgeConnections);
+
+ DebugInFunction
+ << "Return sorted list of boundary connections" << endl;
+
+ return bndEdgeConnections;
+}
+
+
+void Foam::faMesh::setBoundaryConnections
+(
+ const List<Pair<patchTuple>>& bndEdgeConnections
+) const
+{
+ const label nInternalEdges = patch().nInternalEdges();
+ const label nBoundaryEdges = patch().nBoundaryEdges();
+
+ if (bndEdgeConnections.size() != nBoundaryEdges)
+ {
+ FatalErrorInFunction
+ << "Sizing mismatch. Expected " << nBoundaryEdges
+ << " boundary edge connections, but had "
+ << bndEdgeConnections.size() << nl
+ << abort(FatalError);
+ }
+
+ bndConnectPtr_.reset
+ (
+ new List<labelPair>(nBoundaryEdges, labelPair(-1,-1))
+ );
+ auto& bndConnect = *bndConnectPtr_;
+
+ for (const auto& connection : bndEdgeConnections)
+ {
+ const auto& a = connection.first();
+ const auto& b = connection.second();
+
+ if (a.is_finiteArea() && a.is_localProc())
+ {
+ const label bndEdgei = (a.patchEdgei() - nInternalEdges);
+
+ bndConnect[bndEdgei].first() = b.procNo();
+ bndConnect[bndEdgei].second() = b.meshFacei();
+ }
+ else if (b.is_finiteArea() && b.is_localProc())
+ {
+ const label bndEdgei = (b.patchEdgei() - nInternalEdges);
+
+ bndConnect[bndEdgei].first() = a.procNo();
+ bndConnect[bndEdgei].second() = a.meshFacei();
+ }
+ else
+ {
+ FatalErrorInFunction
+ << "Unexpected pairing input " << connection
+ << " ... programming error" << nl
+ << abort(FatalError);
+ }
+ }
+
+ label nInvalid = 0;
+ for (const auto& connection : bndConnect)
+ {
+ if (connection.first() < 0 || connection.second() < 0)
+ {
+ ++nInvalid;
+ }
+ }
+
+ if (Pstream::parRun())
+ {
+ reduce(nInvalid, sumOp<label>());
+ }
+
+ if (nInvalid)
+ {
+ FatalErrorInFunction
+ << "Did not properly match " << nInvalid
+ << " boundary edges" << nl
+ << abort(FatalError);
+ }
+}
+
+
+void Foam::faMesh::calcBoundaryConnections() const
+{
+ setBoundaryConnections(this->getBoundaryEdgeConnections());
+}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+Foam::labelList Foam::faMesh::boundaryProcs() const
+{
+ const auto& connections = this->boundaryConnections();
+
+ labelHashSet procsUsed(2*Pstream::nProcs());
+
+ for (const labelPair& tuple : connections)
+ {
+ procsUsed.insert(tuple.first());
+ }
+
+ procsUsed.erase(-1); // placeholder value
+ procsUsed.erase(Pstream::myProcNo());
+
+ return procsUsed.sortedToc();
+}
+
+
+Foam::List<Foam::labelPair> Foam::faMesh::boundaryProcSizes() const
+{
+ const auto& connections = this->boundaryConnections();
+
+ Map<label> procCount(2*Pstream::nProcs());
+
+ for (const labelPair& tuple : connections)
+ {
+ ++procCount(tuple.first());
+ }
+ procCount.erase(-1); // placeholder value
+ procCount.erase(Pstream::myProcNo());
+
+ // Flatten as list
+ List<labelPair> output(procCount.size());
+ label count = 0;
+ for (const label proci : procCount.sortedToc())
+ {
+ output[count].first() = proci;
+ output[count].second() = procCount[proci]; // size
+ ++count;
+ }
+
+ return output;
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteArea/faMesh/faPatches/faPatch/faPatch.C b/src/finiteArea/faMesh/faPatches/faPatch/faPatch.C
index 6ca7f3db147..8be866c7c03 100644
--- a/src/finiteArea/faMesh/faPatches/faPatch/faPatch.C
+++ b/src/finiteArea/faMesh/faPatches/faPatch/faPatch.C
@@ -116,12 +116,6 @@ Foam::faPatch::~faPatch()
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
-Foam::label Foam::faPatch::ngbPolyPatchIndex() const noexcept
-{
- return nbrPolyPatchId_;
-}
-
-
const Foam::faBoundaryMesh& Foam::faPatch::boundaryMesh() const noexcept
{
return boundaryMesh_;
@@ -134,6 +128,76 @@ Foam::label Foam::faPatch::start() const
}
+Foam::List<Foam::labelPair> Foam::faPatch::boundaryConnections() const
+{
+ const auto& connections = boundaryMesh().mesh().boundaryConnections();
+ const label nInternalEdges = boundaryMesh().mesh().nInternalEdges();
+
+ List<labelPair> output(this->nEdges());
+
+ // Like an IndirectList but removing the nInternalEdges offset
+ label count = 0;
+ for (const label patchEdgei : this->edgeLabels())
+ {
+ const label bndEdgei = (patchEdgei - nInternalEdges);
+ output[count] = connections[bndEdgei];
+ ++count;
+ }
+
+ return output;
+}
+
+
+Foam::labelList Foam::faPatch::boundaryProcs() const
+{
+ const auto& connections = boundaryMesh().mesh().boundaryConnections();
+ const label nInternalEdges = boundaryMesh().mesh().nInternalEdges();
+
+ labelHashSet procsUsed(2*Pstream::nProcs());
+
+ for (const label patchEdgei : this->edgeLabels())
+ {
+ const label bndEdgei = (patchEdgei - nInternalEdges);
+ const label proci = connections[bndEdgei].first();
+ procsUsed.insert(proci);
+ }
+ procsUsed.erase(-1); // placeholder value
+ procsUsed.erase(Pstream::myProcNo());
+
+ return procsUsed.sortedToc();
+}
+
+
+Foam::List<Foam::labelPair> Foam::faPatch::boundaryProcSizes() const
+{
+ const auto& connections = boundaryMesh().mesh().boundaryConnections();
+ const label nInternalEdges = boundaryMesh().mesh().nInternalEdges();
+
+ Map<label> procCount(2*Pstream::nProcs());
+
+ for (const label patchEdgei : this->edgeLabels())
+ {
+ const label bndEdgei = (patchEdgei - nInternalEdges);
+ const label proci = connections[bndEdgei].first();
+ ++procCount(proci);
+ }
+ procCount.erase(-1); // placeholder value
+ procCount.erase(Pstream::myProcNo());
+
+ // Flatten as list
+ List<labelPair> output(procCount.size());
+ label count = 0;
+ for (const label proci : procCount.sortedToc())
+ {
+ output[count].first() = proci;
+ output[count].second() = procCount[proci]; // size
+ ++count;
+ }
+
+ return output;
+}
+
+
const Foam::labelList& Foam::faPatch::pointLabels() const
{
if (!pointLabelsPtr_)
diff --git a/src/finiteArea/faMesh/faPatches/faPatch/faPatch.H b/src/finiteArea/faMesh/faPatches/faPatch/faPatch.H
index 9be7e6153ef..731647bd624 100644
--- a/src/finiteArea/faMesh/faPatches/faPatch/faPatch.H
+++ b/src/finiteArea/faMesh/faPatches/faPatch/faPatch.H
@@ -231,6 +231,7 @@ public:
return static_cast<const labelList&>(*this);
}
+ //- Define new list of edges
void edgeLabels(const UList<label>& newEdgeLabels);
//- Number of patch points
@@ -245,8 +246,11 @@ public:
return labelList::size();
}
- //- Return neighbour polyPatch index
- label ngbPolyPatchIndex() const noexcept;
+ //- The neighbour polyPatch index
+ label ngbPolyPatchIndex() const noexcept
+ {
+ return nbrPolyPatchId_;
+ }
//- Return boundaryMesh reference
const faBoundaryMesh& boundaryMesh() const noexcept;
@@ -284,6 +288,21 @@ public:
virtual void write(Ostream&) const;
+ // Topological information
+
+ //- List of proc/face for the boundary edge neighbours
+ //- in locally reordered edge numbering.
+ List<labelPair> boundaryConnections() const;
+
+ //- Boundary edge neighbour processors
+ //- (does not include own proc)
+ labelList boundaryProcs() const;
+
+ //- List of proc/size for the boundary edge neighbour processors
+ //- (does not include own proc)
+ List<labelPair> boundaryProcSizes() const;
+
+
// Access functions for geometrical data
//- Return patch point labels
diff --git a/src/finiteArea/faMesh/faPatches/faPatch/faPatchData.C b/src/finiteArea/faMesh/faPatches/faPatch/faPatchData.C
index c83150cea4c..00ad7d679cb 100644
--- a/src/finiteArea/faMesh/faPatches/faPatch/faPatchData.C
+++ b/src/finiteArea/faMesh/faPatches/faPatch/faPatchData.C
@@ -28,6 +28,7 @@ License
#include "faPatchData.H"
#include "dictionary.H"
#include "processorFaPatch.H"
+#include "processorPolyPatch.H"
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
@@ -108,6 +109,24 @@ void Foam::faPatchData::assign(const faPatch& fap)
}
+bool Foam::faPatchData::assign_coupled(int ownProci, int neiProci)
+{
+ clear();
+
+ if (ownProci == neiProci)
+ {
+ return false;
+ }
+
+ name_ = processorPolyPatch::newName(ownProci, neiProci);
+ type_ = processorFaPatch::typeName;
+ ownerProcId_ = ownProci;
+ neighProcId_ = neiProci;
+
+ return true;
+}
+
+
int Foam::faPatchData::matchPatchPair
(
const labelPair& patchPair
diff --git a/src/finiteArea/faMesh/faPatches/faPatch/faPatchData.H b/src/finiteArea/faMesh/faPatches/faPatch/faPatchData.H
index e66c91723e6..ee90a80522a 100644
--- a/src/finiteArea/faMesh/faPatches/faPatch/faPatchData.H
+++ b/src/finiteArea/faMesh/faPatches/faPatch/faPatchData.H
@@ -98,6 +98,9 @@ public:
//- Clear and populate with values from finiteArea patch
void assign(const faPatch& fap);
+ //- Set values consistent with a processor coupling
+ bool assign_coupled(int ownProci, int neiProci);
+
//- True if owner/neighbour processor ids are non-equal
bool coupled() const noexcept
{
--
GitLab