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.

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)
         )
     );
 

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,

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] =
+            psf =
+                pSf
+              & (
                     pLambda*vf.boundaryField()[pi].patchInternalField()
-             + (1.0 - pLambda)*vf.boundaryField()[pi].patchNeighbourField();
+                  + (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

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

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()
+        );
+}
+
+
 // ************************************************************************* //