From 9b99bbf312d096dc36d35d69492034cf4b3f25b8 Mon Sep 17 00:00:00 2001
From: Henry <Henry>
Date: Sat, 28 Feb 2015 16:11:57 +0000
Subject: [PATCH] Reynolds stress turbulence models: Add laminar diffusion of R
 and epsilon Also use the innerSqr function where appropriate

---
 .../turbulenceModels/RAS/LRR/LRR.C            | 32 +++++++++++++++--
 .../turbulenceModels/RAS/LRR/LRR.H            |  6 ++++
 .../turbulenceModels/RAS/SSG/SSG.C            | 36 ++++++++++++++++---
 .../turbulenceModels/RAS/SSG/SSG.H            |  6 ++++
 4 files changed, 74 insertions(+), 6 deletions(-)

diff --git a/src/TurbulenceModels/turbulenceModels/RAS/LRR/LRR.C b/src/TurbulenceModels/turbulenceModels/RAS/LRR/LRR.C
index 37391cc37..7dfa45a03 100644
--- a/src/TurbulenceModels/turbulenceModels/RAS/LRR/LRR.C
+++ b/src/TurbulenceModels/turbulenceModels/RAS/LRR/LRR.C
@@ -238,6 +238,34 @@ bool LRR<BasicTurbulenceModel>::read()
 }
 
 
+template<class BasicTurbulenceModel>
+tmp<volSymmTensorField> LRR<BasicTurbulenceModel>::DREff() const
+{
+    return tmp<volSymmTensorField>
+    (
+        new volSymmTensorField
+        (
+            "DREff",
+            (Cs_*(this->k_/this->epsilon_))*this->R_ + I*this->nu()
+        )
+    );
+}
+
+
+template<class BasicTurbulenceModel>
+tmp<volSymmTensorField> LRR<BasicTurbulenceModel>::DepsilonEff() const
+{
+    return tmp<volSymmTensorField>
+    (
+        new volSymmTensorField
+        (
+            "DepsilonEff",
+            (Ceps_*(this->k_/this->epsilon_))*this->R_ + I*this->nu()
+        )
+    );
+}
+
+
 template<class BasicTurbulenceModel>
 void LRR<BasicTurbulenceModel>::correct()
 {
@@ -269,7 +297,7 @@ void LRR<BasicTurbulenceModel>::correct()
     (
         fvm::ddt(alpha, rho, epsilon_)
       + fvm::div(alphaRhoPhi, epsilon_)
-      - fvm::laplacian(Ceps_*alpha*rho*(k_/epsilon_)*R, epsilon_)
+      - fvm::laplacian(alpha*rho*DepsilonEff(), epsilon_)
      ==
         Ceps1_*alpha*rho*G*epsilon_/k_
       - fvm::Sp(Ceps2_*alpha*rho*epsilon_/k_, epsilon_)
@@ -310,7 +338,7 @@ void LRR<BasicTurbulenceModel>::correct()
     (
         fvm::ddt(alpha, rho, R)
       + fvm::div(alphaRhoPhi, R)
-      - fvm::laplacian(Cs_*alpha*rho*(k_/epsilon_)*R, R)
+      - fvm::laplacian(alpha*rho*DREff(), R)
       + fvm::Sp(C1_*alpha*rho*epsilon_/k_, R)
       ==
         alpha*rho*P
diff --git a/src/TurbulenceModels/turbulenceModels/RAS/LRR/LRR.H b/src/TurbulenceModels/turbulenceModels/RAS/LRR/LRR.H
index 2af0da241..68c9bd2d7 100644
--- a/src/TurbulenceModels/turbulenceModels/RAS/LRR/LRR.H
+++ b/src/TurbulenceModels/turbulenceModels/RAS/LRR/LRR.H
@@ -195,6 +195,12 @@ public:
             return epsilon_;
         }
 
+        //- Return the effective diffusivity for R
+        tmp<volSymmTensorField> DREff() const;
+
+        //- Return the effective diffusivity for epsilon
+        tmp<volSymmTensorField> DepsilonEff() const;
+
         //- Solve the turbulence equations and correct eddy-Viscosity and
         //  related properties
         virtual void correct();
diff --git a/src/TurbulenceModels/turbulenceModels/RAS/SSG/SSG.C b/src/TurbulenceModels/turbulenceModels/RAS/SSG/SSG.C
index 720c1f413..4d61b3db6 100644
--- a/src/TurbulenceModels/turbulenceModels/RAS/SSG/SSG.C
+++ b/src/TurbulenceModels/turbulenceModels/RAS/SSG/SSG.C
@@ -248,6 +248,34 @@ bool SSG<BasicTurbulenceModel>::read()
 }
 
 
+template<class BasicTurbulenceModel>
+tmp<volSymmTensorField> SSG<BasicTurbulenceModel>::DREff() const
+{
+    return tmp<volSymmTensorField>
+    (
+        new volSymmTensorField
+        (
+            "DREff",
+            (Cs_*(this->k_/this->epsilon_))*this->R_ + I*this->nu()
+        )
+    );
+}
+
+
+template<class BasicTurbulenceModel>
+tmp<volSymmTensorField> SSG<BasicTurbulenceModel>::DepsilonEff() const
+{
+    return tmp<volSymmTensorField>
+    (
+        new volSymmTensorField
+        (
+            "DepsilonEff",
+            (Ceps_*(this->k_/this->epsilon_))*this->R_ + I*this->nu()
+        )
+    );
+}
+
+
 template<class BasicTurbulenceModel>
 void SSG<BasicTurbulenceModel>::correct()
 {
@@ -279,7 +307,7 @@ void SSG<BasicTurbulenceModel>::correct()
     (
         fvm::ddt(alpha, rho, epsilon_)
       + fvm::div(alphaRhoPhi, epsilon_)
-      - fvm::laplacian(Ceps_*alpha*rho*(k_/epsilon_)*R, epsilon_)
+      - fvm::laplacian(alpha*rho*DepsilonEff(), epsilon_)
      ==
         Ceps1_*alpha*rho*G*epsilon_/k_
       - fvm::Sp(Ceps2_*alpha*rho*epsilon_/k_, epsilon_)
@@ -324,15 +352,15 @@ void SSG<BasicTurbulenceModel>::correct()
     (
         fvm::ddt(alpha, rho, R)
       + fvm::div(alphaRhoPhi, R)
-      - fvm::laplacian(Cs_*alpha*rho*(k_/epsilon_)*R, R)
+      - fvm::laplacian(alpha*rho*DREff(), R)
       + fvm::Sp(((C1_/2)*epsilon_ + (C1s_/2)*G)*alpha*rho/k_, R)
      ==
         alpha*rho*P
       - ((1.0/3.0)*I)*(((2.0 - C1_)*epsilon_ - C1s_*G)*alpha*rho)
-      + (C2_*(alpha*rho*epsilon_))*dev(symm(b&b)) // symm should not be needed
+      + (C2_*(alpha*rho*epsilon_))*dev(innerSqr(b))
       + alpha*rho*k_
        *(
-            (C3_ - C3s_*mag(b))*S
+            (C3_ - C3s_*mag(b))*dev(S)
           + C4_*dev(twoSymm(b&S))
           + C5_*twoSymm(b&Omega)
         )
diff --git a/src/TurbulenceModels/turbulenceModels/RAS/SSG/SSG.H b/src/TurbulenceModels/turbulenceModels/RAS/SSG/SSG.H
index 85e08eacc..a0977fba6 100644
--- a/src/TurbulenceModels/turbulenceModels/RAS/SSG/SSG.H
+++ b/src/TurbulenceModels/turbulenceModels/RAS/SSG/SSG.H
@@ -186,6 +186,12 @@ public:
             return epsilon_;
         }
 
+        //- Return the effective diffusivity for R
+        tmp<volSymmTensorField> DREff() const;
+
+        //- Return the effective diffusivity for epsilon
+        tmp<volSymmTensorField> DepsilonEff() const;
+
         //- Solve the turbulence equations and correct eddy-Viscosity and
         //  related properties
         virtual void correct();