Simple reconstruct method using Gauss rather than least-squares

src/finiteVolume/finiteVolume/fvc/fvcSimpleReconstruct.C

@@ -0,0 +1,155 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 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/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "fvcReconstruct.H"
+#include "fvMesh.H"
+#include "zeroGradientFvPatchFields.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace fvc
+{
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+template<class Type>
+tmp
+<
+    GeometricField
+    <
+        typename outerProduct<vector,Type>::type, fvPatchField, volMesh
+    >
+>
+reconstruct
+(
+    const GeometricField<Type, fvsPatchField, surfaceMesh>& ssf
+)
+{
+    typedef typename outerProduct<vector, Type>::type GradType;
+
+    const fvMesh& mesh = ssf.mesh();
+
+    const labelUList& owner = mesh.owner();
+    const labelUList& neighbour = mesh.neighbour();
+
+    const volVectorField& C = mesh.C();
+    const surfaceVectorField& Cf = mesh.Cf();
+
+    tmp<GeometricField<GradType, fvPatchField, volMesh> > treconField
+    (
+        new GeometricField<GradType, fvPatchField, volMesh>
+        (
+            IOobject
+            (
+                "reconstruct("+ssf.name()+')',
+                ssf.instance(),
+                mesh,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            mesh,
+            dimensioned<GradType>
+            (
+                "0",
+                ssf.dimensions()/dimArea,
+                pTraits<GradType>::zero
+            ),
+            zeroGradientFvPatchField<GradType>::typeName
+        )
+    );
+
+    Field<GradType>& rf = treconField();
+
+    forAll(owner, facei)
+    {
+        label own = owner[facei];
+        label nei = neighbour[facei];
+
+        rf[own] += (Cf[facei] - C[own])*ssf[facei];
+        rf[nei] -= (Cf[facei] - C[nei])*ssf[facei];
+    }
+
+    const typename GeometricField<Type, fvsPatchField, surfaceMesh>::
+    GeometricBoundaryField& bsf = ssf.boundaryField();
+
+    forAll(bsf, patchi)
+    {
+        const fvsPatchField<Type>& psf = bsf[patchi];
+
+        const labelUList& pOwner = mesh.boundary()[patchi].faceCells();
+        const vectorField& pCf = Cf.boundaryField()[patchi];
+
+        forAll(pOwner, pFacei)
+        {
+            label own = pOwner[pFacei];
+            rf[own] += (pCf[pFacei] - C[own])*psf[pFacei];
+        }
+    }
+
+    rf /= mesh.V();
+
+    treconField().correctBoundaryConditions();
+
+    return treconField;
+}
+
+
+template<class Type>
+tmp
+<
+    GeometricField
+    <
+        typename outerProduct<vector, Type>::type, fvPatchField, volMesh
+    >
+>
+reconstruct
+(
+    const tmp<GeometricField<Type, fvsPatchField, surfaceMesh> >& tssf
+)
+{
+    typedef typename outerProduct<vector, Type>::type GradType;
+    tmp<GeometricField<GradType, fvPatchField, volMesh> > tvf
+    (
+        fvc::reconstruct(tssf())
+    );
+    tssf.clear();
+    return tvf;
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace fvc
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// ************************************************************************* //