From 542dae4a6dbad75e296459ca38984dbc1184c3a0 Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Wed, 6 Jan 2021 09:54:37 +0000
Subject: [PATCH] ENH: syncTools: add edge orientation. Fixes #1974.

---
 .../meshes/polyMesh/syncTools/syncTools.H     |  24 ++--
 .../polyMesh/syncTools/syncToolsTemplates.C   | 125 ++++++++++++++++--
 .../snappySnapDriverFeature.C                 |   6 +-
 3 files changed, 130 insertions(+), 25 deletions(-)

diff --git a/src/OpenFOAM/meshes/polyMesh/syncTools/syncTools.H b/src/OpenFOAM/meshes/polyMesh/syncTools/syncTools.H
index 25412ca8ae..a3e8873e82 100644
--- a/src/OpenFOAM/meshes/polyMesh/syncTools/syncTools.H
+++ b/src/OpenFOAM/meshes/polyMesh/syncTools/syncTools.H
@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2016 OpenFOAM Foundation
-    Copyright (C) 2015-2019 OpenCFD Ltd.
+    Copyright (C) 2015-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -138,18 +138,19 @@ public:
             );
 
             //- Synchronize values on all mesh edges.
-            template<class T, class CombineOp, class TransformOp>
+            template<class T, class CombineOp, class TransformOp, class FlipOp>
             static void syncEdgeList
             (
                 const polyMesh& mesh,
                 List<T>& edgeValues,
                 const CombineOp& cop,
                 const T& nullValue,
-                const TransformOp& top
+                const TransformOp& top,
+                const FlipOp& fop
             );
 
             //- Synchronize values on selected mesh edges.
-            template<class T, class CombineOp, class TransformOp>
+            template<class T, class CombineOp, class TransformOp, class FlipOp>
             static void syncEdgeList
             (
                 const polyMesh& mesh,
@@ -157,7 +158,8 @@ public:
                 List<T>& edgeValues,
                 const CombineOp& cop,
                 const T& nullValue,
-                const TransformOp& top
+                const TransformOp& top,
+                const FlipOp& fop
             );
 
             //- Synchronize values on boundary faces only.
@@ -277,7 +279,8 @@ public:
                     edgeValues,
                     cop,
                     nullValue,
-                    mapDistribute::transform()
+                    mapDistribute::transform(),
+                    noOp()
                 );
             }
 
@@ -297,7 +300,8 @@ public:
                     positions,
                     cop,
                     nullValue,
-                    mapDistribute::transformPosition()
+                    mapDistribute::transformPosition(),
+                    noOp()
                 );
             }
 
@@ -319,7 +323,8 @@ public:
                     edgeValues,
                     cop,
                     nullValue,
-                    mapDistribute::transform()
+                    mapDistribute::transform(),
+                    noOp()
                 );
             }
 
@@ -341,7 +346,8 @@ public:
                     positions,
                     cop,
                     nullValue,
-                    mapDistribute::transformPosition()
+                    mapDistribute::transformPosition(),
+                    noOp()
                 );
             }
 
diff --git a/src/OpenFOAM/meshes/polyMesh/syncTools/syncToolsTemplates.C b/src/OpenFOAM/meshes/polyMesh/syncTools/syncToolsTemplates.C
index 133535749f..54be93cee7 100644
--- a/src/OpenFOAM/meshes/polyMesh/syncTools/syncToolsTemplates.C
+++ b/src/OpenFOAM/meshes/polyMesh/syncTools/syncToolsTemplates.C
@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2017 OpenFOAM Foundation
-    Copyright (C) 2015-2020 OpenCFD Ltd.
+    Copyright (C) 2015-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -798,14 +798,15 @@ void Foam::syncTools::syncPointList
 }
 
 
-template<class T, class CombineOp, class TransformOp>
+template<class T, class CombineOp, class TransformOp, class FlipOp>
 void Foam::syncTools::syncEdgeList
 (
     const polyMesh& mesh,
     List<T>& edgeValues,
     const CombineOp& cop,
     const T& nullValue,
-    const TransformOp& top
+    const TransformOp& top,
+    const FlipOp& fop
 )
 {
     if (edgeValues.size() != mesh.nEdges())
@@ -816,12 +817,47 @@ void Foam::syncTools::syncEdgeList
             << mesh.nEdges() << abort(FatalError);
     }
 
+    const edgeList& edges = mesh.edges();
     const globalMeshData& gd = mesh.globalData();
     const labelList& meshEdges = gd.coupledPatchMeshEdges();
+    const indirectPrimitivePatch& cpp = gd.coupledPatch();
+    const edgeList& cppEdges = cpp.edges();
+    const labelList& mp = cpp.meshPoints();
     const globalIndexAndTransform& git = gd.globalTransforms();
     const mapDistribute& edgeMap = gd.globalEdgeSlavesMap();
+    const bitSet& orientation = gd.globalEdgeOrientation();
+
+    List<T> cppFld(meshEdges.size());
+    forAll(meshEdges, i)
+    {
+        const edge& cppE = cppEdges[i];
+        const label meshEdgei = meshEdges[i];
+        const edge& meshE = edges[meshEdgei];
+
+        // 1. is cpp edge oriented as mesh edge
+        // 2. is cpp edge oriented same as master edge
+
+        const int dir = edge::compare(meshE, edge(mp, cppE));
+        if (dir == 0)
+        {
+            FatalErrorInFunction<< "Problem:"
+                << " mesh edge:" << meshE.line(mesh.points())
+                << " coupled edge:" << cppE.line(cpp.localPoints())
+                << exit(FatalError);
+        }
+
+        const bool sameOrientation = ((dir == 1) == orientation[i]);
+
+        if (sameOrientation)
+        {
+            cppFld[i] = edgeValues[meshEdgei];
+        }
+        else
+        {
+            cppFld[i] = fop(edgeValues[meshEdgei]);
+        }
+    }
 
-    List<T> cppFld(UIndirectList<T>(edgeValues, meshEdges));
 
     globalMeshData::syncData
     (
@@ -837,12 +873,29 @@ void Foam::syncTools::syncEdgeList
     // Extract back onto mesh
     forAll(meshEdges, i)
     {
-        edgeValues[meshEdges[i]] = cppFld[i];
+        const edge& cppE = cppEdges[i];
+        const label meshEdgei = meshEdges[i];
+        const edge& meshE = edges[meshEdgei];
+
+        // 1. is cpp edge oriented as mesh edge
+        // 2. is cpp edge oriented same as master edge
+
+        const int dir = edge::compare(meshE, edge(mp, cppE));
+        const bool sameOrientation = ((dir == 1) == orientation[i]);
+
+        if (sameOrientation)
+        {
+            edgeValues[meshEdges[i]] = cppFld[i];
+        }
+        else
+        {
+            edgeValues[meshEdges[i]] =  fop(cppFld[i]);
+        }
     }
 }
 
 
-template<class T, class CombineOp, class TransformOp>
+template<class T, class CombineOp, class TransformOp, class FlipOp>
 void Foam::syncTools::syncEdgeList
 (
     const polyMesh& mesh,
@@ -850,7 +903,8 @@ void Foam::syncTools::syncEdgeList
     List<T>& edgeValues,
     const CombineOp& cop,
     const T& nullValue,
-    const TransformOp& top
+    const TransformOp& top,
+    const FlipOp& fop
 )
 {
     if (edgeValues.size() != meshEdges.size())
@@ -860,19 +914,48 @@ void Foam::syncTools::syncEdgeList
             << " is not equal to the number of meshEdges "
             << meshEdges.size() << abort(FatalError);
     }
+    const edgeList& edges = mesh.edges();
     const globalMeshData& gd = mesh.globalData();
     const indirectPrimitivePatch& cpp = gd.coupledPatch();
+    const edgeList& cppEdges = cpp.edges();
+    const labelList& mp = cpp.meshPoints();
     const Map<label>& mpm = gd.coupledPatchMeshEdgeMap();
+    const bitSet& orientation = gd.globalEdgeOrientation();
 
     List<T> cppFld(cpp.nEdges(), nullValue);
 
     forAll(meshEdges, i)
     {
-        const auto iter = mpm.cfind(meshEdges[i]);
-
+        const label meshEdgei = meshEdges[i];
+        const auto iter = mpm.cfind(meshEdgei);
         if (iter.found())
         {
-            cppFld[*iter] = edgeValues[i];
+            const label cppEdgei = iter();
+            const edge& cppE = cppEdges[cppEdgei];
+            const edge& meshE = edges[meshEdgei];
+
+            // 1. is cpp edge oriented as mesh edge
+            // 2. is cpp edge oriented same as master edge
+
+            const int dir = edge::compare(meshE, edge(mp, cppE));
+            if (dir == 0)
+            {
+                FatalErrorInFunction<< "Problem:"
+                    << " mesh edge:" << meshE.line(mesh.points())
+                    << " coupled edge:" << cppE.line(cpp.localPoints())
+                    << exit(FatalError);
+            }
+
+            const bool sameOrientation = ((dir == 1) == orientation[i]);
+
+            if (sameOrientation)
+            {
+                cppFld[cppEdgei] = edgeValues[i];
+            }
+            else
+            {
+                cppFld[cppEdgei] = fop(edgeValues[i]);
+            }
         }
     }
 
@@ -889,11 +972,25 @@ void Foam::syncTools::syncEdgeList
 
     forAll(meshEdges, i)
     {
-        const auto iter = mpm.cfind(meshEdges[i]);
-
-        if (iter.found())
+        label meshEdgei = meshEdges[i];
+        Map<label>::const_iterator iter = mpm.find(meshEdgei);
+        if (iter != mpm.end())
         {
-            edgeValues[i] = cppFld[*iter];
+            label cppEdgei = iter();
+            const edge& cppE = cppEdges[cppEdgei];
+            const edge& meshE = edges[meshEdgei];
+
+            const int dir = edge::compare(meshE, edge(mp, cppE));
+            const bool sameOrientation = ((dir == 1) == orientation[i]);
+
+            if (sameOrientation)
+            {
+                edgeValues[i] = cppFld[cppEdgei];
+            }
+            else
+            {
+                edgeValues[i] = fop(cppFld[cppEdgei]);
+            }
         }
     }
 }
diff --git a/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappySnapDriverFeature.C b/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappySnapDriverFeature.C
index cbe23e027f..434c1360e3 100644
--- a/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappySnapDriverFeature.C
+++ b/src/mesh/snappyHexMesh/snappyHexMeshDriver/snappySnapDriverFeature.C
@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2015 OpenFOAM Foundation
-    Copyright (C) 2015-2020 OpenCFD Ltd.
+    Copyright (C) 2015-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -2873,6 +2873,7 @@ void Foam::snappySnapDriver::determineBaffleFeatures
         (
             pp.meshEdges(mesh.edges(), mesh.pointEdges())
         );
+        // Collect all coupled edges. Does not filter duplicates/order
         syncTools::syncEdgeList
         (
             mesh,
@@ -2880,7 +2881,8 @@ void Foam::snappySnapDriver::determineBaffleFeatures
             edgeFaceNormals,
             listPlusEqOp<point>(),
             List<point>(),
-            mapDistribute::transform()
+            mapDistribute::transform(),
+            noOp()
         );
     }