Merge branch 'feature-fa-interpolation-scheme' into 'develop'

ENH: leastSquaresEdgeInterpolation: new edge interpolation method

See merge request Development/openfoam!575

src/finiteArea/Make/files

@@ -102,6 +102,7 @@ $(schemes)/linearUpwind/linearUpwindEdgeInterpolationMake.C
 $(schemes)/Gamma/GammaEdgeInterpolationMake.C
 $(schemes)/blended/blendedEdgeInterpolationMake.C
 $(schemes)/skewCorrected/skewCorrectedEdgeInterpolationMake.C
+$(schemes)/leastSquares/leastSquaresEdgeInterpolationMake.C
 
 finiteArea/fa/fa.C
 

src/finiteArea/interpolation/edgeInterpolation/schemes/leastSquares/leastSquaresEdgeInterpolation.H

@@ -0,0 +1,229 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2022 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::leastSquaresEdgeInterpolation
+
+Description
+    Least-squares edge interpolation scheme class, from
+    face centers to points then from points to edges.
+
+    References:
+    \verbatim
+        Governing equations (tag:P):
+            Pesci, C. (2019).
+            Computational analysis of fluid interfaces
+            influenced by soluble surfactant.
+            Darmstadt, Technische Universität. PhD thesis.
+            URI:https://tuprints.ulb.tu-darmstadt.de/id/eprint/9303
+    \endverbatim
+
+SourceFiles
+    leastSquaresEdgeInterpolationMake.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef Foam_leastSquaresEdgeInterpolation_H
+#define Foam_leastSquaresEdgeInterpolation_H
+
+#include "edgeInterpolationScheme.H"
+#include "areaFields.H"
+#include "edgeFields.H"
+#include "leastSquaresFaGrad.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                Class leastSquaresEdgeInterpolation Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class Type>
+class leastSquaresEdgeInterpolation
+:
+    virtual public edgeInterpolationScheme<Type>
+{
+    // Private Member Functions
+
+        //- No copy assignment
+        void operator=(const leastSquaresEdgeInterpolation&) = delete;
+
+
+public:
+
+    //- Runtime type information
+    TypeName("leastSquares");
+
+
+    // Constructors
+
+        //- Construct from mesh
+        leastSquaresEdgeInterpolation(const faMesh& mesh)
+        :
+            edgeInterpolationScheme<Type>(mesh)
+        {}
+
+        //- Construct from Istream
+        leastSquaresEdgeInterpolation(const faMesh& mesh, Istream&)
+        :
+            edgeInterpolationScheme<Type>(mesh)
+        {}
+
+        //- Construct from faceFlux and Istream
+        leastSquaresEdgeInterpolation
+        (
+            const faMesh& mesh,
+            const edgeScalarField&,
+            Istream&
+        )
+        :
+            edgeInterpolationScheme<Type>(mesh)
+        {}
+
+
+    // Member Functions
+
+        //- Return the edge-interpolate of the given face field
+        tmp<GeometricField<Type, faePatchField, edgeMesh>>
+        interpolate
+        (
+            const GeometricField<Type, faPatchField, areaMesh>&
+        ) const;
+
+        //- Return the interpolation weighting factors
+        tmp<edgeScalarField> weights
+        (
+            const GeometricField<Type, faPatchField, areaMesh>&
+        ) const
+        {
+            return this->mesh().edgeInterpolation::weights();
+        }
+};
+
+
+template<class Type>
+tmp<GeometricField<Type, faePatchField, edgeMesh>>
+leastSquaresEdgeInterpolation<Type>::interpolate
+(
+    const GeometricField<Type, faPatchField, areaMesh>& af
+) const
+{
+    typedef typename outerProduct<vector, Type>::type GradType;
+    typedef GeometricField<GradType, faPatchField, areaMesh> GradFieldType;
+    typedef GeometricField<Type, faePatchField, edgeMesh> EdgeFieldType;
+    typedef Field<Type> FieldType;
+
+    const faMesh& mesh = af.mesh();
+
+    const labelList& meshPts = mesh.patch().meshPoints();
+    const labelListList& pointFaceAddr = mesh.patch().pointFaces();
+    const areaVectorField& C = mesh.areaCentres();
+
+    tmp<GradFieldType> tgradAf= fa::leastSquaresFaGrad<Type>(mesh).grad(af);
+    const GradFieldType& gradAf = tgradAf.cref();
+
+    // Field at surface points, Pi (P:p. 29-30)
+    auto tPi = tmp<FieldType>::New(meshPts.size(), Zero);
+    FieldType& Pi = tPi.ref();
+
+    forAll(meshPts, i)
+    {
+        const label nFaces = pointFaceAddr[i].size();
+
+        FieldType Pij(nFaces);
+        vectorField dPC(nFaces);
+
+        for (label facei = 0; facei < nFaces; ++facei)
+        {
+            const label j = pointFaceAddr[i][facei];
+
+            // Euclidean distance between a point and face centres (P:Fig. 3.3)
+            dPC[facei] = mesh.points()[i] - C[j];
+
+            // (P:Eq. 3.14)
+            Pij[facei] = af[j] + (dPC[facei] & gradAf[j]);
+        }
+
+        const scalarField magdPC(max(mag(dPC), SMALL));
+
+        // (P:Eq. 3.15)
+        Pi[i] = sum(Pij/magdPC)/sum(scalar(1)/magdPC);
+    }
+
+    tgradAf.clear();
+
+
+    auto tinterp = tmp<EdgeFieldType>::New
+    (
+        IOobject
+        (
+            "interpolate(" + af.name() + ')',
+            af.instance(),
+            af.db()
+        ),
+        af.mesh(),
+        af.dimensions()
+    );
+    EdgeFieldType& interp = tinterp.ref();
+    FieldType& interpEdges = interp.primitiveFieldRef();
+
+
+    const edgeList& faEdges = mesh.edges();
+    const label nInternalEdges = mesh.nInternalEdges();
+
+    // Internal field
+    for (label edgei = 0; edgei < nInternalEdges; ++edgei)
+    {
+        // (P:Eq. 3.16)
+        interpEdges[edgei] =
+            0.5*
+            (
+                Pi[faEdges[edgei].start()]
+              + Pi[faEdges[edgei].end()]
+            );
+    }
+
+    // Boundary field from boundary condition [fixedValue]
+    // Coupled boundary fields are not handled
+    forAll(interp.boundaryField(), patchi)
+    {
+        interp.boundaryFieldRef()[patchi] = af.boundaryField()[patchi];
+    }
+
+    return tinterp;
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/finiteArea/interpolation/edgeInterpolation/schemes/leastSquares/leastSquaresEdgeInterpolationMake.C

@@ -0,0 +1,38 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2022 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 "leastSquaresEdgeInterpolation.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    makeEdgeInterpolationTypeScheme(leastSquaresEdgeInterpolation, scalar)
+    makeEdgeInterpolationTypeScheme(leastSquaresEdgeInterpolation, vector)
+}
+
+// ************************************************************************* //