From ea92cb82c4e8dfcb773802d37817187ef0d1bd1e Mon Sep 17 00:00:00 2001
From: Mark Olesen <Mark.Olesen@esi-group.com>
Date: Fri, 8 Oct 2021 15:45:19 +0200
Subject: [PATCH] ENH: add boundary halo handling to faMesh

- the finiteArea is typically restricteed to one or more patches on a
  polyMesh but will have an external connectivity to other parts of
  the polyMesh. For proper determination of the point normals (for
  example), the neighbouring information is needed.

  These outside 'halo' faces can be located on different processors,
  but not correspond to an internal processor-processor interface.

  Add a dedicated form of mapDistribute for swapping this type of
  information. Since this is a collective operation, locate
  corresponding methods directly on the faMesh level instead of
  the faPatch level.
---
 src/finiteArea/Make/files                     |   1 +
 src/finiteArea/faMesh/faMesh.C                |  57 ++++-
 src/finiteArea/faMesh/faMesh.H                |  62 +++++-
 src/finiteArea/faMesh/faMeshBoundaryHalo.C    | 194 ++++++++++++++++++
 src/finiteArea/faMesh/faMeshBoundaryHalo.H    | 154 ++++++++++++++
 .../faMesh/faMeshBoundaryHaloTemplates.C      | 107 ++++++++++
 src/finiteArea/faMesh/faMeshTopology.C        | 127 ++++++++++++
 7 files changed, 699 insertions(+), 3 deletions(-)
 create mode 100644 src/finiteArea/faMesh/faMeshBoundaryHalo.C
 create mode 100644 src/finiteArea/faMesh/faMeshBoundaryHalo.H
 create mode 100644 src/finiteArea/faMesh/faMeshBoundaryHaloTemplates.C

diff --git a/src/finiteArea/Make/files b/src/finiteArea/Make/files
index 85d56a5e87..2badfb337d 100644
--- a/src/finiteArea/Make/files
+++ b/src/finiteArea/Make/files
@@ -4,6 +4,7 @@ faMesh/faMeshDemandDrivenData.C
 faMesh/faMeshPatches.C
 faMesh/faMeshTopology.C
 faMesh/faMeshUpdate.C
+faMesh/faMeshBoundaryHalo.C
 faMesh/faBoundaryMesh/faBoundaryMesh.C
 
 faPatches = faMesh/faPatches
diff --git a/src/finiteArea/faMesh/faMesh.C b/src/finiteArea/faMesh/faMesh.C
index 021b83a327..40763bce48 100644
--- a/src/finiteArea/faMesh/faMesh.C
+++ b/src/finiteArea/faMesh/faMesh.C
@@ -27,6 +27,7 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "faMesh.H"
+#include "faMeshBoundaryHalo.H"
 #include "faGlobalMeshData.H"
 #include "Time.H"
 #include "polyMesh.H"
@@ -109,6 +110,36 @@ static labelList selectPatchFaces
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
+void Foam::faMesh::checkBoundaryEdgeLabelRange
+(
+    const labelUList& edgeLabels
+) const
+{
+    label nErrors = 0;
+
+    for (const label edgei : edgeLabels)
+    {
+        if (edgei < nInternalEdges_ || edgei >= nEdges_)
+        {
+            if (!nErrors++)
+            {
+                FatalErrorInFunction
+                    << "Boundary edge label out of range "
+                    << nInternalEdges_ << ".." << (nEdges_-1) << nl
+                    << "   ";
+            }
+
+            FatalError<< ' ' << edgei;
+        }
+    }
+
+    if (nErrors)
+    {
+        FatalError << nl << exit(FatalError);
+    }
+}
+
+
 void Foam::faMesh::initPatch() const
 {
     patchPtr_.reset
@@ -120,6 +151,9 @@ void Foam::faMesh::initPatch() const
         )
     );
     bndConnectPtr_.reset(nullptr);
+    haloMapPtr_.reset(nullptr);
+    haloFaceCentresPtr_.reset(nullptr);
+    haloFaceNormalsPtr_.reset(nullptr);
 }
 
 
@@ -168,10 +202,21 @@ void Foam::faMesh::setPrimitiveMeshData()
 }
 
 
+void Foam::faMesh::clearHalo() const
+{
+    DebugInFunction << "Clearing halo information" << endl;
+
+    haloMapPtr_.reset(nullptr);
+    haloFaceCentresPtr_.reset(nullptr);
+    haloFaceNormalsPtr_.reset(nullptr);
+}
+
+
 void Foam::faMesh::clearGeomNotAreas() const
 {
     DebugInFunction << "Clearing geometry" << endl;
 
+    clearHalo();
     patchPtr_.reset(nullptr);
     bndConnectPtr_.reset(nullptr);
     deleteDemandDrivenData(SPtr_);
@@ -274,7 +319,11 @@ Foam::faMesh::faMesh(const polyMesh& pMesh)
     faceCurvaturesPtr_(nullptr),
     edgeTransformTensorsPtr_(nullptr),
     correctPatchPointNormalsPtr_(nullptr),
-    globalMeshDataPtr_(nullptr)
+    globalMeshDataPtr_(nullptr),
+
+    haloMapPtr_(nullptr),
+    haloFaceCentresPtr_(nullptr),
+    haloFaceNormalsPtr_(nullptr)
 {
     DebugInFunction << "Creating from IOobject" << endl;
 
@@ -368,7 +417,11 @@ Foam::faMesh::faMesh
     faceCurvaturesPtr_(nullptr),
     edgeTransformTensorsPtr_(nullptr),
     correctPatchPointNormalsPtr_(nullptr),
-    globalMeshDataPtr_(nullptr)
+    globalMeshDataPtr_(nullptr),
+
+    haloMapPtr_(nullptr),
+    haloFaceCentresPtr_(nullptr),
+    haloFaceNormalsPtr_(nullptr)
 {}
 
 
diff --git a/src/finiteArea/faMesh/faMesh.H b/src/finiteArea/faMesh/faMesh.H
index 313db0517a..c18674b716 100644
--- a/src/finiteArea/faMesh/faMesh.H
+++ b/src/finiteArea/faMesh/faMesh.H
@@ -71,6 +71,7 @@ namespace Foam
 {
 
 // Forward Declarations
+class faMeshBoundaryHalo;
 class faMeshLduAddressing;
 class faMeshMapper;
 class faPatchData;
@@ -291,11 +292,20 @@ class faMesh
         mutable boolList* correctPatchPointNormalsPtr_;
 
 
-        // Other mesh-related data
+    // Other mesh-related data
 
         //- Parallel info
         mutable autoPtr<faGlobalMeshData> globalMeshDataPtr_;
 
+        //- Mapping/swapping for boundary edge halo neighbours
+        mutable std::unique_ptr<faMeshBoundaryHalo> haloMapPtr_;
+
+        //- Face centres for boundary edge halo neighbours
+        mutable std::unique_ptr<pointField> haloFaceCentresPtr_;
+
+        //- Face normals for boundary edge halo neighbours
+        mutable std::unique_ptr<vectorField> haloFaceNormalsPtr_;
+
 
     // Static Private Data
 
@@ -376,6 +386,9 @@ class faMesh
             //- Clear geometry but not the face areas
             void clearGeomNotAreas() const;
 
+            //- Clear boundary halo information
+            void clearHalo() const;
+
             //- Clear geometry
             void clearGeom() const;
 
@@ -386,6 +399,12 @@ class faMesh
             void clearOut() const;
 
 
+    // Halo handling
+
+        //- Calculate halo centres/normals
+        void calcHaloFaceGeometry() const;
+
+
     // Helpers
 
         //- Create a single patch
@@ -404,6 +423,32 @@ class faMesh
         ) const;
 
 
+        //- Fatal error if edge labels are out of range
+        void checkBoundaryEdgeLabelRange(const labelUList& edgeLabels) const;
+
+        //- Extract list from contiguous (unordered) boundary data
+        //- to the locally sorted order.
+        template<class T>
+        List<T> boundarySubset
+        (
+            const UList<T>& bndField,
+            const labelUList& edgeLabels
+        ) const
+        {
+            #ifdef FULLDEBUG
+            checkBoundaryEdgeLabelRange(edgeLabels);
+            #endif
+            // Like UIndirectList but with an offset
+            List<T> result(edgeLabels.size());
+            forAll(edgeLabels, i)
+            {
+                result[i] = bndField[edgeLabels[i] - nInternalEdges_];
+            }
+            return result;
+        }
+
+
+
 public:
 
     // Public Typedefs
@@ -607,6 +652,21 @@ public:
             //- (does not include own proc)
             List<labelPair> boundaryProcSizes() const;
 
+            //- Mapping/swapping for boundary halo neighbours
+            const faMeshBoundaryHalo& boundaryHaloMap() const;
+
+            //- Face centres of boundary halo neighbours
+            const pointField& haloFaceCentres() const;
+
+            //- Face normals of boundary halo neighbours
+            const vectorField& haloFaceNormals() const;
+
+            //- Face centres of boundary halo neighbours for specified patch
+            tmp<pointField> haloFaceCentres(const label patchi) const;
+
+            //- Face normals of boundary halo neighbours for specified patch
+            tmp<vectorField> haloFaceNormals(const label patchi) const;
+
 
         // Mesh motion and morphing
 
diff --git a/src/finiteArea/faMesh/faMeshBoundaryHalo.C b/src/finiteArea/faMesh/faMeshBoundaryHalo.C
new file mode 100644
index 0000000000..81aa63c117
--- /dev/null
+++ b/src/finiteArea/faMesh/faMeshBoundaryHalo.C
@@ -0,0 +1,194 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  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 "faMeshBoundaryHalo.H"
+#include "faMesh.H"
+#include "globalIndex.H"
+#include "Pstream.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(faMeshBoundaryHalo, 0);
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::faMeshBoundaryHalo::faMeshBoundaryHalo(const label comm)
+:
+    mapDistributeBase(comm),
+    inputMeshFaces_(),
+    boundaryToCompact_()
+{}
+
+
+Foam::faMeshBoundaryHalo::faMeshBoundaryHalo(const faMesh& areaMesh)
+:
+    mapDistributeBase(),
+    inputMeshFaces_(),
+    boundaryToCompact_()
+{
+    reset(areaMesh);
+}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+void Foam::faMeshBoundaryHalo::clear()
+{
+    static_cast<mapDistributeBase&>(*this) = mapDistributeBase();
+
+    inputMeshFaces_.clear();
+    boundaryToCompact_.clear();
+}
+
+
+Foam::label Foam::faMeshBoundaryHalo::haloSize() const
+{
+    if (Pstream::parRun())
+    {
+        return boundaryToCompact_.size();
+    }
+    else
+    {
+        return inputMeshFaces_.size();
+    }
+}
+
+
+void Foam::faMeshBoundaryHalo::reset(const faMesh& areaMesh)
+{
+    inputMeshFaces_.clear();
+    boundaryToCompact_.clear();
+
+    const auto& procConnections = areaMesh.boundaryConnections();
+
+    if (!Pstream::parRun())
+    {
+        // Serial - extract halo numbers directly
+
+        inputMeshFaces_.resize(procConnections.size());
+
+        forAll(procConnections, connecti)
+        {
+            // Connected neighbour, non-parallel = must be local
+            const auto& tuple = procConnections[connecti];
+            // const label nbrProci = tuple.first();
+            const label nbrFacei = tuple.second();
+
+            inputMeshFaces_[connecti] = nbrFacei;
+        }
+
+        return;
+    }
+
+    const label nProcs = Pstream::nProcs(comm_);
+    const label myRank = Pstream::myProcNo(comm_);
+
+    const globalIndex globalFaceNum(areaMesh.mesh().nFaces());
+
+    // Boundary inside faces in polyMesh face ids
+    const labelList insideFaces
+    (
+        UIndirectList<label>
+        (
+            areaMesh.faceLabels(),
+            areaMesh.patch().boundaryFaces()
+        )
+    );
+
+
+    // Slightly circuitous, but allows maximum reuse of mapDistributeBase
+
+    // 1. Construct a connectivity map using global face numbers
+
+    labelListList connectivity(areaMesh.nBoundaryEdges());
+    List<Map<label>> compactMap(nProcs, Map<label>(0));
+
+    // All local mesh faces used
+    labelHashSet localUsed(insideFaces);
+
+    forAll(connectivity, connecti)
+    {
+        labelList& edgeFaces = connectivity[connecti];
+        edgeFaces.resize(2);
+
+        // Owner is the boundary inside face (our side)
+        // Neighbour is the boundary outside face
+
+        // Connected neighbour
+        const auto& tuple = procConnections[connecti];
+        const label nbrProci = tuple.first();
+        const label nbrFacei = tuple.second();
+
+        if (myRank == nbrProci)
+        {
+            // Processor-local connectivity
+            localUsed.insert(nbrFacei);
+        }
+
+        // Global addressing for the connectivity
+        edgeFaces[0] = globalFaceNum.toGlobal(insideFaces[connecti]);
+        edgeFaces[1] = globalFaceNum.toGlobal(nbrProci, nbrFacei);
+    }
+
+    // Create and replace mapping
+    static_cast<mapDistributeBase&>(*this) = mapDistributeBase
+    (
+        globalFaceNum,
+        connectivity,
+        compactMap,
+        Pstream::msgType(),
+        comm_
+    );
+
+    // List of local mesh faces referenced.
+    // Includes inside and locally connected outside faces
+
+    inputMeshFaces_ = localUsed.sortedToc();
+
+    boundaryToCompact_.clear();
+    boundaryToCompact_.resize(connectivity.size());
+
+    // After creating the map, connectivity is localized *and*
+    // uses compact numbering!
+
+    // Extract the neighbour connection (compact numbering)
+    forAll(connectivity, connecti)
+    {
+        const labelList& edgeFaces = connectivity[connecti];
+        // const label face0 = edgeFaces[0];
+        const label face1 = edgeFaces[1];
+
+        boundaryToCompact_[connecti] = face1;
+    }
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteArea/faMesh/faMeshBoundaryHalo.H b/src/finiteArea/faMesh/faMeshBoundaryHalo.H
new file mode 100644
index 0000000000..4a8f58b33e
--- /dev/null
+++ b/src/finiteArea/faMesh/faMeshBoundaryHalo.H
@@ -0,0 +1,154 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  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/>.
+
+Class
+    Foam::faMeshBoundaryHalo
+
+Description
+    Class for obtaining halo face data for the boundary edges.
+    The ordering follows that natural edge ordering of the underlying
+    primitive patch.
+
+Note
+    The halo faces can be located on-processor or off-processor.
+
+SourceFiles
+    faMeshBoundaryHalo.C
+    faMeshBoundaryHaloTemplates.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef faMeshBoundaryHalo_H
+#define faMeshBoundaryHalo_H
+
+#include "mapDistributeBase.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// Forward Declarations
+class faMesh;
+
+/*---------------------------------------------------------------------------*\
+                      Class faMeshBoundaryHalo Declaration
+\*---------------------------------------------------------------------------*/
+
+class faMeshBoundaryHalo
+:
+    public mapDistributeBase
+{
+    // Private Data
+
+        //- List of local input mesh faces required
+        labelList inputMeshFaces_;
+
+        //- Internal mapping from boundary index to compact
+        labelList boundaryToCompact_;
+
+
+public:
+
+    // Declare name of the class and its debug switch
+    ClassName("faMeshBoundaryHalo");
+
+
+    // Constructors
+
+        //- Default construct
+        explicit faMeshBoundaryHalo(const label comm = UPstream::worldComm);
+
+        //- Construct from mesh
+        explicit faMeshBoundaryHalo(const faMesh& mesh);
+
+
+    // Member Functions
+
+        //- Clear out all parameters
+        void clear();
+
+        //- Redefine map and connectivity for a mesh
+        void reset(const faMesh& mesh);
+
+        //- The local data size (output)
+        label haloSize() const;
+
+        //- List of local input mesh faces required.
+        //  \note will not correspond exactly to the boundary inside faces.
+        //  Duplicates have been removed and it also contains the
+        //  processor-local neighbour faces, which would otherwise not be
+        //  handled by the distribute method.
+        const labelList& inputMeshFaces() const noexcept
+        {
+            return inputMeshFaces_;
+        }
+
+
+    // Other
+
+        //- Distribute sparse data.
+        //  On output it is adjusted.
+        template<class Type>
+        void distributeSparse
+        (
+            List<Type>& fld,
+            const labelUList& sparseInputLocations,
+            const labelUList& compactOutputMapping
+        ) const;
+
+        //- Distribute sparse data.
+        //  On output it is adjusted.
+        template<class Type>
+        void distributeSparse
+        (
+            List<Type>& fld,
+            const labelUList& sparseInputLocations
+        ) const;
+
+        //- Distribute sparse data.
+        //  The input field one enty per sparse id (inputMeshFaces).
+        //  On output it will have for the input sparse
+        //  The input field contains  location.
+        template<class Type>
+        void distributeSparse(List<Type>& fld) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+    #include "faMeshBoundaryHaloTemplates.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/finiteArea/faMesh/faMeshBoundaryHaloTemplates.C b/src/finiteArea/faMesh/faMeshBoundaryHaloTemplates.C
new file mode 100644
index 0000000000..69f2214b04
--- /dev/null
+++ b/src/finiteArea/faMesh/faMeshBoundaryHaloTemplates.C
@@ -0,0 +1,107 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  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 "faMeshBoundaryHalo.H"
+#include "UIndirectList.H"
+#include "Pstream.H"
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+template<class Type>
+void Foam::faMeshBoundaryHalo::distributeSparse
+(
+    List<Type>& fld,
+    const labelUList& sparseInputLocations,
+    const labelUList& compactOutputMapping
+) const
+{
+    if (!Pstream::parRun())
+    {
+        return;  // No-op in serial
+    }
+
+    // Construct data in compact addressing
+    List<Type> compactFld(constructSize_, Zero);
+
+    if (sparseInputLocations.empty())
+    {
+        // Copy in as dense field
+        forAll(fld, i)
+        {
+            compactFld[i] = fld[i];
+        }
+    }
+    else
+    {
+        if (fld.size() != sparseInputLocations.size())
+        {
+            FatalErrorInFunction
+                << "Input field size (" << fld.size()
+                << " != sparse ids size ("
+                << sparseInputLocations.size() << ")\n"
+                << exit(FatalError);
+        }
+
+        // Copy in sparse locations
+        forAll(sparseInputLocations, i)
+        {
+            const label idx = sparseInputLocations[i];
+            if (idx != -1)
+            {
+                compactFld[idx] = fld[i];
+            }
+        }
+    }
+
+
+    // Pull all data
+    mapDistributeBase::distribute<Type>(compactFld);
+
+    // Rewrite to output
+    fld = UIndirectList<Type>(compactFld, compactOutputMapping);
+}
+
+
+template<class Type>
+void Foam::faMeshBoundaryHalo::distributeSparse
+(
+    List<Type>& fld,
+    const labelUList& sparseInputLocations
+) const
+{
+    this->distributeSparse(fld, sparseInputLocations, boundaryToCompact_);
+}
+
+
+template<class Type>
+void Foam::faMeshBoundaryHalo::distributeSparse(List<Type>& fld) const
+{
+    this->distributeSparse(fld, inputMeshFaces_, boundaryToCompact_);
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteArea/faMesh/faMeshTopology.C b/src/finiteArea/faMesh/faMeshTopology.C
index 179054ac0b..938dbe9dbd 100644
--- a/src/finiteArea/faMesh/faMeshTopology.C
+++ b/src/finiteArea/faMesh/faMeshTopology.C
@@ -26,6 +26,7 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "faMesh.H"
+#include "faMeshBoundaryHalo.H"
 #include "globalMeshData.H"
 #include "indirectPrimitivePatch.H"
 #include "edgeHashes.H"
@@ -800,4 +801,130 @@ Foam::List<Foam::labelPair> Foam::faMesh::boundaryProcSizes() const
 }
 
 
+const Foam::faMeshBoundaryHalo& Foam::faMesh::boundaryHaloMap() const
+{
+    if (!haloMapPtr_)
+    {
+        haloMapPtr_.reset(new faMeshBoundaryHalo(*this));
+    }
+
+    return *haloMapPtr_;
+}
+
+
+void Foam::faMesh::calcHaloFaceGeometry() const
+{
+    if (haloFaceCentresPtr_ || haloFaceNormalsPtr_)
+    {
+        FatalErrorInFunction
+            << "Halo centres/normals already calculated"
+            << exit(FatalError);
+    }
+
+    DebugInFunction
+        << "Calculating halo face centres/normals" << endl;
+
+    const faceList& faces = mesh().faces();
+    const pointField& points = mesh().points();
+
+    const faMeshBoundaryHalo& halo = boundaryHaloMap();
+
+    const labelList& inputFaceIds = halo.inputMeshFaces();
+
+    haloFaceCentresPtr_.reset(new pointField);
+    haloFaceNormalsPtr_.reset(new vectorField);
+
+    auto& centres = *haloFaceCentresPtr_;
+    auto& normals = *haloFaceNormalsPtr_;
+
+    centres.resize(inputFaceIds.size());
+    normals.resize(inputFaceIds.size());
+
+    // My values
+    forAll(inputFaceIds, i)
+    {
+        const face& f = faces[inputFaceIds[i]];
+
+        centres[i] = f.centre(points);
+        normals[i] = f.unitNormal(points);
+    }
+
+    // Swap information and resize
+    halo.distributeSparse(centres);
+    halo.distributeSparse(normals);
+}
+
+
+const Foam::pointField& Foam::faMesh::haloFaceCentres() const
+{
+    if (!haloFaceCentresPtr_ || !haloFaceNormalsPtr_)
+    {
+        calcHaloFaceGeometry();
+    }
+
+    return *haloFaceCentresPtr_;
+}
+
+
+const Foam::vectorField& Foam::faMesh::haloFaceNormals() const
+{
+    if (!haloFaceCentresPtr_ || !haloFaceNormalsPtr_)
+    {
+        calcHaloFaceGeometry();
+    }
+
+    return *haloFaceNormalsPtr_;
+}
+
+
+Foam::tmp<Foam::pointField>
+Foam::faMesh::haloFaceCentres(const label patchi) const
+{
+    if (patchi < 0 || patchi >= boundary().size())
+    {
+        FatalErrorInFunction
+            << "Patch " << patchi << " is out-of-range 0.."
+            << (boundary().size()-1) << nl
+            << exit(FatalError);
+    }
+
+    return tmp<pointField>::New
+    (
+        List<point>
+        (
+            boundarySubset
+            (
+                this->haloFaceCentres(),
+                boundary()[patchi].edgeLabels()
+            )
+        )
+    );
+}
+
+
+Foam::tmp<Foam::vectorField>
+Foam::faMesh::haloFaceNormals(const label patchi) const
+{
+    if (patchi < 0 || patchi >= boundary().size())
+    {
+        FatalErrorInFunction
+            << "Patch " << patchi << " is out-of-range 0.."
+            << (boundary().size()-1) << nl
+            << exit(FatalError);
+    }
+
+    return tmp<vectorField>::New
+    (
+        List<vector>
+        (
+            boundarySubset
+            (
+                this->haloFaceNormals(),
+                boundary()[patchi].edgeLabels()
+            )
+        )
+    );
+}
+
+
 // ************************************************************************* //