From b162c4818b9bbdb824db9e75baffef6ac7c5a7a0 Mon Sep 17 00:00:00 2001
From: Henry Weller <http://openfoam.org>
Date: Thu, 21 Feb 2019 13:54:41 +0000
Subject: [PATCH] ddtScheme::fvcDdtPhiCoeff: Reverted to previous
 flux-normalised scheme

The current Courant number limited formulation is preferable for cases running
with VERY small Courant numbers but there are stability issues with it on some
cases for which the previous flux-normalised scheme worked well.  Overall the
previous scheme has proved more reliable and fvcDdtPhiCoeff has been reverted to
it pending further research and development in this area which will require
funding:

https://openfoam.org/news/funding-2019/
https://openfoam.org/maintenance/
---
 .../ddtSchemes/ddtScheme/ddtScheme.C           | 18 +++++++++++++++---
 .../localEulerDdtScheme/localEulerDdtScheme.C  |  6 +++++-
 .../localEulerDdtScheme/localEulerDdtScheme.H  | 16 ++++++++--------
 3 files changed, 28 insertions(+), 12 deletions(-)

diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/ddtScheme/ddtScheme.C b/src/finiteVolume/finiteVolume/ddtSchemes/ddtScheme/ddtScheme.C
index 28ec0eeee7..eb69175ac4 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     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2019 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -156,6 +156,8 @@ tmp<surfaceScalarField> ddtScheme<Type>::fvcDdtPhiCoeff
     const fluxFieldType& phiCorr
 )
 {
+    // Courant number limited formulation
+    /*
     tmp<surfaceScalarField> tddtCouplingCoeff = scalar(1)
       - min
         (
@@ -163,6 +165,16 @@ tmp<surfaceScalarField> ddtScheme<Type>::fvcDdtPhiCoeff
            *mesh().time().deltaT()*mesh().deltaCoeffs()/mesh().magSf(),
             scalar(1)
         );
+    */
+
+    // Flux normalised formulation
+    tmp<surfaceScalarField> tddtCouplingCoeff = scalar(1)
+      - min
+        (
+            mag(phiCorr)
+           /(mag(phi) + dimensionedScalar("small", phi.dimensions(), SMALL)),
+            scalar(1)
+        );
 
     surfaceScalarField& ddtCouplingCoeff = tddtCouplingCoeff.ref();
 
@@ -203,7 +215,7 @@ tmp<surfaceScalarField> ddtScheme<Type>::fvcDdtPhiCoeff
     const volScalarField& rho
 )
 {
-    return fvcDdtPhiCoeff(U, phi, phiCorr/fvc::interpolate(rho));
+    return fvcDdtPhiCoeff(U, phi, phiCorr);
 }
 
 
@@ -230,7 +242,7 @@ tmp<surfaceScalarField> ddtScheme<Type>::fvcDdtPhiCoeff
     (
         U,
         phi,
-        (phi - fvc::dotInterpolate(mesh().Sf(), U))/fvc::interpolate(rho)
+        (phi - fvc::dotInterpolate(mesh().Sf(), U))
     );
 }
 
diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/localEulerDdtScheme/localEulerDdtScheme.C b/src/finiteVolume/finiteVolume/ddtSchemes/localEulerDdtScheme/localEulerDdtScheme.C
index 52e77d20d5..75ef442deb 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     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2019 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -347,6 +347,8 @@ localEulerDdtScheme<Type>::fvmDdt
 }
 
 
+/*
+// Courant number limited formulation
 template<class Type>
 tmp<surfaceScalarField> localEulerDdtScheme<Type>::fvcDdtPhiCoeff
 (
@@ -355,6 +357,7 @@ tmp<surfaceScalarField> localEulerDdtScheme<Type>::fvcDdtPhiCoeff
     const fluxFieldType& phiCorr
 )
 {
+    // Courant number limited formulation
     tmp<surfaceScalarField> tddtCouplingCoeff = scalar(1)
       - min
         (
@@ -391,6 +394,7 @@ tmp<surfaceScalarField> localEulerDdtScheme<Type>::fvcDdtPhiCoeff
 
     return tddtCouplingCoeff;
 }
+*/
 
 
 template<class Type>
diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/localEulerDdtScheme/localEulerDdtScheme.H b/src/finiteVolume/finiteVolume/ddtSchemes/localEulerDdtScheme/localEulerDdtScheme.H
index edb911607e..c5f3ded8f8 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     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2019 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -215,14 +215,14 @@ public:
 
         typedef typename ddtScheme<Type>::fluxFieldType fluxFieldType;
 
-        using ddtScheme<Type>::fvcDdtPhiCoeff;
+        // using ddtScheme<Type>::fvcDdtPhiCoeff;
 
-        virtual tmp<surfaceScalarField> fvcDdtPhiCoeff
-        (
-            const GeometricField<Type, fvPatchField, volMesh>& U,
-            const fluxFieldType& phi,
-            const fluxFieldType& phiCorr
-        );
+        // virtual tmp<surfaceScalarField> fvcDdtPhiCoeff
+        // (
+        //     const GeometricField<Type, fvPatchField, volMesh>& U,
+        //     const fluxFieldType& phi,
+        //     const fluxFieldType& phiCorr
+        // );
 
         virtual tmp<fluxFieldType> fvcDdtUfCorr
         (