From d52546b3ce776dbfcfc24eb3ffc30e80064d3b4c Mon Sep 17 00:00:00 2001
From: Henry Weller <http://cfd.direct>
Date: Wed, 6 Apr 2016 14:23:18 +0100
Subject: [PATCH] surfaceInterpolationScheme: Added dotInterpolate
 member-function

dotInterpolate interpolates the field and "dots" the resulting
face-values with the vector field provided which removes the need to
create a temporary field for the interpolate.  This reduces the peak
storage of OpenFOAM caused by the divergence of the gradient of vector
fields, improves memory management and under some conditions decreases
run-time.

This development is based on a patch contributed by Paul Edwards, Intel.
---
 .../gaussDivScheme/gaussDivScheme.C           |   2 +-
 .../correctedSnGrad/correctedSnGrad.C         |   4 +-
 .../surfaceInterpolationScheme.C              | 132 ++++++++++++++++--
 .../surfaceInterpolationScheme.H              |  72 +++++++++-
 .../surfaceInterpolationSchemes.C             |  42 +++++-
 5 files changed, 235 insertions(+), 17 deletions(-)

diff --git a/src/finiteVolume/finiteVolume/divSchemes/gaussDivScheme/gaussDivScheme.C b/src/finiteVolume/finiteVolume/divSchemes/gaussDivScheme/gaussDivScheme.C
index 9721a04eca..308f0daff7 100644
--- a/src/finiteVolume/finiteVolume/divSchemes/gaussDivScheme/gaussDivScheme.C
+++ b/src/finiteVolume/finiteVolume/divSchemes/gaussDivScheme/gaussDivScheme.C
@@ -58,7 +58,7 @@ gaussDivScheme<Type>::fvcDiv
     (
         fvc::surfaceIntegrate
         (
-            this->mesh_.Sf() & this->tinterpScheme_().interpolate(vf)
+            this->tinterpScheme_().dotInterpolate(this->mesh_.Sf(), vf)
         )
     );
 
diff --git a/src/finiteVolume/finiteVolume/snGradSchemes/correctedSnGrad/correctedSnGrad.C b/src/finiteVolume/finiteVolume/snGradSchemes/correctedSnGrad/correctedSnGrad.C
index ea58ee0249..305743076c 100644
--- a/src/finiteVolume/finiteVolume/snGradSchemes/correctedSnGrad/correctedSnGrad.C
+++ b/src/finiteVolume/finiteVolume/snGradSchemes/correctedSnGrad/correctedSnGrad.C
@@ -50,9 +50,9 @@ Foam::fv::correctedSnGrad<Type>::fullGradCorrection
 
     // construct GeometricField<Type, fvsPatchField, surfaceMesh>
     tmp<GeometricField<Type, fvsPatchField, surfaceMesh>> tssf =
-        mesh.nonOrthCorrectionVectors()
-      & linear<typename outerProduct<vector, Type>::type>(mesh).interpolate
+        linear<typename outerProduct<vector, Type>::type>(mesh).dotInterpolate
         (
+            mesh.nonOrthCorrectionVectors(),
             gradScheme<Type>::New
             (
                 mesh,
diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationScheme.C b/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationScheme.C
index d871edb6e2..4f6a79c9a2 100644
--- a/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationScheme.C
+++ b/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationScheme.C
@@ -26,6 +26,7 @@ License
 #include "surfaceInterpolationScheme.H"
 #include "volFields.H"
 #include "surfaceFields.H"
+#include "geometricOneField.H"
 #include "coupledFvPatchField.H"
 
 // * * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * //
@@ -215,9 +216,19 @@ Foam::surfaceInterpolationScheme<Type>::interpolate
 
 
 template<class Type>
-Foam::tmp<Foam::GeometricField<Type, Foam::fvsPatchField, Foam::surfaceMesh>>
-Foam::surfaceInterpolationScheme<Type>::interpolate
+template<class SFType>
+Foam::tmp
+<
+    Foam::GeometricField
+    <
+        typename Foam::innerProduct<typename SFType::value_type, Type>::type,
+        Foam::fvsPatchField,
+        Foam::surfaceMesh
+    >
+>
+Foam::surfaceInterpolationScheme<Type>::dotInterpolate
 (
+    const SFType& Sf,
     const GeometricField<Type, fvPatchField, volMesh>& vf,
     const tmp<surfaceScalarField>& tlambdas
 )
@@ -233,6 +244,9 @@ Foam::surfaceInterpolationScheme<Type>::interpolate
             << endl;
     }
 
+    typedef typename Foam::innerProduct<typename SFType::value_type, Type>::type
+        RetType;
+
     const surfaceScalarField& lambdas = tlambdas();
 
     const Field<Type>& vfi = vf.internalField();
@@ -242,9 +256,9 @@ Foam::surfaceInterpolationScheme<Type>::interpolate
     const labelUList& P = mesh.owner();
     const labelUList& N = mesh.neighbour();
 
-    tmp<GeometricField<Type, fvsPatchField, surfaceMesh>> tsf
+    tmp<GeometricField<RetType, fvsPatchField, surfaceMesh>> tsf
     (
-        new GeometricField<Type, fvsPatchField, surfaceMesh>
+        new GeometricField<RetType, fvsPatchField, surfaceMesh>
         (
             IOobject
             (
@@ -253,16 +267,18 @@ Foam::surfaceInterpolationScheme<Type>::interpolate
                 vf.db()
             ),
             mesh,
-            vf.dimensions()
+            Sf.dimensions()*vf.dimensions()
         )
     );
-    GeometricField<Type, fvsPatchField, surfaceMesh>& sf = tsf.ref();
+    GeometricField<RetType, fvsPatchField, surfaceMesh>& sf = tsf.ref();
 
-    Field<Type>& sfi = sf.internalField();
+    Field<RetType>& sfi = sf.internalField();
+
+    const typename SFType::InternalField& Sfi = Sf.internalField();
 
     for (label fi=0; fi<P.size(); fi++)
     {
-        sfi[fi] = lambda[fi]*(vfi[P[fi]] - vfi[N[fi]]) + vfi[N[fi]];
+        sfi[fi] = Sfi[fi] & (lambda[fi]*(vfi[P[fi]] - vfi[N[fi]]) + vfi[N[fi]]);
     }
 
     // Interpolate across coupled patches using given lambdas
@@ -270,16 +286,21 @@ Foam::surfaceInterpolationScheme<Type>::interpolate
     forAll(lambdas.boundaryField(), pi)
     {
         const fvsPatchScalarField& pLambda = lambdas.boundaryField()[pi];
+        const typename SFType::PatchFieldType& pSf = Sf.boundaryField()[pi];
+        fvsPatchField<RetType>& psf = sf.boundaryField()[pi];
 
         if (vf.boundaryField()[pi].coupled())
         {
-            tsf.ref().boundaryField()[pi] =
-                pLambda*vf.boundaryField()[pi].patchInternalField()
-             + (1.0 - pLambda)*vf.boundaryField()[pi].patchNeighbourField();
+            psf =
+                pSf
+              & (
+                    pLambda*vf.boundaryField()[pi].patchInternalField()
+                  + (1.0 - pLambda)*vf.boundaryField()[pi].patchNeighbourField()
+                );
         }
         else
         {
-            sf.boundaryField()[pi] = vf.boundaryField()[pi];
+            psf = pSf & vf.boundaryField()[pi];
         }
     }
 
@@ -289,6 +310,93 @@ Foam::surfaceInterpolationScheme<Type>::interpolate
 }
 
 
+template<class Type>
+Foam::tmp<Foam::GeometricField<Type, Foam::fvsPatchField, Foam::surfaceMesh>>
+Foam::surfaceInterpolationScheme<Type>::interpolate
+(
+    const GeometricField<Type, fvPatchField, volMesh>& vf,
+    const tmp<surfaceScalarField>& tlambdas
+)
+{
+    return dotInterpolate(geometricOneField(), vf, tlambdas);
+}
+
+
+template<class Type>
+Foam::tmp
+<
+    Foam::GeometricField
+    <
+        typename Foam::innerProduct<Foam::vector, Type>::type,
+        Foam::fvsPatchField,
+        Foam::surfaceMesh
+    >
+>
+Foam::surfaceInterpolationScheme<Type>::dotInterpolate
+(
+    const surfaceVectorField& Sf,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+) const
+{
+    if (surfaceInterpolation::debug)
+    {
+        InfoInFunction
+            << "Interpolating "
+            << vf.type() << " "
+            << vf.name()
+            << " from cells to faces"
+            << endl;
+    }
+
+    tmp
+    <
+        GeometricField
+        <
+            typename Foam::innerProduct<Foam::vector, Type>::type,
+            fvsPatchField,
+            surfaceMesh
+        >
+    > tsf = dotInterpolate(Sf, vf, weights(vf));
+
+    if (corrected())
+    {
+        tsf.ref() += Sf & correction(vf);
+    }
+
+    return tsf;
+}
+
+
+template<class Type>
+Foam::tmp
+<
+    Foam::GeometricField
+    <
+        typename Foam::innerProduct<Foam::vector, Type>::type,
+        Foam::fvsPatchField,
+        Foam::surfaceMesh
+    >
+>
+Foam::surfaceInterpolationScheme<Type>::dotInterpolate
+(
+    const surfaceVectorField& Sf,
+    const tmp<GeometricField<Type, fvPatchField, volMesh>>& tvf
+) const
+{
+    tmp
+    <
+        GeometricField
+        <
+            typename Foam::innerProduct<Foam::vector, Type>::type,
+            fvsPatchField,
+            surfaceMesh
+        >
+    > tSfDotinterpVf = dotInterpolate(Sf, tvf());
+    tvf.clear();
+    return tSfDotinterpVf;
+}
+
+
 template<class Type>
 Foam::tmp<Foam::GeometricField<Type, Foam::fvsPatchField, Foam::surfaceMesh>>
 Foam::surfaceInterpolationScheme<Type>::interpolate
diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationScheme.H b/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationScheme.H
index 7780f9d8bb..da6a64fbff 100644
--- a/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationScheme.H
+++ b/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationScheme.H
@@ -156,6 +156,25 @@ public:
             const tmp<surfaceScalarField>&
         );
 
+        //- Return the face-interpolate of the given cell field
+        //  with the given weighting factors dotted with given field Sf
+        template<class SFType>
+        static tmp
+        <
+            GeometricField
+            <
+                typename innerProduct<typename SFType::value_type, Type>::type,
+                fvsPatchField,
+                surfaceMesh
+            >
+        >
+        dotInterpolate
+        (
+            const SFType& Sf,
+            const GeometricField<Type, fvPatchField, volMesh>& vf,
+            const tmp<surfaceScalarField>& tlambdas
+        );
+
         //- Return the face-interpolate of the given cell field
         //  with the given weighting factors
         static tmp<GeometricField<Type, fvsPatchField, surfaceMesh>>
@@ -165,7 +184,6 @@ public:
             const tmp<surfaceScalarField>&
         );
 
-
         //- Return the interpolation weighting factors for the given field
         virtual tmp<surfaceScalarField> weights
         (
@@ -186,6 +204,41 @@ public:
             return tmp<GeometricField<Type, fvsPatchField, surfaceMesh>>(NULL);
         }
 
+        //- Return the face-interpolate of the given cell field
+        //  with explicit correction dotted with given field Sf
+        virtual
+        tmp
+        <
+            GeometricField
+            <
+                typename innerProduct<vector, Type>::type,
+                fvsPatchField,
+                surfaceMesh
+            >
+        >
+        dotInterpolate
+        (
+            const surfaceVectorField& Sf,
+            const GeometricField<Type, fvPatchField, volMesh>& vf
+        ) const;
+
+        //- Return the face-interpolate of the given tmp cell field
+        //  with explicit correction dotted with given field Sf
+        tmp
+        <
+            GeometricField
+            <
+                typename innerProduct<vector, Type>::type,
+                fvsPatchField,
+                surfaceMesh
+            >
+        >
+        dotInterpolate
+        (
+            const surfaceVectorField& Sf,
+            const tmp<GeometricField<Type, fvPatchField, volMesh>>&
+        ) const;
+
         //- Return the face-interpolate of the given cell field
         //  with explicit correction
         virtual tmp<GeometricField<Type, fvsPatchField, surfaceMesh>>
@@ -201,6 +254,23 @@ public:
 };
 
 
+template<>
+tmp
+<
+    GeometricField
+    <
+        typename innerProduct<vector, scalar>::type,
+        fvsPatchField,
+        surfaceMesh
+    >
+>
+surfaceInterpolationScheme<scalar>::dotInterpolate
+(
+    const surfaceVectorField& Sf,
+    const GeometricField<scalar, fvPatchField, volMesh>&
+) const;
+
+
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 } // End namespace Foam
diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationSchemes.C b/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationSchemes.C
index 7332b7cbde..a6da1c8678 100644
--- a/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationSchemes.C
+++ b/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationSchemes.C
@@ -54,9 +54,49 @@ makeBaseSurfaceInterpolationScheme(sphericalTensor)
 makeBaseSurfaceInterpolationScheme(symmTensor)
 makeBaseSurfaceInterpolationScheme(tensor)
 
-
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 } // End namespace Foam
 
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+template<>
+Foam::tmp
+<
+    Foam::GeometricField
+    <
+        typename Foam::innerProduct<Foam::vector, Foam::scalar>::type,
+        Foam::fvsPatchField,
+        Foam::surfaceMesh
+    >
+>
+Foam::surfaceInterpolationScheme<Foam::scalar>::dotInterpolate
+(
+    const surfaceVectorField& Sf,
+    const GeometricField<scalar, fvPatchField, volMesh>&
+) const
+{
+    NotImplemented;
+
+    return
+        tmp
+        <
+            GeometricField
+            <
+                typename innerProduct<vector, scalar>::type,
+                fvsPatchField,
+                surfaceMesh
+            >
+        >
+        (
+            GeometricField
+            <
+                typename innerProduct<vector, scalar>::type,
+                fvsPatchField,
+                surfaceMesh
+            >::null()
+        );
+}
+
+
 // ************************************************************************* //
-- 
GitLab