From b7fcf23f690bd53aef7890cccf42dececf0b5f31 Mon Sep 17 00:00:00 2001
From: Henry <Henry>
Date: Tue, 29 Apr 2014 14:54:18 +0100
Subject: [PATCH] finiteVolume: added support for ddt(alpha, rho, psi)

---
 .../gaussConvectionScheme.H                   |   2 +-
 .../CoEulerDdtScheme/CoEulerDdtScheme.C       | 113 +++++++-
 .../CoEulerDdtScheme/CoEulerDdtScheme.H       |  16 +-
 .../CrankNicolsonDdtScheme.C                  | 243 +++++++++++++++++-
 .../CrankNicolsonDdtScheme.H                  |  16 +-
 .../EulerDdtScheme/EulerDdtScheme.C           | 113 +++++++-
 .../EulerDdtScheme/EulerDdtScheme.H           |  16 +-
 .../ddtSchemes/SLTSDdtScheme/SLTSDdtScheme.C  | 113 +++++++-
 .../ddtSchemes/SLTSDdtScheme/SLTSDdtScheme.H  |  16 +-
 .../backwardDdtScheme/backwardDdtScheme.C     | 163 +++++++++++-
 .../backwardDdtScheme/backwardDdtScheme.H     |  16 +-
 .../boundedDdtScheme/boundedDdtScheme.C       |  28 +-
 .../boundedDdtScheme/boundedDdtScheme.H       |  16 +-
 .../ddtSchemes/ddtScheme/ddtScheme.C          |  42 ++-
 .../ddtSchemes/ddtScheme/ddtScheme.H          |  15 +-
 .../localEulerDdtScheme/localEulerDdtScheme.C | 113 +++++++-
 .../localEulerDdtScheme/localEulerDdtScheme.H |  16 +-
 .../steadyStateDdtScheme.C                    |  55 +++-
 .../steadyStateDdtScheme.H                    |  16 +-
 src/finiteVolume/finiteVolume/fvc/fvcDdt.C    |  25 +-
 src/finiteVolume/finiteVolume/fvc/fvcDdt.H    |  10 +-
 src/finiteVolume/finiteVolume/fvm/fvmDdt.C    |  25 +-
 src/finiteVolume/finiteVolume/fvm/fvmDdt.H    |  10 +-
 23 files changed, 1175 insertions(+), 23 deletions(-)

diff --git a/src/finiteVolume/finiteVolume/convectionSchemes/gaussConvectionScheme/gaussConvectionScheme.H b/src/finiteVolume/finiteVolume/convectionSchemes/gaussConvectionScheme/gaussConvectionScheme.H
index 37761441281..158ff21a13b 100644
--- a/src/finiteVolume/finiteVolume/convectionSchemes/gaussConvectionScheme/gaussConvectionScheme.H
+++ b/src/finiteVolume/finiteVolume/convectionSchemes/gaussConvectionScheme/gaussConvectionScheme.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/CoEulerDdtScheme/CoEulerDdtScheme.C b/src/finiteVolume/finiteVolume/ddtSchemes/CoEulerDdtScheme/CoEulerDdtScheme.C
index ec0bcdca259..8607d49fc1e 100644
--- a/src/finiteVolume/finiteVolume/ddtSchemes/CoEulerDdtScheme/CoEulerDdtScheme.C
+++ b/src/finiteVolume/finiteVolume/ddtSchemes/CoEulerDdtScheme/CoEulerDdtScheme.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -368,6 +368,75 @@ CoEulerDdtScheme<Type>::fvcDdt
 }
 
 
+template<class Type>
+tmp<GeometricField<Type, fvPatchField, volMesh> >
+CoEulerDdtScheme<Type>::fvcDdt
+(
+    const volScalarField& alpha,
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    const volScalarField rDeltaT(CorDeltaT());
+
+    IOobject ddtIOobject
+    (
+        "ddt("+alpha.name()+','+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()
+               *alpha.dimensions()*rho.dimensions()*vf.dimensions(),
+                rDeltaT.internalField()*
+                (
+                    alpha.internalField()
+                   *rho.internalField()
+                   *vf.internalField()
+
+                  - alpha.oldTime().internalField()
+                   *rho.oldTime().internalField()
+                   *vf.oldTime().internalField()*mesh().Vsc0()/mesh().Vsc()
+                ),
+                rDeltaT.boundaryField()*
+                (
+                    alpha.boundaryField()
+                   *rho.boundaryField()
+                   *vf.boundaryField()
+
+                  - alpha.oldTime().boundaryField()
+                   *rho.oldTime().boundaryField()
+                   *vf.oldTime().boundaryField()
+                )
+            )
+        );
+    }
+    else
+    {
+        return tmp<GeometricField<Type, fvPatchField, volMesh> >
+        (
+            new GeometricField<Type, fvPatchField, volMesh>
+            (
+                ddtIOobject,
+                rDeltaT
+               *(
+                   alpha*rho*vf
+                 - alpha.oldTime()*rho.oldTime()*vf.oldTime()
+                )
+            )
+        );
+    }
+}
+
+
 template<class Type>
 tmp<fvMatrix<Type> >
 CoEulerDdtScheme<Type>::fvmDdt
@@ -479,6 +548,48 @@ CoEulerDdtScheme<Type>::fvmDdt
 }
 
 
+template<class Type>
+tmp<fvMatrix<Type> >
+CoEulerDdtScheme<Type>::fvmDdt
+(
+    const volScalarField& alpha,
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    tmp<fvMatrix<Type> > tfvm
+    (
+        new fvMatrix<Type>
+        (
+            vf,
+            alpha.dimensions()*rho.dimensions()*vf.dimensions()*dimVol/dimTime
+        )
+    );
+    fvMatrix<Type>& fvm = tfvm();
+
+    scalarField rDeltaT(CorDeltaT()().internalField());
+
+    fvm.diag() = rDeltaT*alpha.internalField()*rho.internalField()*mesh().Vsc();
+
+    if (mesh().moving())
+    {
+        fvm.source() = rDeltaT
+            *alpha.oldTime().internalField()
+            *rho.oldTime().internalField()
+            *vf.oldTime().internalField()*mesh().Vsc0();
+    }
+    else
+    {
+        fvm.source() = rDeltaT
+            *alpha.oldTime().internalField()
+            *rho.oldTime().internalField()
+            *vf.oldTime().internalField()*mesh().Vsc();
+    }
+
+    return tfvm;
+}
+
+
 template<class Type>
 tmp<typename CoEulerDdtScheme<Type>::fluxFieldType>
 CoEulerDdtScheme<Type>::fvcDdtUfCorr
diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/CoEulerDdtScheme/CoEulerDdtScheme.H b/src/finiteVolume/finiteVolume/ddtSchemes/CoEulerDdtScheme/CoEulerDdtScheme.H
index 204a7935871..e743a0cb403 100644
--- a/src/finiteVolume/finiteVolume/ddtSchemes/CoEulerDdtScheme/CoEulerDdtScheme.H
+++ b/src/finiteVolume/finiteVolume/ddtSchemes/CoEulerDdtScheme/CoEulerDdtScheme.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -139,6 +139,13 @@ public:
             const GeometricField<Type, fvPatchField, volMesh>&
         );
 
+        tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDdt
+        (
+            const volScalarField& alpha,
+            const volScalarField& rho,
+            const GeometricField<Type, fvPatchField, volMesh>& psi
+        );
+
         tmp<fvMatrix<Type> > fvmDdt
         (
             const GeometricField<Type, fvPatchField, volMesh>&
@@ -156,6 +163,13 @@ public:
             const GeometricField<Type, fvPatchField, volMesh>&
         );
 
+        tmp<fvMatrix<Type> > fvmDdt
+        (
+            const volScalarField& alpha,
+            const volScalarField& rho,
+            const GeometricField<Type, fvPatchField, volMesh>& psi
+        );
+
         typedef typename ddtScheme<Type>::fluxFieldType fluxFieldType;
 
         tmp<fluxFieldType> fvcDdtUfCorr
diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/CrankNicolsonDdtScheme/CrankNicolsonDdtScheme.C b/src/finiteVolume/finiteVolume/ddtSchemes/CrankNicolsonDdtScheme/CrankNicolsonDdtScheme.C
index 739a5488555..c9b095fac0d 100644
--- a/src/finiteVolume/finiteVolume/ddtSchemes/CrankNicolsonDdtScheme/CrankNicolsonDdtScheme.C
+++ b/src/finiteVolume/finiteVolume/ddtSchemes/CrankNicolsonDdtScheme/CrankNicolsonDdtScheme.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -620,6 +620,136 @@ CrankNicolsonDdtScheme<Type>::fvcDdt
 }
 
 
+template<class Type>
+tmp<GeometricField<Type, fvPatchField, volMesh> >
+CrankNicolsonDdtScheme<Type>::fvcDdt
+(
+    const volScalarField& alpha,
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    DDt0Field<GeometricField<Type, fvPatchField, volMesh> >& ddt0 =
+        ddt0_<GeometricField<Type, fvPatchField, volMesh> >
+        (
+            "ddt0(" + alpha.name() + ',' + rho.name() + ',' + vf.name() + ')',
+            alpha.dimensions()*rho.dimensions()*vf.dimensions()
+        );
+
+    IOobject ddtIOobject
+    (
+        "ddt(" + alpha.name() + ',' + rho.name() + ',' + vf.name() + ')',
+        mesh().time().timeName(),
+        mesh()
+    );
+
+    dimensionedScalar rDtCoef = rDtCoef_(ddt0);
+
+    if (mesh().moving())
+    {
+        if (evaluate(ddt0))
+        {
+            scalar rDtCoef0 = rDtCoef0_(ddt0).value();
+
+            ddt0.internalField() =
+            (
+                rDtCoef0*
+                (
+                    mesh().V0()
+                   *alpha.oldTime().internalField()
+                   *rho.oldTime().internalField()
+                   *vf.oldTime().internalField()
+
+                  - mesh().V00()
+                   *alpha.oldTime().oldTime().internalField()
+                   *rho.oldTime().oldTime().internalField()
+                   *vf.oldTime().oldTime().internalField()
+                ) - mesh().V00()*offCentre_(ddt0.internalField())
+            )/mesh().V0();
+
+            ddt0.boundaryField() =
+            (
+                rDtCoef0*
+                (
+                    alpha.oldTime().boundaryField()
+                   *rho.oldTime().boundaryField()
+                   *vf.oldTime().boundaryField()
+
+                  - alpha.oldTime().oldTime().boundaryField()
+                   *rho.oldTime().oldTime().boundaryField()
+                   *vf.oldTime().oldTime().boundaryField()
+                ) - offCentre_(ff(ddt0.boundaryField()))
+            );
+        }
+
+        return tmp<GeometricField<Type, fvPatchField, volMesh> >
+        (
+            new GeometricField<Type, fvPatchField, volMesh>
+            (
+                ddtIOobject,
+                mesh(),
+                rDtCoef.dimensions()
+               *alpha.dimensions()*rho.dimensions()*vf.dimensions(),
+                (
+                    rDtCoef.value()*
+                    (
+                        mesh().V()
+                       *alpha.internalField()
+                       *rho.internalField()
+                       *vf.internalField()
+
+                      - mesh().V0()
+                       *alpha.oldTime().internalField()
+                       *rho.oldTime().internalField()
+                       *vf.oldTime().internalField()
+                    ) - mesh().V00()*offCentre_(ddt0.internalField())
+                )/mesh().V(),
+                rDtCoef.value()*
+                (
+                    alpha.boundaryField()
+                   *rho.boundaryField()
+                   *vf.boundaryField()
+
+                  - alpha.oldTime().boundaryField()
+                   *rho.oldTime().boundaryField()
+                   *vf.oldTime().boundaryField()
+                ) - offCentre_(ff(ddt0.boundaryField()))
+            )
+        );
+    }
+    else
+    {
+        if (evaluate(ddt0))
+        {
+            ddt0 = rDtCoef0_(ddt0)*
+            (
+                alpha.oldTime()
+               *rho.oldTime()
+               *vf.oldTime()
+
+              - alpha.oldTime().oldTime()
+               *rho.oldTime().oldTime()
+               *vf.oldTime().oldTime()
+            ) - offCentre_(ddt0());
+        }
+
+        return tmp<GeometricField<Type, fvPatchField, volMesh> >
+        (
+            new GeometricField<Type, fvPatchField, volMesh>
+            (
+                ddtIOobject,
+                rDtCoef
+               *(
+                   alpha*rho*vf
+                 - alpha.oldTime()*rho.oldTime()*vf.oldTime()
+                )
+              - offCentre_(ddt0())
+            )
+        );
+    }
+}
+
+
 template<class Type>
 tmp<fvMatrix<Type> >
 CrankNicolsonDdtScheme<Type>::fvmDdt
@@ -874,6 +1004,117 @@ CrankNicolsonDdtScheme<Type>::fvmDdt
 }
 
 
+template<class Type>
+tmp<fvMatrix<Type> >
+CrankNicolsonDdtScheme<Type>::fvmDdt
+(
+    const volScalarField& alpha,
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    DDt0Field<GeometricField<Type, fvPatchField, volMesh> >& ddt0 =
+        ddt0_<GeometricField<Type, fvPatchField, volMesh> >
+        (
+            "ddt0(" + alpha.name() + ',' + rho.name() + ',' + vf.name() + ')',
+            alpha.dimensions()*rho.dimensions()*vf.dimensions()
+        );
+
+    tmp<fvMatrix<Type> > tfvm
+    (
+        new fvMatrix<Type>
+        (
+            vf,
+            alpha.dimensions()*rho.dimensions()*vf.dimensions()*dimVol/dimTime
+        )
+    );
+    fvMatrix<Type>& fvm = tfvm();
+
+    scalar rDtCoef = rDtCoef_(ddt0).value();
+    fvm.diag() = rDtCoef*alpha.internalField()*rho.internalField()*mesh().V();
+
+    vf.oldTime().oldTime();
+    alpha.oldTime().oldTime();
+    rho.oldTime().oldTime();
+
+    if (mesh().moving())
+    {
+        if (evaluate(ddt0))
+        {
+            scalar rDtCoef0 = rDtCoef0_(ddt0).value();
+
+            ddt0.internalField() =
+            (
+                rDtCoef0*
+                (
+                    mesh().V0()
+                   *alpha.oldTime().internalField()
+                   *rho.oldTime().internalField()
+                   *vf.oldTime().internalField()
+
+                  - mesh().V00()
+                   *alpha.oldTime().oldTime().internalField()
+                   *rho.oldTime().oldTime().internalField()
+                   *vf.oldTime().oldTime().internalField()
+                )
+              - mesh().V00()*offCentre_(ddt0.internalField())
+            )/mesh().V0();
+
+            ddt0.boundaryField() =
+            (
+                rDtCoef0*
+                (
+                    alpha.oldTime().boundaryField()
+                   *rho.oldTime().boundaryField()
+                   *vf.oldTime().boundaryField()
+
+                  - alpha.oldTime().oldTime().boundaryField()
+                   *rho.oldTime().oldTime().boundaryField()
+                   *vf.oldTime().oldTime().boundaryField()
+                )
+              - offCentre_(ff(ddt0.boundaryField()))
+            );
+        }
+
+        fvm.source() =
+        (
+            rDtCoef
+           *alpha.internalField()
+           *rho.internalField()
+           *vf.oldTime().internalField()
+          + offCentre_(ddt0.internalField())
+        )*mesh().V0();
+    }
+    else
+    {
+        if (evaluate(ddt0))
+        {
+            ddt0 = rDtCoef0_(ddt0)*
+            (
+                alpha.oldTime()
+               *rho.oldTime()
+               *vf.oldTime()
+
+              - alpha.oldTime().oldTime()
+               *rho.oldTime().oldTime()
+               *vf.oldTime().oldTime()
+            ) - offCentre_(ddt0());
+        }
+
+        fvm.source() =
+        (
+            rDtCoef
+           *alpha.oldTime().internalField()
+           *rho.oldTime().internalField()
+           *vf.oldTime().internalField()
+          + offCentre_(ddt0.internalField())
+        )*mesh().V();
+    }
+
+    return tfvm;
+}
+
+
 template<class Type>
 tmp<typename CrankNicolsonDdtScheme<Type>::fluxFieldType>
 CrankNicolsonDdtScheme<Type>::fvcDdtUfCorr
diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/CrankNicolsonDdtScheme/CrankNicolsonDdtScheme.H b/src/finiteVolume/finiteVolume/ddtSchemes/CrankNicolsonDdtScheme/CrankNicolsonDdtScheme.H
index 94c34b8ea82..e94c8c2c0bf 100644
--- a/src/finiteVolume/finiteVolume/ddtSchemes/CrankNicolsonDdtScheme/CrankNicolsonDdtScheme.H
+++ b/src/finiteVolume/finiteVolume/ddtSchemes/CrankNicolsonDdtScheme/CrankNicolsonDdtScheme.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -210,6 +210,13 @@ public:
             const GeometricField<Type, fvPatchField, volMesh>&
         );
 
+        tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDdt
+        (
+            const volScalarField& alpha,
+            const volScalarField& rho,
+            const GeometricField<Type, fvPatchField, volMesh>& psi
+        );
+
         tmp<fvMatrix<Type> > fvmDdt
         (
             const GeometricField<Type, fvPatchField, volMesh>&
@@ -227,6 +234,13 @@ public:
             const GeometricField<Type, fvPatchField, volMesh>&
         );
 
+        tmp<fvMatrix<Type> > fvmDdt
+        (
+            const volScalarField& alpha,
+            const volScalarField& rho,
+            const GeometricField<Type, fvPatchField, volMesh>& psi
+        );
+
         typedef typename ddtScheme<Type>::fluxFieldType fluxFieldType;
 
         tmp<fluxFieldType> fvcDdtUfCorr
diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/EulerDdtScheme/EulerDdtScheme.C b/src/finiteVolume/finiteVolume/ddtSchemes/EulerDdtScheme/EulerDdtScheme.C
index 957524992b5..a9ba06b6e29 100644
--- a/src/finiteVolume/finiteVolume/ddtSchemes/EulerDdtScheme/EulerDdtScheme.C
+++ b/src/finiteVolume/finiteVolume/ddtSchemes/EulerDdtScheme/EulerDdtScheme.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -257,6 +257,75 @@ EulerDdtScheme<Type>::fvcDdt
 }
 
 
+template<class Type>
+tmp<GeometricField<Type, fvPatchField, volMesh> >
+EulerDdtScheme<Type>::fvcDdt
+(
+    const volScalarField& alpha,
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
+
+    IOobject ddtIOobject
+    (
+        "ddt("+alpha.name()+','+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()
+               *alpha.dimensions()*rho.dimensions()*vf.dimensions(),
+                rDeltaT.value()*
+                (
+                    alpha.internalField()
+                   *rho.internalField()
+                   *vf.internalField()
+
+                  - alpha.oldTime().internalField()
+                   *rho.oldTime().internalField()
+                   *vf.oldTime().internalField()*mesh().Vsc0()/mesh().Vsc()
+                ),
+                rDeltaT.value()*
+                (
+                    alpha.boundaryField()
+                   *rho.boundaryField()
+                   *vf.boundaryField()
+
+                  - alpha.oldTime().boundaryField()
+                   *rho.oldTime().boundaryField()
+                   *vf.oldTime().boundaryField()
+                )
+            )
+        );
+    }
+    else
+    {
+        return tmp<GeometricField<Type, fvPatchField, volMesh> >
+        (
+            new GeometricField<Type, fvPatchField, volMesh>
+            (
+                ddtIOobject,
+                rDeltaT
+               *(
+                   alpha*rho*vf
+                 - alpha.oldTime()*rho.oldTime()*vf.oldTime()
+                )
+            )
+        );
+    }
+}
+
+
 template<class Type>
 tmp<fvMatrix<Type> >
 EulerDdtScheme<Type>::fvmDdt
@@ -368,6 +437,48 @@ EulerDdtScheme<Type>::fvmDdt
 }
 
 
+template<class Type>
+tmp<fvMatrix<Type> >
+EulerDdtScheme<Type>::fvmDdt
+(
+    const volScalarField& alpha,
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    tmp<fvMatrix<Type> > tfvm
+    (
+        new fvMatrix<Type>
+        (
+            vf,
+            alpha.dimensions()*rho.dimensions()*vf.dimensions()*dimVol/dimTime
+        )
+    );
+    fvMatrix<Type>& fvm = tfvm();
+
+    scalar rDeltaT = 1.0/mesh().time().deltaTValue();
+
+    fvm.diag() = rDeltaT*alpha.internalField()*rho.internalField()*mesh().Vsc();
+
+    if (mesh().moving())
+    {
+        fvm.source() = rDeltaT
+            *alpha.oldTime().internalField()
+            *rho.oldTime().internalField()
+            *vf.oldTime().internalField()*mesh().Vsc0();
+    }
+    else
+    {
+        fvm.source() = rDeltaT
+            *alpha.oldTime().internalField()
+            *rho.oldTime().internalField()
+            *vf.oldTime().internalField()*mesh().Vsc();
+    }
+
+    return tfvm;
+}
+
+
 template<class Type>
 tmp<typename EulerDdtScheme<Type>::fluxFieldType>
 EulerDdtScheme<Type>::fvcDdtUfCorr
diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/EulerDdtScheme/EulerDdtScheme.H b/src/finiteVolume/finiteVolume/ddtSchemes/EulerDdtScheme/EulerDdtScheme.H
index 1e2b840a136..2818fa123f9 100644
--- a/src/finiteVolume/finiteVolume/ddtSchemes/EulerDdtScheme/EulerDdtScheme.H
+++ b/src/finiteVolume/finiteVolume/ddtSchemes/EulerDdtScheme/EulerDdtScheme.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -117,6 +117,13 @@ public:
             const GeometricField<Type, fvPatchField, volMesh>&
         );
 
+        tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDdt
+        (
+            const volScalarField& alpha,
+            const volScalarField& rho,
+            const GeometricField<Type, fvPatchField, volMesh>& psi
+        );
+
         tmp<fvMatrix<Type> > fvmDdt
         (
             const GeometricField<Type, fvPatchField, volMesh>&
@@ -134,6 +141,13 @@ public:
             const GeometricField<Type, fvPatchField, volMesh>&
         );
 
+        tmp<fvMatrix<Type> > fvmDdt
+        (
+            const volScalarField& alpha,
+            const volScalarField& rho,
+            const GeometricField<Type, fvPatchField, volMesh>& psi
+        );
+
         typedef typename ddtScheme<Type>::fluxFieldType fluxFieldType;
 
         tmp<fluxFieldType> fvcDdtUfCorr
diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/SLTSDdtScheme/SLTSDdtScheme.C b/src/finiteVolume/finiteVolume/ddtSchemes/SLTSDdtScheme/SLTSDdtScheme.C
index 95a786c5acc..955bba08071 100644
--- a/src/finiteVolume/finiteVolume/ddtSchemes/SLTSDdtScheme/SLTSDdtScheme.C
+++ b/src/finiteVolume/finiteVolume/ddtSchemes/SLTSDdtScheme/SLTSDdtScheme.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -367,6 +367,75 @@ SLTSDdtScheme<Type>::fvcDdt
 }
 
 
+template<class Type>
+tmp<GeometricField<Type, fvPatchField, volMesh> >
+SLTSDdtScheme<Type>::fvcDdt
+(
+    const volScalarField& alpha,
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    const volScalarField rDeltaT(SLrDeltaT());
+
+    IOobject ddtIOobject
+    (
+        "ddt("+alpha.name()+','+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()
+               *alpha.dimensions()*rho.dimensions()*vf.dimensions(),
+                rDeltaT.internalField()*
+                (
+                    alpha.internalField()
+                   *rho.internalField()
+                   *vf.internalField()
+
+                  - alpha.oldTime().internalField()
+                   *rho.oldTime().internalField()
+                   *vf.oldTime().internalField()*mesh().Vsc0()/mesh().Vsc()
+                ),
+                rDeltaT.boundaryField()*
+                (
+                    alpha.boundaryField()
+                   *rho.boundaryField()
+                   *vf.boundaryField()
+
+                  - alpha.oldTime().boundaryField()
+                   *rho.oldTime().boundaryField()
+                   *vf.oldTime().boundaryField()
+                )
+            )
+        );
+    }
+    else
+    {
+        return tmp<GeometricField<Type, fvPatchField, volMesh> >
+        (
+            new GeometricField<Type, fvPatchField, volMesh>
+            (
+                ddtIOobject,
+                rDeltaT
+               *(
+                   alpha*rho*vf
+                 - alpha.oldTime()*rho.oldTime()*vf.oldTime()
+                )
+            )
+        );
+    }
+}
+
+
 template<class Type>
 tmp<fvMatrix<Type> >
 SLTSDdtScheme<Type>::fvmDdt
@@ -481,6 +550,48 @@ SLTSDdtScheme<Type>::fvmDdt
 }
 
 
+template<class Type>
+tmp<fvMatrix<Type> >
+SLTSDdtScheme<Type>::fvmDdt
+(
+    const volScalarField& alpha,
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    tmp<fvMatrix<Type> > tfvm
+    (
+        new fvMatrix<Type>
+        (
+            vf,
+            alpha.dimensions()*rho.dimensions()*vf.dimensions()*dimVol/dimTime
+        )
+    );
+    fvMatrix<Type>& fvm = tfvm();
+
+    scalarField rDeltaT(SLrDeltaT()().internalField());
+
+    fvm.diag() = rDeltaT*alpha.internalField()*rho.internalField()*mesh().Vsc();
+
+    if (mesh().moving())
+    {
+        fvm.source() = rDeltaT
+            *alpha.oldTime().internalField()
+            *rho.oldTime().internalField()
+            *vf.oldTime().internalField()*mesh().Vsc0();
+    }
+    else
+    {
+        fvm.source() = rDeltaT
+            *alpha.oldTime().internalField()
+            *rho.oldTime().internalField()
+            *vf.oldTime().internalField()*mesh().Vsc();
+    }
+
+    return tfvm;
+}
+
+
 template<class Type>
 tmp<typename SLTSDdtScheme<Type>::fluxFieldType>
 SLTSDdtScheme<Type>::fvcDdtUfCorr
diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/SLTSDdtScheme/SLTSDdtScheme.H b/src/finiteVolume/finiteVolume/ddtSchemes/SLTSDdtScheme/SLTSDdtScheme.H
index b420f2b1614..ea6fdd987e5 100644
--- a/src/finiteVolume/finiteVolume/ddtSchemes/SLTSDdtScheme/SLTSDdtScheme.H
+++ b/src/finiteVolume/finiteVolume/ddtSchemes/SLTSDdtScheme/SLTSDdtScheme.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -140,6 +140,13 @@ public:
             const GeometricField<Type, fvPatchField, volMesh>&
         );
 
+        tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDdt
+        (
+            const volScalarField& alpha,
+            const volScalarField& rho,
+            const GeometricField<Type, fvPatchField, volMesh>& psi
+        );
+
         tmp<fvMatrix<Type> > fvmDdt
         (
             const GeometricField<Type, fvPatchField, volMesh>&
@@ -157,6 +164,13 @@ public:
             const GeometricField<Type, fvPatchField, volMesh>&
         );
 
+        tmp<fvMatrix<Type> > fvmDdt
+        (
+            const volScalarField& alpha,
+            const volScalarField& rho,
+            const GeometricField<Type, fvPatchField, volMesh>& psi
+        );
+
         typedef typename ddtScheme<Type>::fluxFieldType fluxFieldType;
 
         tmp<fluxFieldType> fvcDdtUfCorr
diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/backwardDdtScheme/backwardDdtScheme.C b/src/finiteVolume/finiteVolume/ddtSchemes/backwardDdtScheme/backwardDdtScheme.C
index 06116eb6edd..1474b8c89d9 100644
--- a/src/finiteVolume/finiteVolume/ddtSchemes/backwardDdtScheme/backwardDdtScheme.C
+++ b/src/finiteVolume/finiteVolume/ddtSchemes/backwardDdtScheme/backwardDdtScheme.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -281,6 +281,7 @@ backwardDdtScheme<Type>::fvcDdt
     }
 }
 
+
 template<class Type>
 tmp<GeometricField<Type, fvPatchField, volMesh> >
 backwardDdtScheme<Type>::fvcDdt
@@ -356,6 +357,100 @@ backwardDdtScheme<Type>::fvcDdt
 }
 
 
+template<class Type>
+tmp<GeometricField<Type, fvPatchField, volMesh> >
+backwardDdtScheme<Type>::fvcDdt
+(
+    const volScalarField& alpha,
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
+
+    IOobject ddtIOobject
+    (
+        "ddt("+alpha.name()+','+rho.name()+','+vf.name()+')',
+        mesh().time().timeName(),
+        mesh()
+    );
+
+    scalar deltaT = deltaT_();
+    scalar deltaT0 = deltaT0_(vf);
+
+    scalar coefft   = 1 + deltaT/(deltaT + deltaT0);
+    scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
+    scalar coefft0  = coefft + coefft00;
+
+    if (mesh().moving())
+    {
+        return tmp<GeometricField<Type, fvPatchField, volMesh> >
+        (
+            new GeometricField<Type, fvPatchField, volMesh>
+            (
+                ddtIOobject,
+                mesh(),
+                rDeltaT.dimensions()
+               *alpha.dimensions()*rho.dimensions()*vf.dimensions(),
+                rDeltaT.value()*
+                (
+                    coefft
+                   *alpha.internalField()
+                   *rho.internalField()
+                   *vf.internalField() -
+                    (
+                        coefft0
+                       *alpha.oldTime().internalField()
+                       *rho.oldTime().internalField()
+                       *vf.oldTime().internalField()*mesh().V0()
+
+                      - coefft00
+                       *alpha.oldTime().oldTime().internalField()
+                       *rho.oldTime().oldTime().internalField()
+                       *vf.oldTime().oldTime().internalField()*mesh().V00()
+                    )/mesh().V()
+                ),
+                rDeltaT.value()*
+                (
+                    coefft
+                   *alpha.boundaryField()
+                   *rho.boundaryField()
+                   *vf.boundaryField() -
+                    (
+                        coefft0
+                       *alpha.oldTime().boundaryField()
+                       *rho.oldTime().boundaryField()
+                       *vf.oldTime().boundaryField()
+
+                      - coefft00
+                       *alpha.oldTime().oldTime().boundaryField()
+                       *rho.oldTime().oldTime().boundaryField()
+                       *vf.oldTime().oldTime().boundaryField()
+                    )
+                )
+            )
+        );
+    }
+    else
+    {
+        return tmp<GeometricField<Type, fvPatchField, volMesh> >
+        (
+            new GeometricField<Type, fvPatchField, volMesh>
+            (
+                ddtIOobject,
+                rDeltaT*
+                (
+                    coefft*alpha*rho*vf
+                  - coefft0*alpha.oldTime()*rho.oldTime()*vf.oldTime()
+                  + coefft00*alpha.oldTime().oldTime()
+                   *rho.oldTime().oldTime()*vf.oldTime().oldTime()
+                )
+            )
+        );
+    }
+}
+
+
 template<class Type>
 tmp<fvMatrix<Type> >
 backwardDdtScheme<Type>::fvmDdt
@@ -512,6 +607,72 @@ backwardDdtScheme<Type>::fvmDdt
 }
 
 
+template<class Type>
+tmp<fvMatrix<Type> >
+backwardDdtScheme<Type>::fvmDdt
+(
+    const volScalarField& alpha,
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    tmp<fvMatrix<Type> > tfvm
+    (
+        new fvMatrix<Type>
+        (
+            vf,
+            alpha.dimensions()*rho.dimensions()*vf.dimensions()*dimVol/dimTime
+        )
+    );
+    fvMatrix<Type>& fvm = tfvm();
+
+    scalar rDeltaT = 1.0/deltaT_();
+
+    scalar deltaT = deltaT_();
+    scalar deltaT0 = deltaT0_(vf);
+
+    scalar coefft   = 1 + deltaT/(deltaT + deltaT0);
+    scalar coefft00 = deltaT*deltaT/(deltaT0*(deltaT + deltaT0));
+    scalar coefft0  = coefft + coefft00;
+
+    fvm.diag() =
+        (coefft*rDeltaT)*alpha.internalField()*rho.internalField()*mesh().V();
+
+    if (mesh().moving())
+    {
+        fvm.source() = rDeltaT*
+        (
+            coefft0
+           *alpha.oldTime().internalField()
+           *rho.oldTime().internalField()
+           *vf.oldTime().internalField()*mesh().V0()
+
+          - coefft00
+           *alpha.oldTime().oldTime().internalField()
+           *rho.oldTime().oldTime().internalField()
+           *vf.oldTime().oldTime().internalField()*mesh().V00()
+        );
+    }
+    else
+    {
+        fvm.source() = rDeltaT*mesh().V()*
+        (
+            coefft0
+           *alpha.oldTime().internalField()
+           *rho.oldTime().internalField()
+           *vf.oldTime().internalField()
+
+          - coefft00
+           *alpha.oldTime().oldTime().internalField()
+           *rho.oldTime().oldTime().internalField()
+           *vf.oldTime().oldTime().internalField()
+        );
+    }
+
+    return tfvm;
+}
+
+
 template<class Type>
 tmp<typename backwardDdtScheme<Type>::fluxFieldType>
 backwardDdtScheme<Type>::fvcDdtUfCorr
diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/backwardDdtScheme/backwardDdtScheme.H b/src/finiteVolume/finiteVolume/ddtSchemes/backwardDdtScheme/backwardDdtScheme.H
index 425e78ece65..9855ad55997 100644
--- a/src/finiteVolume/finiteVolume/ddtSchemes/backwardDdtScheme/backwardDdtScheme.H
+++ b/src/finiteVolume/finiteVolume/ddtSchemes/backwardDdtScheme/backwardDdtScheme.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -128,6 +128,13 @@ public:
             const GeometricField<Type, fvPatchField, volMesh>&
         );
 
+        tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDdt
+        (
+            const volScalarField& alpha,
+            const volScalarField& rho,
+            const GeometricField<Type, fvPatchField, volMesh>& psi
+        );
+
         tmp<fvMatrix<Type> > fvmDdt
         (
             const GeometricField<Type, fvPatchField, volMesh>&
@@ -145,6 +152,13 @@ public:
             const GeometricField<Type, fvPatchField, volMesh>&
         );
 
+        tmp<fvMatrix<Type> > fvmDdt
+        (
+            const volScalarField& alpha,
+            const volScalarField& rho,
+            const GeometricField<Type, fvPatchField, volMesh>& psi
+        );
+
         typedef typename ddtScheme<Type>::fluxFieldType fluxFieldType;
 
         tmp<fluxFieldType> fvcDdtUfCorr
diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/boundedDdtScheme/boundedDdtScheme.C b/src/finiteVolume/finiteVolume/ddtSchemes/boundedDdtScheme/boundedDdtScheme.C
index 2edefd19b29..d5ea7dbcbb1 100644
--- a/src/finiteVolume/finiteVolume/ddtSchemes/boundedDdtScheme/boundedDdtScheme.C
+++ b/src/finiteVolume/finiteVolume/ddtSchemes/boundedDdtScheme/boundedDdtScheme.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -87,6 +87,19 @@ boundedDdtScheme<Type>::fvcDdt
 }
 
 
+template<class Type>
+tmp<GeometricField<Type, fvPatchField, volMesh> >
+boundedDdtScheme<Type>::fvcDdt
+(
+    const volScalarField& alpha,
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    return scheme_().fvcDdt(alpha, rho, vf) - fvc::ddt(alpha, rho)*vf;
+}
+
+
 template<class Type>
 tmp<fvMatrix<Type> >
 boundedDdtScheme<Type>::fvmDdt
@@ -122,6 +135,19 @@ boundedDdtScheme<Type>::fvmDdt
 }
 
 
+template<class Type>
+tmp<fvMatrix<Type> >
+boundedDdtScheme<Type>::fvmDdt
+(
+    const volScalarField& alpha,
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    return scheme_().fvmDdt(alpha, rho, vf) - fvm::Sp(fvc::ddt(alpha, rho), vf);
+}
+
+
 template<class Type>
 tmp<typename boundedDdtScheme<Type>::fluxFieldType>
 boundedDdtScheme<Type>::fvcDdtUfCorr
diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/boundedDdtScheme/boundedDdtScheme.H b/src/finiteVolume/finiteVolume/ddtSchemes/boundedDdtScheme/boundedDdtScheme.H
index db309c8c815..741e7c1e0fa 100644
--- a/src/finiteVolume/finiteVolume/ddtSchemes/boundedDdtScheme/boundedDdtScheme.H
+++ b/src/finiteVolume/finiteVolume/ddtSchemes/boundedDdtScheme/boundedDdtScheme.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -126,6 +126,13 @@ public:
             const GeometricField<Type, fvPatchField, volMesh>&
         );
 
+        tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDdt
+        (
+            const volScalarField& alpha,
+            const volScalarField& rho,
+            const GeometricField<Type, fvPatchField, volMesh>& psi
+        );
+
         tmp<fvMatrix<Type> > fvmDdt
         (
             const GeometricField<Type, fvPatchField, volMesh>&
@@ -143,6 +150,13 @@ public:
             const GeometricField<Type, fvPatchField, volMesh>&
         );
 
+        tmp<fvMatrix<Type> > fvmDdt
+        (
+            const volScalarField& alpha,
+            const volScalarField& rho,
+            const GeometricField<Type, fvPatchField, volMesh>& psi
+        );
+
         typedef typename ddtScheme<Type>::fluxFieldType fluxFieldType;
 
         tmp<fluxFieldType> fvcDdtUfCorr
diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/ddtScheme/ddtScheme.C b/src/finiteVolume/finiteVolume/ddtSchemes/ddtScheme/ddtScheme.C
index ccfa56538e7..c8e0ce6fa6f 100644
--- a/src/finiteVolume/finiteVolume/ddtSchemes/ddtScheme/ddtScheme.C
+++ b/src/finiteVolume/finiteVolume/ddtSchemes/ddtScheme/ddtScheme.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -26,6 +26,7 @@ License
 #include "fv.H"
 #include "HashTable.H"
 #include "surfaceInterpolate.H"
+#include "fvMatrix.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -95,6 +96,45 @@ ddtScheme<Type>::~ddtScheme()
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
+template<class Type>
+tmp<GeometricField<Type, fvPatchField, volMesh> > ddtScheme<Type>::fvcDdt
+(
+    const volScalarField& alpha,
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    notImplemented("fvmDdt(alpha, rho, psi");
+
+    return tmp<GeometricField<Type, fvPatchField, volMesh> >
+    (
+        GeometricField<Type, fvPatchField, volMesh>::null()
+    );
+}
+
+
+template<class Type>
+tmp<fvMatrix<Type> > ddtScheme<Type>::fvmDdt
+(
+    const volScalarField& alpha,
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    notImplemented("fvmDdt(alpha, rho, psi");
+
+    return tmp<fvMatrix<Type> >
+    (
+        new fvMatrix<Type>
+        (
+            vf,
+            alpha.dimensions()*rho.dimensions()
+            *vf.dimensions()*dimVol/dimTime
+        )
+    );
+}
+
+
 template<class Type>
 tmp<surfaceScalarField> ddtScheme<Type>::fvcDdtPhiCoeff
 (
diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/ddtScheme/ddtScheme.H b/src/finiteVolume/finiteVolume/ddtSchemes/ddtScheme/ddtScheme.H
index ceb6b1d63ac..dccf7f9619d 100644
--- a/src/finiteVolume/finiteVolume/ddtSchemes/ddtScheme/ddtScheme.H
+++ b/src/finiteVolume/finiteVolume/ddtSchemes/ddtScheme/ddtScheme.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -160,6 +160,13 @@ public:
             const GeometricField<Type, fvPatchField, volMesh>&
         ) = 0;
 
+        virtual tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDdt
+        (
+            const volScalarField& alpha,
+            const volScalarField& rho,
+            const GeometricField<Type, fvPatchField, volMesh>&
+        ) = 0;
+
         virtual tmp<fvMatrix<Type> > fvmDdt
         (
             const GeometricField<Type, fvPatchField, volMesh>&
@@ -177,6 +184,12 @@ public:
             const GeometricField<Type, fvPatchField, volMesh>&
         ) = 0;
 
+        virtual tmp<fvMatrix<Type> > fvmDdt
+        (
+            const volScalarField& alpha,
+            const volScalarField& rho,
+            const GeometricField<Type, fvPatchField, volMesh>& vf
+        ) = 0;
 
         typedef GeometricField
         <
diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/localEulerDdtScheme/localEulerDdtScheme.C b/src/finiteVolume/finiteVolume/ddtSchemes/localEulerDdtScheme/localEulerDdtScheme.C
index ed5c3f0f349..e1c0ea5e482 100644
--- a/src/finiteVolume/finiteVolume/ddtSchemes/localEulerDdtScheme/localEulerDdtScheme.C
+++ b/src/finiteVolume/finiteVolume/ddtSchemes/localEulerDdtScheme/localEulerDdtScheme.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -265,6 +265,75 @@ localEulerDdtScheme<Type>::fvcDdt
 }
 
 
+template<class Type>
+tmp<GeometricField<Type, fvPatchField, volMesh> >
+localEulerDdtScheme<Type>::fvcDdt
+(
+    const volScalarField& alpha,
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    const volScalarField& rDeltaT = localRDeltaT();
+
+    IOobject ddtIOobject
+    (
+        "ddt("+alpha.name()+','+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()
+               *alpha.dimensions()*rho.dimensions()*vf.dimensions(),
+                rDeltaT.internalField()*
+                (
+                    alpha.internalField()
+                   *rho.internalField()
+                   *vf.internalField()
+
+                  - alpha.oldTime().internalField()
+                   *rho.oldTime().internalField()
+                   *vf.oldTime().internalField()*mesh().Vsc0()/mesh().Vsc()
+                ),
+                rDeltaT.boundaryField()*
+                (
+                    alpha.boundaryField()
+                   *rho.boundaryField()
+                   *vf.boundaryField()
+
+                  - alpha.oldTime().boundaryField()
+                   *rho.oldTime().boundaryField()
+                   *vf.oldTime().boundaryField()
+                )
+            )
+        );
+    }
+    else
+    {
+        return tmp<GeometricField<Type, fvPatchField, volMesh> >
+        (
+            new GeometricField<Type, fvPatchField, volMesh>
+            (
+                ddtIOobject,
+                rDeltaT
+               *(
+                   alpha*rho*vf
+                 - alpha.oldTime()*rho.oldTime()*vf.oldTime()
+                )
+            )
+        );
+    }
+}
+
+
 template<class Type>
 tmp<fvMatrix<Type> >
 localEulerDdtScheme<Type>::fvmDdt
@@ -376,6 +445,48 @@ localEulerDdtScheme<Type>::fvmDdt
 }
 
 
+template<class Type>
+tmp<fvMatrix<Type> >
+localEulerDdtScheme<Type>::fvmDdt
+(
+    const volScalarField& alpha,
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    tmp<fvMatrix<Type> > tfvm
+    (
+        new fvMatrix<Type>
+        (
+            vf,
+            alpha.dimensions()*rho.dimensions()*vf.dimensions()*dimVol/dimTime
+        )
+    );
+    fvMatrix<Type>& fvm = tfvm();
+
+    const scalarField& rDeltaT = localRDeltaT().internalField();
+
+    fvm.diag() = rDeltaT*alpha.internalField()*rho.internalField()*mesh().Vsc();
+
+    if (mesh().moving())
+    {
+        fvm.source() = rDeltaT
+            *alpha.oldTime().internalField()
+            *rho.oldTime().internalField()
+            *vf.oldTime().internalField()*mesh().Vsc0();
+    }
+    else
+    {
+        fvm.source() = rDeltaT
+            *alpha.oldTime().internalField()
+            *rho.oldTime().internalField()
+            *vf.oldTime().internalField()*mesh().Vsc();
+    }
+
+    return tfvm;
+}
+
+
 template<class Type>
 tmp<typename localEulerDdtScheme<Type>::fluxFieldType>
 localEulerDdtScheme<Type>::fvcDdtUfCorr
diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/localEulerDdtScheme/localEulerDdtScheme.H b/src/finiteVolume/finiteVolume/ddtSchemes/localEulerDdtScheme/localEulerDdtScheme.H
index 50aed99d523..dc7bb0c79ba 100644
--- a/src/finiteVolume/finiteVolume/ddtSchemes/localEulerDdtScheme/localEulerDdtScheme.H
+++ b/src/finiteVolume/finiteVolume/ddtSchemes/localEulerDdtScheme/localEulerDdtScheme.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -136,6 +136,13 @@ public:
             const GeometricField<Type, fvPatchField, volMesh>&
         );
 
+        tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDdt
+        (
+            const volScalarField& alpha,
+            const volScalarField& rho,
+            const GeometricField<Type, fvPatchField, volMesh>& psi
+        );
+
         tmp<fvMatrix<Type> > fvmDdt
         (
             const GeometricField<Type, fvPatchField, volMesh>&
@@ -153,6 +160,13 @@ public:
             const GeometricField<Type, fvPatchField, volMesh>&
         );
 
+        tmp<fvMatrix<Type> > fvmDdt
+        (
+            const volScalarField& alpha,
+            const volScalarField& rho,
+            const GeometricField<Type, fvPatchField, volMesh>& psi
+        );
+
         typedef typename ddtScheme<Type>::fluxFieldType fluxFieldType;
 
         tmp<fluxFieldType> fvcDdtUfCorr
diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/steadyStateDdtScheme/steadyStateDdtScheme.C b/src/finiteVolume/finiteVolume/ddtSchemes/steadyStateDdtScheme/steadyStateDdtScheme.C
index f21a6a0b6f5..fb695f1b4d8 100644
--- a/src/finiteVolume/finiteVolume/ddtSchemes/steadyStateDdtScheme/steadyStateDdtScheme.C
+++ b/src/finiteVolume/finiteVolume/ddtSchemes/steadyStateDdtScheme/steadyStateDdtScheme.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -157,6 +157,37 @@ steadyStateDdtScheme<Type>::fvcDdt
 }
 
 
+template<class Type>
+tmp<GeometricField<Type, fvPatchField, volMesh> >
+steadyStateDdtScheme<Type>::fvcDdt
+(
+    const volScalarField& alpha,
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    return tmp<GeometricField<Type, fvPatchField, volMesh> >
+    (
+        new GeometricField<Type, fvPatchField, volMesh>
+        (
+            IOobject
+            (
+                "ddt("+alpha.name()+','+rho.name()+','+vf.name()+')',
+                mesh().time().timeName(),
+                mesh()
+            ),
+            mesh(),
+            dimensioned<Type>
+            (
+                "0",
+                rho.dimensions()*vf.dimensions()/dimTime,
+                pTraits<Type>::zero
+            )
+        )
+    );
+}
+
+
 template<class Type>
 tmp<fvMatrix<Type> >
 steadyStateDdtScheme<Type>::fvmDdt
@@ -219,6 +250,28 @@ steadyStateDdtScheme<Type>::fvmDdt
 }
 
 
+template<class Type>
+tmp<fvMatrix<Type> >
+steadyStateDdtScheme<Type>::fvmDdt
+(
+    const volScalarField& alpha,
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    tmp<fvMatrix<Type> > tfvm
+    (
+        new fvMatrix<Type>
+        (
+            vf,
+            alpha.dimensions()*rho.dimensions()*vf.dimensions()*dimVol/dimTime
+        )
+    );
+
+    return tfvm;
+}
+
+
 template<class Type>
 tmp<typename steadyStateDdtScheme<Type>::fluxFieldType>
 steadyStateDdtScheme<Type>::fvcDdtUfCorr
diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/steadyStateDdtScheme/steadyStateDdtScheme.H b/src/finiteVolume/finiteVolume/ddtSchemes/steadyStateDdtScheme/steadyStateDdtScheme.H
index 4512d053102..7d990746f2e 100644
--- a/src/finiteVolume/finiteVolume/ddtSchemes/steadyStateDdtScheme/steadyStateDdtScheme.H
+++ b/src/finiteVolume/finiteVolume/ddtSchemes/steadyStateDdtScheme/steadyStateDdtScheme.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -116,6 +116,13 @@ public:
             const GeometricField<Type, fvPatchField, volMesh>&
         );
 
+        tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDdt
+        (
+            const volScalarField& alpha,
+            const volScalarField& rho,
+            const GeometricField<Type, fvPatchField, volMesh>& psi
+        );
+
         tmp<fvMatrix<Type> > fvmDdt
         (
             const GeometricField<Type, fvPatchField, volMesh>&
@@ -133,6 +140,13 @@ public:
             const GeometricField<Type, fvPatchField, volMesh>&
         );
 
+        tmp<fvMatrix<Type> > fvmDdt
+        (
+            const volScalarField& alpha,
+            const volScalarField& rho,
+            const GeometricField<Type, fvPatchField, volMesh>& psi
+        );
+
         typedef typename ddtScheme<Type>::fluxFieldType fluxFieldType;
 
         tmp<fluxFieldType> fvcDdtUfCorr
diff --git a/src/finiteVolume/finiteVolume/fvc/fvcDdt.C b/src/finiteVolume/finiteVolume/fvc/fvcDdt.C
index 4026ceeb9d0..e300c072d49 100644
--- a/src/finiteVolume/finiteVolume/fvc/fvcDdt.C
+++ b/src/finiteVolume/finiteVolume/fvc/fvcDdt.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -102,6 +102,29 @@ ddt
 }
 
 
+template<class Type>
+tmp<GeometricField<Type, fvPatchField, volMesh> >
+ddt
+(
+    const volScalarField& alpha,
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    return fv::ddtScheme<Type>::New
+    (
+        vf.mesh(),
+        vf.mesh().ddtScheme
+        (
+            "ddt("
+          + alpha.name() + ','
+          + rho.name() + ','
+          + vf.name() + ')'
+        )
+    )().fvcDdt(alpha, rho, vf);
+}
+
+
 template<class Type>
 tmp<GeometricField<Type, fvPatchField, volMesh> >
 ddt
diff --git a/src/finiteVolume/finiteVolume/fvc/fvcDdt.H b/src/finiteVolume/finiteVolume/fvc/fvcDdt.H
index 5aca4724690..1f7047cdec7 100644
--- a/src/finiteVolume/finiteVolume/fvc/fvcDdt.H
+++ b/src/finiteVolume/finiteVolume/fvc/fvcDdt.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -80,6 +80,14 @@ namespace fvc
         const GeometricField<Type, fvPatchField, volMesh>&
     );
 
+    template<class Type>
+    tmp<GeometricField<Type, fvPatchField, volMesh> > ddt
+    (
+        const volScalarField&,
+        const volScalarField&,
+        const GeometricField<Type, fvPatchField, volMesh>&
+    );
+
     template<class Type>
     tmp<GeometricField<Type, fvPatchField, volMesh> > ddt
     (
diff --git a/src/finiteVolume/finiteVolume/fvm/fvmDdt.C b/src/finiteVolume/finiteVolume/fvm/fvmDdt.C
index 9a73ae36958..a2f75da0242 100644
--- a/src/finiteVolume/finiteVolume/fvm/fvmDdt.C
+++ b/src/finiteVolume/finiteVolume/fvm/fvmDdt.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -99,6 +99,29 @@ ddt
 }
 
 
+template<class Type>
+tmp<fvMatrix<Type> >
+ddt
+(
+    const volScalarField& alpha,
+    const volScalarField& rho,
+    const GeometricField<Type, fvPatchField, volMesh>& vf
+)
+{
+    return fv::ddtScheme<Type>::New
+    (
+        vf.mesh(),
+        vf.mesh().ddtScheme
+        (
+            "ddt("
+          + alpha.name() + ','
+          + rho.name() + ','
+          + vf.name() + ')'
+        )
+    )().fvmDdt(alpha, rho, vf);
+}
+
+
 template<class Type>
 tmp<fvMatrix<Type> >
 ddt
diff --git a/src/finiteVolume/finiteVolume/fvm/fvmDdt.H b/src/finiteVolume/finiteVolume/fvm/fvmDdt.H
index efc28a2cca2..4091c5be6b1 100644
--- a/src/finiteVolume/finiteVolume/fvm/fvmDdt.H
+++ b/src/finiteVolume/finiteVolume/fvm/fvmDdt.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -77,6 +77,14 @@ namespace fvm
         const GeometricField<Type, fvPatchField, volMesh>&
     );
 
+    template<class Type>
+    tmp<fvMatrix<Type> > ddt
+    (
+        const volScalarField&,
+        const volScalarField&,
+        const GeometricField<Type, fvPatchField, volMesh>&
+    );
+
     template<class Type>
     tmp<fvMatrix<Type> > ddt
     (
-- 
GitLab