fvMeshStitcher: Fixes for cases with multiple NCC-s

A number of bugs have been fixed relating to cases in which multiple NCC
patches are locally edge-connected, and/or are edge-connected across a
processor boundary.

src/finiteVolume/fvMesh/fvMeshStitchers/fvMeshStitcher/fvMeshStitcher.C

@@ -1122,7 +1122,6 @@ Foam::fvMeshStitcher::calculateOwnerOrigBoundaryEdgeParts
 }
 
 
-
 void Foam::fvMeshStitcher::applyOwnerOrigBoundaryEdgeParts
 (
     surfaceVectorField& SfSf,
@@ -1159,7 +1158,8 @@ void Foam::fvMeshStitcher::applyOwnerOrigBoundaryEdgeParts
 
             const label patchi =
                 mesh_.isInternalFace(facei)
-              ? -1 : pbMesh.patchIndices()[facei - mesh_.nInternalFaces()];
+              ? -1
+              : pbMesh.patchIndices()[facei - mesh_.nInternalFaces()];
 
             if (patchi != -1 && patchIsOwnerOrig[patchi])
             {
@@ -1428,6 +1428,8 @@ void Foam::fvMeshStitcher::intersect
 
     const nonConformalBoundary& ncb = nonConformalBoundary::New(mesh_);
     const labelList ownerOrigPatchIndices = ncb.ownerOrigPatchIndices();
+    const edgeList& ownerOrigBoundaryMeshEdges =
+        ncb.ownerOrigBoundaryMeshEdges();
 
     // Alias the boundary geometry fields
     surfaceVectorField::Boundary& SfBf = SfSf.boundaryFieldRef();
@@ -1562,21 +1564,37 @@ void Foam::fvMeshStitcher::intersect
             const label ownerOrigBoundaryEdgei =
                 origPatchEdgeOwnerOrigBoundaryEdges[origPatchEdgei];
 
+            const label sign =
+                edge::compare
+                (
+                    meshEdge(origPatch, origPatchEdgei),
+                    ownerOrigBoundaryMeshEdges[ownerOrigBoundaryEdgei]
+                );
+
             part errorP =
                 patchEdgeParts[origPatchi][origPatchEdgei];
-            errorP -= ownerOrigBoundaryEdgeParts[ownerOrigBoundaryEdgei];
+            errorP -=
+                sign > 0
+              ? ownerOrigBoundaryEdgeParts[ownerOrigBoundaryEdgei]
+              : -ownerOrigBoundaryEdgeParts[ownerOrigBoundaryEdgei];
 
             forAll(origPatch.edgeFaces()[origPatchEdgei], patchEdgeFacei)
             {
                 const label patchFacei =
                     origPatch.edgeFaces()[origPatchEdgei][patchEdgeFacei];
 
+                const label sign =
+                    origPatch.localFaces()[patchFacei].edgeDirection
+                    (
+                        origPatch.edges()[origPatchEdgei]
+                    );
+
                 part p
                 (
                     SfBf[origPatchi][patchFacei],
                     CfBf[origPatchi][patchFacei]
                 );
-                p += errorP;
+                p += sign > 0 ? errorP : -errorP;
 
                 SfBf[origPatchi][patchFacei] = p.area;
                 CfBf[origPatchi][patchFacei] = p.centre;

src/finiteVolume/fvMesh/fvPatches/constraint/processor/processorFvPatch.C

@@ -25,7 +25,8 @@ License
 
 #include "processorFvPatch.H"
 #include "addToRunTimeSelectionTable.H"
-#include "transformField.H"
+#include "volFields.H"
+#include "surfaceFields.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -42,13 +43,26 @@ void Foam::processorFvPatch::makeWeights(scalarField& w) const
 {
     if (Pstream::parRun())
     {
-        coupledFvPatch::makeWeights
-        (
-            w,
-            procPolyPatch_.neighbFaceAreas(),
-            procPolyPatch_.neighbFaceCentres()
-          - procPolyPatch_.neighbFaceCellCentres()
-        );
+        if (!boundaryMesh().mesh().conformal())
+        {
+            coupledFvPatch::makeWeights
+            (
+                w,
+              - boundaryMesh().mesh().Sf().boundaryField()[index()],
+                boundaryMesh().mesh().Cf().boundaryField()[index()]
+              - boundaryMesh().mesh().C().boundaryField()[index()]
+            );
+        }
+        else
+        {
+            coupledFvPatch::makeWeights
+            (
+                w,
+                procPolyPatch_.neighbFaceAreas(),
+                procPolyPatch_.neighbFaceCentres()
+              - procPolyPatch_.neighbFaceCellCentres()
+            );
+        }
     }
     else
     {
@@ -61,12 +75,24 @@ Foam::tmp<Foam::vectorField> Foam::processorFvPatch::delta() const
 {
     if (Pstream::parRun())
     {
-        return
-            coupledFvPatch::delta
-            (
-                procPolyPatch_.neighbFaceCentres()
-              - procPolyPatch_.neighbFaceCellCentres()
-            );
+        if (!boundaryMesh().mesh().conformal())
+        {
+            return
+                coupledFvPatch::delta
+                (
+                    boundaryMesh().mesh().Cf().boundaryField()[index()]
+                  - boundaryMesh().mesh().C().boundaryField()[index()]
+                );
+        }
+        else
+        {
+            return
+                coupledFvPatch::delta
+                (
+                    procPolyPatch_.neighbFaceCentres()
+                  - procPolyPatch_.neighbFaceCellCentres()
+                );
+        }
     }
     else
     {

src/meshTools/nonConformal/boundary/nonConformalBoundary.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2022-2023 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2022-2024 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -179,28 +179,6 @@ void Foam::nonConformalBoundary::nonConformalOtherPatchIndices
 {}
 
 
-const Foam::labelList&
-Foam::nonConformalBoundary::meshPointOwnerOrigBoundaryPoint() const
-{
-    if (!meshPointOwnerOrigBoundaryPointPtr_.valid())
-    {
-        meshPointOwnerOrigBoundaryPointPtr_.set
-        (
-            new labelList(mesh().nPoints(), -1)
-        );
-
-        forAll(ownerOrigBoundary_.meshPoints(), ownerOrigBoundaryPointi)
-        {
-            meshPointOwnerOrigBoundaryPointPtr_()
-                [ownerOrigBoundary_.meshPoints()[ownerOrigBoundaryPointi]] =
-                ownerOrigBoundaryPointi;
-        }
-    }
-
-    return meshPointOwnerOrigBoundaryPointPtr_();
-}
-
-
 const Foam::vectorField&
 Foam::nonConformalBoundary::ownerOrigBoundaryPointNormals() const
 {
@@ -373,20 +351,30 @@ Foam::nonConformalBoundary::ownerOrigBoundaryPointMeshPoint() const
 {
     if (!ownerOrigBoundaryPointMeshPointPtr_.valid())
     {
+        // Construct the local maps using the owner-orig primitive patch
         ownerOrigBoundaryPointMeshPointPtr_.set
         (
             new labelList(ownerOrigBoundary_.meshPoints())
         );
 
-        // ...
-        meshPointOwnerOrigBoundaryPoint();
-        labelList& map = meshPointOwnerOrigBoundaryPointPtr_();
+        meshPointOwnerOrigBoundaryPointPtr_.set
+        (
+            new labelList(mesh().nPoints(), -1)
+        );
 
-        // ...
-        label ownerOrigBoundaryPointi = ownerOrigBoundary_.nPoints();
-        DynamicList<label> remotePoints;
+        labelList& meshPointOwnerOrigBoundaryPoint =
+            meshPointOwnerOrigBoundaryPointPtr_();
+
+        forAll(ownerOrigBoundary_.meshPoints(), ownerOrigBoundaryPointi)
+        {
+            meshPointOwnerOrigBoundaryPoint
+                [ownerOrigBoundary_.meshPoints()[ownerOrigBoundaryPointi]] =
+                ownerOrigBoundaryPointi;
+        }
 
-        // ...
+        // Construct the remote map by enumerating newly identified points
+        label ownerOrigBoundaryPointi = ownerOrigBoundary_.nPoints();
+        DynamicList<label> remoteMeshPoints;
         for
         (
             label ownerOrigBoundaryEdgei = ownerOrigBoundary_.nEdges();
@@ -403,16 +391,17 @@ Foam::nonConformalBoundary::ownerOrigBoundaryPointMeshPoint() const
             {
                 const label meshPointi = e[i];
 
-                if (map[meshPointi] == -1)
+                if (meshPointOwnerOrigBoundaryPoint[meshPointi] == -1)
                 {
-                    map[meshPointi] = ownerOrigBoundaryPointi ++;
-                    remotePoints.append(meshPointi);
+                    meshPointOwnerOrigBoundaryPoint[meshPointi] =
+                        ownerOrigBoundaryPointi ++;
+                    remoteMeshPoints.append(meshPointi);
                 }
             }
         }
 
-        // ...
-        ownerOrigBoundaryPointMeshPointPtr_->append(remotePoints);
+        // Append to the point-mesh-point map
+        ownerOrigBoundaryPointMeshPointPtr_->append(remoteMeshPoints);
     }
 
     return ownerOrigBoundaryPointMeshPointPtr_();
@@ -528,10 +517,11 @@ Foam::nonConformalBoundary::ownerOrigBoundaryEdges() const
             new edgeList(ownerOrigBoundary_.edges())
         );
 
-        const labelList& map = meshPointOwnerOrigBoundaryPoint();
+        ownerOrigBoundaryPointMeshPoint();
+        const labelList& meshPointOwnerOrigBoundaryPoint =
+            meshPointOwnerOrigBoundaryPointPtr_();
 
         DynamicList<edge> remoteEdges;
-
         for
         (
             label ownerOrigBoundaryEdgei = ownerOrigBoundary_.nEdges();
@@ -544,7 +534,14 @@ Foam::nonConformalBoundary::ownerOrigBoundaryEdges() const
 
             const edge& e = mesh().edges()[meshEdgei];
 
-            remoteEdges.append(edge(map[e.start()], map[e.end()]));
+            remoteEdges.append
+            (
+                edge
+                (
+                    meshPointOwnerOrigBoundaryPoint[e.start()],
+                    meshPointOwnerOrigBoundaryPoint[e.end()]
+                )
+            );
         }
 
         ownerOrigBoundaryEdgesPtr_->append(remoteEdges);
@@ -594,11 +591,15 @@ Foam::nonConformalBoundary::patchPointOwnerOrigBoundaryPoints
     {
         const polyPatch& pp = mesh().boundaryMesh()[patchi];
 
+        ownerOrigBoundaryPointMeshPoint();
+        const labelList& meshPointOwnerOrigBoundaryPoint =
+            meshPointOwnerOrigBoundaryPointPtr_();
+
         const faceList patchOwnerOrigBoundaryLocalFaces
         (
             renumber
             (
-                meshPointOwnerOrigBoundaryPoint(),
+                meshPointOwnerOrigBoundaryPoint,
                 static_cast<const faceList&>(pp)
             )
         );

src/meshTools/nonConformal/boundary/nonConformalBoundary.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2022-2023 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2022-2024 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -65,7 +65,7 @@ class nonConformalBoundary
         //- Primitive patch of the owner-orig boundary
         indirectPrimitivePatch ownerOrigBoundary_;
 
-        //- A map from owner-orig boundary point to mesh point
+        //- A map from mesh point to owner-orig boundary point
         mutable autoPtr<labelList> meshPointOwnerOrigBoundaryPointPtr_;
 
         //- A map from owner-orig boundary point to mesh point
@@ -135,10 +135,6 @@ class nonConformalBoundary
             const label side
         ) const;
 
-        //- Get a map from mesh point to owner-orig boundary points. Non
-        //  owner-orig boundary points are indicated by the value -1.
-        const labelList& meshPointOwnerOrigBoundaryPoint() const;
-
         //- Get point normals for the owner-orig boundary
         const vectorField& ownerOrigBoundaryPointNormals() const;