ENH: Added localEulerDdtScheme time scheme from 1.7.x

src/finiteVolume/Make/files

@@ -287,6 +287,7 @@ $(ddtSchemes)/steadyStateDdtScheme/steadyStateDdtSchemes.C
 $(ddtSchemes)/EulerDdtScheme/EulerDdtSchemes.C
 $(ddtSchemes)/CoEulerDdtScheme/CoEulerDdtSchemes.C
 $(ddtSchemes)/SLTSDdtScheme/SLTSDdtSchemes.C
+$(ddtSchemes)/localEulerDdtScheme/localEulerDdtSchemes.C
 $(ddtSchemes)/backwardDdtScheme/backwardDdtSchemes.C
 $(ddtSchemes)/boundedBackwardDdtScheme/boundedBackwardDdtScheme.C
 $(ddtSchemes)/boundedBackwardDdtScheme/boundedBackwardDdtSchemes.C

src/finiteVolume/finiteVolume/ddtSchemes/localEulerDdtScheme/localEulerDdtScheme.C

+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 1991-2010 OpenCFD Ltd.
+     \\/     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 "localEulerDdtScheme.H"
+#include "surfaceInterpolate.H"
+#include "fvMatrices.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace fv
+{
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+template<class Type>
+const volScalarField& localEulerDdtScheme<Type>::localRDeltaT() const
+{
+    return mesh().objectRegistry::lookupObject<volScalarField>(rDeltaTName_);
+}
+
+
+template<class Type>
+tmp<GeometricField<Type, fvPatchField, volMesh> >
+localEulerDdtScheme<Type>::fvcDdt
+(
+    const dimensioned<Type>& dt
+)
+{
+    const volScalarField& rDeltaT = localRDeltaT();
+
+    IOobject ddtIOobject
+    (
+        "ddt(" + dt.name() + ')',
+        mesh().time().timeName(),
+        mesh()
+    );
+
+    if (mesh().moving())
+    {
+        tmp<GeometricField<Type, fvPatchField, volMesh> > tdtdt
+        (
+            new GeometricField<Type, fvPatchField, volMesh>
+            (
+                ddtIOobject,
+                mesh(),
+                dimensioned<Type>
+                (
+                    "0",
+                    dt.dimensions()/dimTime,
+                    pTraits<Type>::zero
+                )
+            )
+        );
+
+        tdtdt().internalField() =
+            rDeltaT.internalField()*dt.value()*(1.0 - mesh().V0()/mesh().V());
+
+        return tdtdt;
+    }
+    else
+    {
+        return tmp<GeometricField<Type, fvPatchField, volMesh> >
+        (
+            new GeometricField<Type, fvPatchField, volMesh>
+            (
+                ddtIOobject,
+                mesh(),
+                dimensioned<Type>
+                (
+                    "0",
+                    dt.dimensions()/dimTime,
+                    pTraits<Type>::zero
+                ),
+                calculatedFvPatchField<Type>::typeName
+            )
+        );
+    }
+}
+
+
+template<class Type>
+tmp<GeometricField<Type, fvPatchField, volMesh> >
+localEulerDdtScheme<Type>::fvcDdt
+(
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    const volScalarField& rDeltaT = localRDeltaT();
+
+    IOobject ddtIOobject
+    (
+        "ddt(" + vf.name() + ')',
+        mesh().time().timeName(),
+        mesh()
+    );
+
+    if (mesh().moving())
+    {
+        return tmp<GeometricField<Type, fvPatchField, volMesh> >
+        (
+            new GeometricField<Type, fvPatchField, volMesh>
+            (
+                ddtIOobject,
+                mesh(),
+                rDeltaT.dimensions()*vf.dimensions(),
+                rDeltaT.internalField()*
+                (
+                    vf.internalField()
+                  - vf.oldTime().internalField()*mesh().V0()/mesh().V()
+                ),
+                rDeltaT.boundaryField()*
+                (
+                    vf.boundaryField() - vf.oldTime().boundaryField()
+                )
+            )
+        );
+    }
+    else
+    {
+        return tmp<GeometricField<Type, fvPatchField, volMesh> >
+        (
+            new GeometricField<Type, fvPatchField, volMesh>
+            (
+                ddtIOobject,
+                rDeltaT*(vf - vf.oldTime())
+            )
+        );
+    }
+}
+
+
+template<class Type>
+tmp<GeometricField<Type, fvPatchField, volMesh> >
+localEulerDdtScheme<Type>::fvcDdt
+(
+    const dimensionedScalar& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    const volScalarField& rDeltaT = localRDeltaT();
+
+    IOobject ddtIOobject
+    (
+        "ddt(" + rho.name() + ',' + vf.name() + ')',
+        mesh().time().timeName(),
+        mesh()
+    );
+
+    if (mesh().moving())
+    {
+        return tmp<GeometricField<Type, fvPatchField, volMesh> >
+        (
+            new GeometricField<Type, fvPatchField, volMesh>
+            (
+                ddtIOobject,
+                mesh(),
+                rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
+                rDeltaT.internalField()*rho.value()*
+                (
+                    vf.internalField()
+                  - vf.oldTime().internalField()*mesh().V0()/mesh().V()
+                ),
+                rDeltaT.boundaryField()*rho.value()*
+                (
+                    vf.boundaryField() - vf.oldTime().boundaryField()
+                )
+            )
+        );
+    }
+    else
+    {
+        return tmp<GeometricField<Type, fvPatchField, volMesh> >
+        (
+            new GeometricField<Type, fvPatchField, volMesh>
+            (
+                ddtIOobject,
+                rDeltaT*rho*(vf - vf.oldTime())
+            )
+        );
+    }
+}
+
+
+template<class Type>
+tmp<GeometricField<Type, fvPatchField, volMesh> >
+localEulerDdtScheme<Type>::fvcDdt
+(
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    const volScalarField& rDeltaT = localRDeltaT();
+
+    IOobject ddtIOobject
+    (
+        "ddt(" + rho.name() + ',' + vf.name() + ')',
+        mesh().time().timeName(),
+        mesh()
+    );
+
+    if (mesh().moving())
+    {
+        return tmp<GeometricField<Type, fvPatchField, volMesh> >
+        (
+            new GeometricField<Type, fvPatchField, volMesh>
+            (
+                ddtIOobject,
+                mesh(),
+                rDeltaT.dimensions()*rho.dimensions()*vf.dimensions(),
+                rDeltaT.internalField()*
+                (
+                    rho.internalField()*vf.internalField()
+                  - rho.oldTime().internalField()
+                   *vf.oldTime().internalField()*mesh().V0()/mesh().V()
+                ),
+                rDeltaT.boundaryField()*
+                (
+                    rho.boundaryField()*vf.boundaryField()
+                  - rho.oldTime().boundaryField()
+                   *vf.oldTime().boundaryField()
+                )
+            )
+        );
+    }
+    else
+    {
+        return tmp<GeometricField<Type, fvPatchField, volMesh> >
+        (
+            new GeometricField<Type, fvPatchField, volMesh>
+            (
+                ddtIOobject,
+                rDeltaT*(rho*vf - rho.oldTime()*vf.oldTime())
+            )
+        );
+    }
+}
+
+
+template<class Type>
+tmp<fvMatrix<Type> >
+localEulerDdtScheme<Type>::fvmDdt
+(
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    tmp<fvMatrix<Type> > tfvm
+    (
+        new fvMatrix<Type>
+        (
+            vf,
+            vf.dimensions()*dimVol/dimTime
+        )
+    );
+
+    fvMatrix<Type>& fvm = tfvm();
+
+    const scalarField& rDeltaT = localRDeltaT().internalField();
+
+    fvm.diag() = rDeltaT*mesh().V();
+
+    if (mesh().moving())
+    {
+        fvm.source() = rDeltaT*vf.oldTime().internalField()*mesh().V0();
+    }
+    else
+    {
+        fvm.source() = rDeltaT*vf.oldTime().internalField()*mesh().V();
+    }
+
+    return tfvm;
+}
+
+
+template<class Type>
+tmp<fvMatrix<Type> >
+localEulerDdtScheme<Type>::fvmDdt
+(
+    const dimensionedScalar& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    tmp<fvMatrix<Type> > tfvm
+    (
+        new fvMatrix<Type>
+        (
+            vf,
+            rho.dimensions()*vf.dimensions()*dimVol/dimTime
+        )
+    );
+    fvMatrix<Type>& fvm = tfvm();
+
+    const scalarField& rDeltaT = localRDeltaT().internalField();
+
+    fvm.diag() = rDeltaT*rho.value()*mesh().V();
+
+    if (mesh().moving())
+    {
+        fvm.source() = rDeltaT
+            *rho.value()*vf.oldTime().internalField()*mesh().V0();
+    }
+    else
+    {
+        fvm.source() = rDeltaT
+            *rho.value()*vf.oldTime().internalField()*mesh().V();
+    }
+
+    return tfvm;
+}
+
+
+template<class Type>
+tmp<fvMatrix<Type> >
+localEulerDdtScheme<Type>::fvmDdt
+(
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    tmp<fvMatrix<Type> > tfvm
+    (
+        new fvMatrix<Type>
+        (
+            vf,
+            rho.dimensions()*vf.dimensions()*dimVol/dimTime
+        )
+    );
+    fvMatrix<Type>& fvm = tfvm();
+
+    const scalarField& rDeltaT = localRDeltaT().internalField();
+
+    fvm.diag() = rDeltaT*rho.internalField()*mesh().V();
+
+    if (mesh().moving())
+    {
+        fvm.source() = rDeltaT
+            *rho.oldTime().internalField()
+            *vf.oldTime().internalField()*mesh().V0();
+    }
+    else
+    {
+        fvm.source() = rDeltaT
+            *rho.oldTime().internalField()
+            *vf.oldTime().internalField()*mesh().V();
+    }
+
+    return tfvm;
+}
+
+
+template<class Type>
+tmp<typename localEulerDdtScheme<Type>::fluxFieldType>
+localEulerDdtScheme<Type>::fvcDdtPhiCorr
+(
+    const volScalarField& rA,
+    const GeometricField<Type, fvPatchField, volMesh>& U,
+    const fluxFieldType& phi
+)
+{
+    IOobject ddtIOobject
+    (
+        "ddtPhiCorr(" + rA.name() + ',' + U.name() + ',' + phi.name() + ')',
+        mesh().time().timeName(),
+        mesh()
+    );
+
+    if (mesh().moving())
+    {
+        return tmp<fluxFieldType>
+        (
+            new fluxFieldType
+            (
+                ddtIOobject,
+                mesh(),
+                dimensioned<typename flux<Type>::type>
+                (
+                    "0",
+                    rA.dimensions()*phi.dimensions()/dimTime,
+                    pTraits<typename flux<Type>::type>::zero
+                )
+            )
+        );
+    }
+    else
+    {
+        const volScalarField& rDeltaT = localRDeltaT();
+
+        return tmp<fluxFieldType>
+        (
+            new fluxFieldType
+            (
+                ddtIOobject,
+                fvcDdtPhiCoeff(U.oldTime(), phi.oldTime())*
+                (
+                    fvc::interpolate(rDeltaT*rA)*phi.oldTime()
+                  - (fvc::interpolate(rDeltaT*rA*U.oldTime()) & mesh().Sf())
+                )
+            )
+        );
+    }
+}
+
+
+template<class Type>
+tmp<typename localEulerDdtScheme<Type>::fluxFieldType>
+localEulerDdtScheme<Type>::fvcDdtPhiCorr
+(
+    const volScalarField& rA,
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& U,
+    const fluxFieldType& phi
+)
+{
+    IOobject ddtIOobject
+    (
+        "ddtPhiCorr("
+      + rA.name() + ',' + rho.name() + ',' + U.name() + ',' + phi.name() + ')',
+        mesh().time().timeName(),
+        mesh()
+    );
+
+    if (mesh().moving())
+    {
+        return tmp<fluxFieldType>
+        (
+            new fluxFieldType
+            (
+                ddtIOobject,
+                mesh(),
+                dimensioned<typename flux<Type>::type>
+                (
+                    "0",
+                    rA.dimensions()*rho.dimensions()*phi.dimensions()/dimTime,
+                    pTraits<typename flux<Type>::type>::zero
+                )
+            )
+        );
+    }
+    else
+    {
+        const volScalarField& rDeltaT = localRDeltaT();
+
+        if
+        (
+            U.dimensions() == dimVelocity
+         && phi.dimensions() == dimVelocity*dimArea
+        )
+        {
+            return tmp<fluxFieldType>
+            (
+                new fluxFieldType
+                (
+                    ddtIOobject,
+                    fvcDdtPhiCoeff(U.oldTime(), phi.oldTime())
+                   *(
+                        fvc::interpolate(rDeltaT*rA*rho.oldTime())*phi.oldTime()
+                      - (fvc::interpolate(rDeltaT*rA*rho.oldTime()*U.oldTime())
+                      & mesh().Sf())
+                    )
+                )
+            );
+        }
+        else if
+        (
+            U.dimensions() == dimVelocity
+         && phi.dimensions() == dimDensity*dimVelocity*dimArea
+        )
+        {
+            return tmp<fluxFieldType>
+            (
+                new fluxFieldType
+                (
+                    ddtIOobject,
+                    fvcDdtPhiCoeff
+                    (
+                        U.oldTime(),
+                        phi.oldTime()/fvc::interpolate(rho.oldTime())
+                    )
+                   *(
+                        fvc::interpolate(rDeltaT*rA*rho.oldTime())
+                       *phi.oldTime()/fvc::interpolate(rho.oldTime())
+                      - (
+                            fvc::interpolate
+                            (
+                                rDeltaT*rA*rho.oldTime()*U.oldTime()
+                            ) & mesh().Sf()
+                        )
+                    )
+                )
+            );
+        }
+        else if
+        (
+            U.dimensions() == dimDensity*dimVelocity
+         && phi.dimensions() == dimDensity*dimVelocity*dimArea
+        )
+        {
+            return tmp<fluxFieldType>
+            (
+                new fluxFieldType
+                (
+                    ddtIOobject,
+                    fvcDdtPhiCoeff(rho.oldTime(), U.oldTime(), phi.oldTime())
+                   *(
+                        fvc::interpolate(rDeltaT*rA)*phi.oldTime()
+                      - (
+                            fvc::interpolate(rDeltaT*rA*U.oldTime())&mesh().Sf()
+                        )
+                    )
+                )
+            );
+        }
+        else
+        {
+            FatalErrorIn
+            (
+                "localEulerDdtScheme<Type>::fvcDdtPhiCorr"
+            )   << "dimensions of phi are not correct"
+                << abort(FatalError);
+
+            return fluxFieldType::null();
+        }
+    }
+}
+
+
+template<class Type>
+tmp<surfaceScalarField> localEulerDdtScheme<Type>::meshPhi
+(
+    const GeometricField<Type, fvPatchField, volMesh>&
+)
+{
+    return tmp<surfaceScalarField>
+    (
+        new surfaceScalarField
+        (
+            IOobject
+            (
+                "meshPhi",
+                mesh().time().timeName(),
+                mesh()
+            ),
+            mesh(),
+            dimensionedScalar("0", dimVolume/dimTime, 0.0)
+        )
+    );
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //