From efa306a033b7d9f4f19d0b16ba656e8aa7082228 Mon Sep 17 00:00:00 2001
From: Henry Weller <http://cfd.direct>
Date: Fri, 10 Mar 2017 11:03:42 +0000
Subject: [PATCH] turbulenceModels::LES::dynamicLagrangian: Converted advection
 to compressible convection form

Formally this is equivalent to the previous formulation but more convenient to
use given that for compressible flow the mass flux rather than the volume flux
is available.
---
 .../LES/dynamicLagrangian/dynamicLagrangian.C | 19 ++++++++++---------
 1 file changed, 10 insertions(+), 9 deletions(-)

diff --git a/src/TurbulenceModels/turbulenceModels/LES/dynamicLagrangian/dynamicLagrangian.C b/src/TurbulenceModels/turbulenceModels/LES/dynamicLagrangian/dynamicLagrangian.C
index 00154232e..28ad1eb75 100644
--- a/src/TurbulenceModels/turbulenceModels/LES/dynamicLagrangian/dynamicLagrangian.C
+++ b/src/TurbulenceModels/turbulenceModels/LES/dynamicLagrangian/dynamicLagrangian.C
@@ -159,8 +159,9 @@ void dynamicLagrangian<BasicTurbulenceModel>::correct()
     }
 
     // Local references
-    tmp<surfaceScalarField> tphi = this->phi();
-    const surfaceScalarField& phi = tphi();
+    const alphaField& alpha = this->alpha_;
+    const rhoField& rho = this->rho_;
+    const surfaceScalarField& alphaRhoPhi = this->alphaRhoPhi_;
     const volVectorField& U = this->U_;
     fv::options& fvOptions(fv::options::New(this->mesh_));
 
@@ -184,19 +185,19 @@ void dynamicLagrangian<BasicTurbulenceModel>::correct()
 
     volScalarField invT
     (
-        (1.0/(theta_.value()*this->delta()))*pow(flm_*fmm_, 1.0/8.0)
+        alpha*rho*(1.0/(theta_.value()*this->delta()))*pow(flm_*fmm_, 1.0/8.0)
     );
 
     volScalarField LM(L && M);
 
     fvScalarMatrix flmEqn
     (
-        fvm::ddt(flm_)
-      + fvm::div(phi, flm_)
+        fvm::ddt(alpha, rho, flm_)
+      + fvm::div(alphaRhoPhi, flm_)
      ==
         invT*LM
       - fvm::Sp(invT, flm_)
-      + fvOptions(flm_)
+      + fvOptions(alpha, rho, flm_)
     );
 
     flmEqn.relax();
@@ -209,12 +210,12 @@ void dynamicLagrangian<BasicTurbulenceModel>::correct()
 
     fvScalarMatrix fmmEqn
     (
-        fvm::ddt(fmm_)
-      + fvm::div(phi, fmm_)
+        fvm::ddt(alpha, rho, fmm_)
+      + fvm::div(alphaRhoPhi, fmm_)
      ==
         invT*MM
       - fvm::Sp(invT, fmm_)
-      + fvOptions(fmm_)
+      + fvOptions(alpha, rho, fmm_)
     );
 
     fmmEqn.relax();