From d67a1df92eacc0f69912a44a817cb36f72ad2e0d Mon Sep 17 00:00:00 2001
From: Henry Weller <http://cfd.direct>
Date: Tue, 21 Jun 2016 14:16:18 +0100
Subject: [PATCH] mapFieldsPar: updated to enable mapping from source patches
 (instead of recreating)

  - patchFields now get mapped (instead of created)
  - with -consistent it now maps all patches except for processor ones (they are
    the only ones that are processor-local)
  - all constraint patches get evaluated after mapping to bring them up to date.

Patch contributed by Mattijs Janssens
---
 .../preProcessing/mapFieldsPar/MapVolFields.H |  95 +++++++++-
 .../preProcessing/mapFieldsPar/mapFieldsPar.C |  10 +-
 .../preProcessing/mapFieldsPar/setTimeIndex.H |   8 +
 .../fvPatchFields/fvPatchField/fvPatchField.C |  49 ++---
 .../calcMethod/direct/directMethod.C          |  24 +--
 .../distributedWeightedFvPatchFieldMapper.H   | 170 ++++++++++++++++++
 src/sampling/meshToMesh/meshToMesh.C          | 116 +++++++++++-
 src/sampling/meshToMesh/meshToMesh.H          | 107 +++++++++++
 src/sampling/meshToMesh/meshToMeshTemplates.C | 148 +++++++++------
 9 files changed, 617 insertions(+), 110 deletions(-)
 create mode 100644 src/sampling/meshToMesh/distributedWeightedFvPatchFieldMapper.H

diff --git a/applications/utilities/preProcessing/mapFieldsPar/MapVolFields.H b/applications/utilities/preProcessing/mapFieldsPar/MapVolFields.H
index ac5c18b077..19791dac16 100644
--- a/applications/utilities/preProcessing/mapFieldsPar/MapVolFields.H
+++ b/applications/utilities/preProcessing/mapFieldsPar/MapVolFields.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -35,6 +35,88 @@ License
 namespace Foam
 {
 
+template<class Type>
+void evaluateConstraintTypes(GeometricField<Type, fvPatchField, volMesh>& fld)
+{
+    typename GeometricField<Type, fvPatchField, volMesh>::
+        Boundary& fldBf = fld.boundaryFieldRef();
+
+    if
+    (
+        Pstream::defaultCommsType == Pstream::blocking
+     || Pstream::defaultCommsType == Pstream::nonBlocking
+    )
+    {
+        label nReq = Pstream::nRequests();
+
+        forAll(fldBf, patchi)
+        {
+            fvPatchField<Type>& tgtField = fldBf[patchi];
+
+            if
+            (
+                tgtField.type() == tgtField.patch().patch().type()
+             && polyPatch::constraintType(tgtField.patch().patch().type())
+            )
+            {
+                tgtField.initEvaluate(Pstream::defaultCommsType);
+            }
+        }
+
+        // Block for any outstanding requests
+        if
+        (
+            Pstream::parRun()
+         && Pstream::defaultCommsType == Pstream::nonBlocking
+        )
+        {
+            Pstream::waitRequests(nReq);
+        }
+
+        forAll(fldBf, patchi)
+        {
+            fvPatchField<Type>& tgtField = fldBf[patchi];
+
+            if
+            (
+                tgtField.type() == tgtField.patch().patch().type()
+             && polyPatch::constraintType(tgtField.patch().patch().type())
+            )
+            {
+                tgtField.evaluate(Pstream::defaultCommsType);
+            }
+        }
+    }
+    else if (Pstream::defaultCommsType == Pstream::scheduled)
+    {
+        const lduSchedule& patchSchedule =
+            fld.mesh().globalData().patchSchedule();
+
+        forAll(patchSchedule, patchEvali)
+        {
+            label patchi = patchSchedule[patchEvali].patch;
+            fvPatchField<Type>& tgtField = fldBf[patchi];
+
+            if
+            (
+                tgtField.type() == tgtField.patch().patch().type()
+             && polyPatch::constraintType(tgtField.patch().patch().type())
+            )
+            {
+                if (patchSchedule[patchEvali].init)
+                {
+                    tgtField.initEvaluate(Pstream::scheduled);
+                }
+                else
+                {
+                    tgtField.evaluate(Pstream::scheduled);
+                }
+            }
+        }
+    }
+}
+
+
 template<class Type, class CombineOp>
 void MapVolFields
 (
@@ -57,8 +139,6 @@ void MapVolFields
 
         if (selectedFields.empty() || selectedFields.found(fieldName))
         {
-            Info<< "    interpolating " << fieldName << endl;
-
             const fieldType fieldSource(*fieldIter(), meshSource);
 
             IOobject targetIO
@@ -71,14 +151,21 @@ void MapVolFields
 
             if (targetIO.headerOk())
             {
+                Info<< "    interpolating onto existing field "
+                    << fieldName << endl;
                 fieldType fieldTarget(targetIO, meshTarget);
 
                 interp.mapSrcToTgt(fieldSource, cop, fieldTarget);
 
+                evaluateConstraintTypes(fieldTarget);
+
                 fieldTarget.write();
             }
             else
             {
+                Info<< "    creating new field "
+                    << fieldName << endl;
+
                 targetIO.readOpt() = IOobject::NO_READ;
 
                 tmp<fieldType>
@@ -86,6 +173,8 @@ void MapVolFields
 
                 fieldType fieldTarget(targetIO, tfieldTarget);
 
+                evaluateConstraintTypes(fieldTarget);
+
                 fieldTarget.write();
             }
         }
diff --git a/applications/utilities/preProcessing/mapFieldsPar/mapFieldsPar.C b/applications/utilities/preProcessing/mapFieldsPar/mapFieldsPar.C
index ddadaffb4f..fc10a7f8ea 100644
--- a/applications/utilities/preProcessing/mapFieldsPar/mapFieldsPar.C
+++ b/applications/utilities/preProcessing/mapFieldsPar/mapFieldsPar.C
@@ -263,11 +263,11 @@ int main(int argc, char *argv[])
 
     if (!consistent)
     {
-        IOdictionary mapFieldsParDict
+        IOdictionary mapFieldsDict
         (
             IOobject
             (
-                "mapFieldsParDict",
+                "mapFieldsDict",
                 runTimeTarget.system(),
                 runTimeTarget,
                 IOobject::MUST_READ_IF_MODIFIED,
@@ -276,8 +276,8 @@ int main(int argc, char *argv[])
             )
         );
 
-        mapFieldsParDict.lookup("patchMap") >> patchMap;
-        mapFieldsParDict.lookup("cuttingPatches") >>  cuttingPatches;
+        mapFieldsDict.lookup("patchMap") >> patchMap;
+        mapFieldsDict.lookup("cuttingPatches") >>  cuttingPatches;
     }
 
     #include "setTimeIndex.H"
@@ -326,7 +326,7 @@ int main(int argc, char *argv[])
             meshSource,
             meshTarget,
             patchMap,
-            addProcessorPatches(meshTarget, cuttingPatches),
+            cuttingPatches,
             mapMethod,
             subtract,
             selectedFields,
diff --git a/applications/utilities/preProcessing/mapFieldsPar/setTimeIndex.H b/applications/utilities/preProcessing/mapFieldsPar/setTimeIndex.H
index a468d50648..e1ca52e097 100644
--- a/applications/utilities/preProcessing/mapFieldsPar/setTimeIndex.H
+++ b/applications/utilities/preProcessing/mapFieldsPar/setTimeIndex.H
@@ -1,5 +1,13 @@
 {
     instantList sourceTimes = runTimeSource.times();
+
+    if (sourceTimes.empty())
+    {
+        FatalErrorInFunction << "No result times in source "
+            << runTimeSource.caseName()
+            << exit(FatalError);
+    }
+
     label sourceTimeIndex = runTimeSource.timeIndex();
     if (args.optionFound("sourceTime"))
     {
diff --git a/src/finiteVolume/fields/fvPatchFields/fvPatchField/fvPatchField.C b/src/finiteVolume/fields/fvPatchFields/fvPatchField/fvPatchField.C
index 5ecbd2b1fa..ea2935e73d 100644
--- a/src/finiteVolume/fields/fvPatchFields/fvPatchField/fvPatchField.C
+++ b/src/finiteVolume/fields/fvPatchFields/fvPatchField/fvPatchField.C
@@ -153,7 +153,7 @@ Foam::fvPatchField<Type>::fvPatchField
     patchType_(ptf.patchType_)
 {
     // For unmapped faces set to internal field value (zero-gradient)
-    if (notNull(iF) && iF.size())
+    if (notNull(iF) && mapper.hasUnmapped())
     {
         fvPatchField<Type>::operator=(this->patchInternalField());
     }
@@ -243,7 +243,7 @@ void Foam::fvPatchField<Type>::autoMap
 {
     Field<Type>& f = *this;
 
-    if (!this->size())
+    if (!this->size() && !mapper.distributed())
     {
         f.setSize(mapper.size());
         if (f.size())
@@ -257,38 +257,39 @@ void Foam::fvPatchField<Type>::autoMap
         Field<Type>::autoMap(mapper);
 
         // For unmapped faces set to internal field value (zero-gradient)
-        if
-        (
-            mapper.direct()
-         && notNull(mapper.directAddressing())
-         && mapper.directAddressing().size()
-        )
+        if (mapper.hasUnmapped())
         {
             Field<Type> pif(this->patchInternalField());
 
-            const labelList& mapAddressing = mapper.directAddressing();
-
-            forAll(mapAddressing, i)
+            if
+            (
+                mapper.direct()
+             && notNull(mapper.directAddressing())
+             && mapper.directAddressing().size()
+            )
             {
-                if (mapAddressing[i] < 0)
+                const labelList& mapAddressing = mapper.directAddressing();
+
+                forAll(mapAddressing, i)
                 {
-                    f[i] = pif[i];
+                    if (mapAddressing[i] < 0)
+                    {
+                        f[i] = pif[i];
+                    }
                 }
             }
-        }
-        else if (!mapper.direct() && mapper.addressing().size())
-        {
-            Field<Type> pif(this->patchInternalField());
-
-            const labelListList& mapAddressing = mapper.addressing();
-
-            forAll(mapAddressing, i)
+            else if (!mapper.direct() && mapper.addressing().size())
             {
-                const labelList& localAddrs = mapAddressing[i];
+                const labelListList& mapAddressing = mapper.addressing();
 
-                if (!localAddrs.size())
+                forAll(mapAddressing, i)
                 {
-                    f[i] = pif[i];
+                    const labelList& localAddrs = mapAddressing[i];
+
+                    if (!localAddrs.size())
+                    {
+                        f[i] = pif[i];
+                    }
                 }
             }
         }
diff --git a/src/sampling/meshToMesh/calcMethod/direct/directMethod.C b/src/sampling/meshToMesh/calcMethod/direct/directMethod.C
index 46d5706b75..4fc9157b8d 100644
--- a/src/sampling/meshToMesh/calcMethod/direct/directMethod.C
+++ b/src/sampling/meshToMesh/calcMethod/direct/directMethod.C
@@ -72,21 +72,15 @@ bool Foam::directMethod::findInitialSeeds
 
         if (mapFlag[srcI])
         {
-            const pointField
-                pts(srcCells[srcI].points(srcFaces, srcPts).xfer());
-
-            forAll(pts, ptI)
-            {
-                const point& pt = pts[ptI];
-                label tgtI = tgt_.cellTree().findInside(pt);
+            const point srcCtr(srcCells[srcI].centre(srcPts, srcFaces));
+            label tgtI = tgt_.cellTree().findInside(srcCtr);
 
                 if (tgtI != -1 && intersect(srcI, tgtI))
                 {
                     srcSeedI = srcI;
                     tgtSeedI = tgtI;
 
-                    return true;
-                }
+                return true;
             }
         }
     }
@@ -178,8 +172,6 @@ void Foam::directMethod::appendToDirectSeeds
     const labelList& srcNbr = src_.cellCells()[srcSeedI];
     const labelList& tgtNbr = tgt_.cellCells()[tgtSeedI];
 
-    const vectorField& srcCentre = src_.cellCentres();
-
     forAll(srcNbr, i)
     {
         label srcI = srcNbr[i];
@@ -194,15 +186,7 @@ void Foam::directMethod::appendToDirectSeeds
             {
                 label tgtI = tgtNbr[j];
 
-                if
-                (
-                    tgt_.pointInCell
-                    (
-                        srcCentre[srcI],
-                        tgtI,
-                        polyMesh::FACE_PLANES
-                    )
-                )
+                if (intersect(srcI, tgtI))
                 {
                     // new match - append to lists
                     found = true;
diff --git a/src/sampling/meshToMesh/distributedWeightedFvPatchFieldMapper.H b/src/sampling/meshToMesh/distributedWeightedFvPatchFieldMapper.H
new file mode 100644
index 0000000000..537f7d8de5
--- /dev/null
+++ b/src/sampling/meshToMesh/distributedWeightedFvPatchFieldMapper.H
@@ -0,0 +1,170 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2015-2016 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+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::distributedWeightedFvPatchFieldMapper
+
+Description
+    FieldMapper with weighted mapping from (optionally remote) quantities.
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef distributedWeightedFvPatchFieldMapper_H
+#define distributedWeightedFvPatchFieldMapper_H
+
+#include "fvPatchFieldMapper.H"
+#include "mapDistributeBase.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+           Class distributedWeightedFvPatchFieldMapper Declaration
+\*---------------------------------------------------------------------------*/
+
+class distributedWeightedFvPatchFieldMapper
+:
+    public fvPatchFieldMapper
+{
+    const label singlePatchProc_;
+
+    const mapDistributeBase* distMapPtr_;
+
+    const labelListList& addressing_;
+
+    const scalarListList& weights_;
+
+    bool hasUnmapped_;
+
+public:
+
+    // Constructors
+
+        //- Construct given addressing
+        distributedWeightedFvPatchFieldMapper
+        (
+            const label singlePatchProc,
+            const mapDistributeBase* distMapPtr,
+            const labelListList& addressing,
+            const scalarListList& weights
+        )
+        :
+            singlePatchProc_(singlePatchProc),
+            distMapPtr_(distMapPtr),
+            addressing_(addressing),
+            weights_(weights),
+            hasUnmapped_(false)
+        {
+            forAll(addressing_, i)
+            {
+                if (addressing_[i].size() == 0)
+                {
+                    hasUnmapped_ = true;
+                }
+            }
+
+            if ((singlePatchProc_ == -1) != (distMapPtr_ != NULL))
+            {
+                FatalErrorIn
+                (
+                    "distributedWeightedFvPatchFieldMapper::"
+                    "distributedWeightedFvPatchFieldMapper(..)"
+                )   << "Supply a mapDistributeBase if and only if "
+                    << "singlePatchProc is -1"
+                    << " singlePatchProc_:" << singlePatchProc_
+                    << " distMapPtr_:" << (distMapPtr_ != NULL)
+                    << exit(FatalError);
+            }
+        }
+
+    //- Destructor
+    virtual ~distributedWeightedFvPatchFieldMapper()
+    {}
+
+
+    // Member Functions
+
+        virtual label size() const
+        {
+            if (distributed())
+            {
+                return distributeMap().constructSize();
+            }
+            else
+            {
+                return addressing().size();
+            }
+        }
+
+        virtual bool direct() const
+        {
+            return false;
+        }
+
+        virtual bool distributed() const
+        {
+            return singlePatchProc_ == -1;
+        }
+
+        virtual const mapDistributeBase& distributeMap() const
+        {
+            if (!distMapPtr_)
+            {
+                FatalErrorIn
+                (
+                    "distributedWeightedFvPatchFieldMapper::"
+                    "distributeMap()"
+                )   << "Cannot ask for distributeMap on a non-distributed"
+                    << " mapper" << exit(FatalError);
+            }
+            return *distMapPtr_;
+        }
+
+        virtual bool hasUnmapped() const
+        {
+            return hasUnmapped_;
+        }
+
+        virtual const labelListList& addressing() const
+        {
+            return addressing_;
+        }
+
+        virtual const scalarListList& weights() const
+        {
+            return weights_;
+        }
+
+};
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/sampling/meshToMesh/meshToMesh.C b/src/sampling/meshToMesh/meshToMesh.C
index b49ae2a7b0..0c2efaaeed 100644
--- a/src/sampling/meshToMesh/meshToMesh.C
+++ b/src/sampling/meshToMesh/meshToMesh.C
@@ -53,6 +53,116 @@ namespace Foam
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
+template<>
+void Foam::meshToMesh::mapAndOpSrcToTgt
+(
+    const AMIPatchToPatchInterpolation& AMI,
+    const Field<scalar>& srcField,
+    Field<scalar>& tgtField,
+    const plusEqOp<scalar>& cop
+) const
+{}
+
+
+template<>
+void Foam::meshToMesh::mapAndOpSrcToTgt
+(
+    const AMIPatchToPatchInterpolation& AMI,
+    const Field<vector>& srcField,
+    Field<vector>& tgtField,
+    const plusEqOp<vector>& cop
+) const
+{}
+
+
+template<>
+void Foam::meshToMesh::mapAndOpSrcToTgt
+(
+    const AMIPatchToPatchInterpolation& AMI,
+    const Field<sphericalTensor>& srcField,
+    Field<sphericalTensor>& tgtField,
+    const plusEqOp<sphericalTensor>& cop
+) const
+{}
+
+
+template<>
+void Foam::meshToMesh::mapAndOpSrcToTgt
+(
+    const AMIPatchToPatchInterpolation& AMI,
+    const Field<symmTensor>& srcField,
+    Field<symmTensor>& tgtField,
+    const plusEqOp<symmTensor>& cop
+) const
+{}
+
+
+template<>
+void Foam::meshToMesh::mapAndOpSrcToTgt
+(
+    const AMIPatchToPatchInterpolation& AMI,
+    const Field<tensor>& srcField,
+    Field<tensor>& tgtField,
+    const plusEqOp<tensor>& cop
+) const
+{}
+
+
+template<>
+void Foam::meshToMesh::mapAndOpTgtToSrc
+(
+    const AMIPatchToPatchInterpolation& AMI,
+    Field<scalar>& srcField,
+    const Field<scalar>& tgtField,
+    const plusEqOp<scalar>& cop
+) const
+{}
+
+
+template<>
+void Foam::meshToMesh::mapAndOpTgtToSrc
+(
+    const AMIPatchToPatchInterpolation& AMI,
+    Field<vector>& srcField,
+    const Field<vector>& tgtField,
+    const plusEqOp<vector>& cop
+) const
+{}
+
+
+template<>
+void Foam::meshToMesh::mapAndOpTgtToSrc
+(
+    const AMIPatchToPatchInterpolation& AMI,
+    Field<sphericalTensor>& srcField,
+    const Field<sphericalTensor>& tgtField,
+    const plusEqOp<sphericalTensor>& cop
+) const
+{}
+
+
+template<>
+void Foam::meshToMesh::mapAndOpTgtToSrc
+(
+    const AMIPatchToPatchInterpolation& AMI,
+    Field<symmTensor>& srcField,
+    const Field<symmTensor>& tgtField,
+    const plusEqOp<symmTensor>& cop
+) const
+{}
+
+
+template<>
+void Foam::meshToMesh::mapAndOpTgtToSrc
+(
+    const AMIPatchToPatchInterpolation& AMI,
+    Field<tensor>& srcField,
+    const Field<tensor>& tgtField,
+    const plusEqOp<tensor>& cop
+) const
+{}
+
+
 Foam::labelList Foam::meshToMesh::maskCells
 (
     const polyMesh& src,
@@ -443,7 +553,11 @@ void Foam::meshToMesh::constructNoCuttingPatches
         forAll(srcBM, patchi)
         {
             const polyPatch& pp = srcBM[patchi];
-            if (!polyPatch::constraintType(pp.type()))
+
+            // We want to map all the global patches, including constraint
+            // patches (since they might have mappable properties, e.g.
+            // jumpCyclic). We'll fix the value afterwards.
+            if (!isA<processorPolyPatch>(pp))
             {
                 srcPatchID.append(pp.index());
 
diff --git a/src/sampling/meshToMesh/meshToMesh.H b/src/sampling/meshToMesh/meshToMesh.H
index 6d5e5cd058..1c114e34b5 100644
--- a/src/sampling/meshToMesh/meshToMesh.H
+++ b/src/sampling/meshToMesh/meshToMesh.H
@@ -129,6 +129,28 @@ private:
         template<class Type>
         void add(UList<Type>& fld, const label offset) const;
 
+        //- Helper function to interpolate patch field. Template
+        //  specialisations  below
+        template<class Type, class CombineOp>
+        void mapAndOpSrcToTgt
+        (
+            const AMIPatchToPatchInterpolation& AMI,
+            const Field<Type>& srcField,
+            Field<Type>& tgtField,
+            const CombineOp& cop
+        ) const;
+
+        //- Helper function to interpolate patch field. Template
+        //  specialisations  below
+        template<class Type, class CombineOp>
+        void mapAndOpTgtToSrc
+        (
+            const AMIPatchToPatchInterpolation& AMI,
+            Field<Type>& srcField,
+            const Field<Type>& tgtField,
+            const CombineOp& cop
+        ) const;
+
         //- Return src cell IDs for the overlap region
         labelList maskCells(const polyMesh& src, const polyMesh& tgt) const;
 
@@ -524,6 +546,91 @@ public:
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
+// Disable fvPatchField value override after rmap
+
+template<>
+void meshToMesh::mapAndOpSrcToTgt
+(
+    const AMIPatchToPatchInterpolation& AMI,
+    const Field<scalar>& srcField,
+    Field<scalar>& tgtField,
+    const plusEqOp<scalar>& cop
+) const;
+template<>
+void meshToMesh::mapAndOpSrcToTgt
+(
+    const AMIPatchToPatchInterpolation& AMI,
+    const Field<vector>& srcField,
+    Field<vector>& tgtField,
+    const plusEqOp<vector>& cop
+) const;
+template<>
+void meshToMesh::mapAndOpSrcToTgt
+(
+    const AMIPatchToPatchInterpolation& AMI,
+    const Field<sphericalTensor>& srcField,
+    Field<sphericalTensor>& tgtField,
+    const plusEqOp<sphericalTensor>& cop
+) const;
+template<>
+void meshToMesh::mapAndOpSrcToTgt
+(
+    const AMIPatchToPatchInterpolation& AMI,
+    const Field<symmTensor>& srcField,
+    Field<symmTensor>& tgtField,
+    const plusEqOp<symmTensor>& cop
+) const;
+template<>
+void meshToMesh::mapAndOpSrcToTgt
+(
+    const AMIPatchToPatchInterpolation& AMI,
+    const Field<tensor>& srcField,
+    Field<tensor>& tgtField,
+    const plusEqOp<tensor>& cop
+) const;
+
+
+template<>
+void meshToMesh::mapAndOpTgtToSrc
+(
+    const AMIPatchToPatchInterpolation& AMI,
+    Field<scalar>& srcField,
+    const Field<scalar>& tgtField,
+    const plusEqOp<scalar>& cop
+) const;
+template<>
+void meshToMesh::mapAndOpTgtToSrc
+(
+    const AMIPatchToPatchInterpolation& AMI,
+    Field<vector>& srcField,
+    const Field<vector>& tgtField,
+    const plusEqOp<vector>& cop
+) const;
+template<>
+void meshToMesh::mapAndOpTgtToSrc
+(
+    const AMIPatchToPatchInterpolation& AMI,
+    Field<sphericalTensor>& srcField,
+    const Field<sphericalTensor>& tgtField,
+    const plusEqOp<sphericalTensor>& cop
+) const;
+template<>
+void meshToMesh::mapAndOpTgtToSrc
+(
+    const AMIPatchToPatchInterpolation& AMI,
+    Field<symmTensor>& srcField,
+    const Field<symmTensor>& tgtField,
+    const plusEqOp<symmTensor>& cop
+) const;
+template<>
+void meshToMesh::mapAndOpTgtToSrc
+(
+    const AMIPatchToPatchInterpolation& AMI,
+    Field<tensor>& srcField,
+    const Field<tensor>& tgtField,
+    const plusEqOp<tensor>& cop
+) const;
+
 } // End namespace Foam
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
diff --git a/src/sampling/meshToMesh/meshToMeshTemplates.C b/src/sampling/meshToMesh/meshToMeshTemplates.C
index 446b46ecf3..7f2038a371 100644
--- a/src/sampling/meshToMesh/meshToMeshTemplates.C
+++ b/src/sampling/meshToMesh/meshToMeshTemplates.C
@@ -27,7 +27,7 @@ License
 #include "volFields.H"
 #include "directFvPatchFieldMapper.H"
 #include "calculatedFvPatchField.H"
-#include "weightedFvPatchFieldMapper.H"
+#include "distributedWeightedFvPatchFieldMapper.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -317,6 +317,27 @@ Foam::tmp<Foam::Field<Type>> Foam::meshToMesh::mapTgtToSrc
 }
 
 
+template<class Type, class CombineOp>
+void Foam::meshToMesh::mapAndOpSrcToTgt
+(
+    const AMIPatchToPatchInterpolation& AMI,
+    const Field<Type>& srcField,
+    Field<Type>& tgtField,
+    const CombineOp& cop
+) const
+{
+    tgtField = pTraits<Type>::zero;
+
+    AMI.interpolateToTarget
+    (
+        srcField,
+        multiplyWeightedOp<Type, CombineOp>(cop),
+        tgtField,
+        UList<Type>::null()
+    );
+}
+
+
 template<class Type, class CombineOp>
 void Foam::meshToMesh::mapSrcToTgt
 (
@@ -340,40 +361,35 @@ void Foam::meshToMesh::mapSrcToTgt
         const fvPatchField<Type>& srcField = field.boundaryField()[srcPatchi];
         fvPatchField<Type>& tgtField = resultBf[tgtPatchi];
 
-        // 2.3 does not do distributed mapping yet so only do if
-        // running on single processor
-        if (AMIList[i].singlePatchProc() != -1)
-        {
-            // Clone and map (since rmap does not do general mapping)
-            tmp<fvPatchField<Type>> tnewTgt
+        // Clone and map (since rmap does not do general mapping)
+        tmp<fvPatchField<Type>> tnewTgt
+        (
+            fvPatchField<Type>::New
             (
-                fvPatchField<Type>::New
+                srcField,
+                tgtField.patch(),
+                result(),
+                distributedWeightedFvPatchFieldMapper
                 (
-                    srcField,
-                    tgtField.patch(),
-                    result(),
-                    weightedFvPatchFieldMapper
+                    AMIList[i].singlePatchProc(),
                     (
-                        AMIList[i].tgtAddress(),
-                        AMIList[i].tgtWeights()
-                    )
+                        AMIList[i].singlePatchProc() == -1
+                      ? &AMIList[i].srcMap()
+                      : NULL
+                    ),
+                    AMIList[i].tgtAddress(),
+                    AMIList[i].tgtWeights()
                 )
-            );
-
-            // Transfer all mapped quantities (value and e.g. gradient) onto
-            // tgtField. Value will get overwritten below.
-            tgtField.rmap(tnewTgt(), identity(tgtField.size()));
-        }
+            )
+        );
 
-        tgtField == Type(Zero);
+        // Transfer all mapped quantities (value and e.g. gradient) onto
+        // tgtField. Value will get overwritten below.
+        tgtField.rmap(tnewTgt(), identity(tgtField.size()));
 
-        AMIList[i].interpolateToTarget
-        (
-            srcField,
-            multiplyWeightedOp<Type, CombineOp>(cop),
-            tgtField,
-            UList<Type>::null()
-        );
+        // Override value to account for CombineOp (note: is dummy template
+        // specialisation for plusEqOp)
+        mapAndOpSrcToTgt(AMIList[i], srcField, tgtField, cop);
     }
 
     forAll(cuttingPatches_, i)
@@ -403,7 +419,7 @@ Foam::meshToMesh::mapSrcToTgt
 
     PtrList<fvPatchField<Type>> tgtPatchFields(tgtBm.size());
 
-    // constuct tgt boundary patch types as copy of 'field' boundary types
+    // construct tgt boundary patch types as copy of 'field' boundary types
     // note: this will provide place holders for fields with additional
     // entries, but these values will need to be reset
     forAll(tgtPatchID_, i)
@@ -509,6 +525,27 @@ Foam::meshToMesh::mapSrcToTgt
 }
 
 
+template<class Type, class CombineOp>
+void Foam::meshToMesh::mapAndOpTgtToSrc
+(
+    const AMIPatchToPatchInterpolation& AMI,
+    Field<Type>& srcField,
+    const Field<Type>& tgtField,
+    const CombineOp& cop
+) const
+{
+    srcField = pTraits<Type>::zero;
+
+    AMI.interpolateToSource
+    (
+        tgtField,
+        multiplyWeightedOp<Type, CombineOp>(cop),
+        srcField,
+        UList<Type>::null()
+    );
+}
+
+
 template<class Type, class CombineOp>
 void Foam::meshToMesh::mapTgtToSrc
 (
@@ -529,40 +566,37 @@ void Foam::meshToMesh::mapTgtToSrc
         fvPatchField<Type>& srcField = result.boundaryField()[srcPatchi];
         const fvPatchField<Type>& tgtField = field.boundaryField()[tgtPatchi];
 
-        // 2.3 does not do distributed mapping yet so only do if
-        // running on single processor
-        if (AMIList[i].singlePatchProc() != -1)
-        {
-            // Clone and map (since rmap does not do general mapping)
-            tmp<fvPatchField<Type>> tnewSrc
+
+        // Clone and map (since rmap does not do general mapping)
+        tmp<fvPatchField<Type>> tnewSrc
+        (
+            fvPatchField<Type>::New
             (
-                fvPatchField<Type>::New
+                tgtField,
+                srcField.patch(),
+                result(),
+                distributedWeightedFvPatchFieldMapper
                 (
-                    tgtField,
-                    srcField.patch(),
-                    result(),
-                    weightedFvPatchFieldMapper
+                    AMIList[i].singlePatchProc(),
                     (
-                        AMIList[i].srcAddress(),
-                        AMIList[i].srcWeights()
-                    )
+                        AMIList[i].singlePatchProc() == -1
+                      ? &AMIList[i].tgtMap()
+                      : NULL
+                    ),
+                    AMIList[i].srcAddress(),
+                    AMIList[i].srcWeights()
                 )
-            );
+            )
+        );
 
-            // Transfer all mapped quantities (value and e.g. gradient) onto
-            // srcField. Value will get overwritten below
-            srcField.rmap(tnewSrc(), identity(srcField.size()));
-        }
+        // Transfer all mapped quantities (value and e.g. gradient) onto
+        // srcField. Value will get overwritten below
+        srcField.rmap(tnewSrc(), identity(srcField.size()));
 
-        srcField == Type(Zero);
 
-        AMIList[i].interpolateToSource
-        (
-            tgtField,
-            multiplyWeightedOp<Type, CombineOp>(cop),
-            srcField,
-            UList<Type>::null()
-        );
+        // Override value to account for CombineOp (could be dummy for
+        // plusEqOp)
+        mapAndOpTgtToSrc(AMIList[i], srcField, tgtField, cop);
     }
 
     forAll(cuttingPatches_, i)
-- 
GitLab