From 48844d3f12ec8e567e00848ecd4b7abdf78304a8 Mon Sep 17 00:00:00 2001
From: ignis <ignis@postech.ac.kr>
Date: Wed, 22 Apr 2020 13:58:28 +0900
Subject: [PATCH] turbulent diffusivity simply added

---
 .../cantera/transport/TurbulentTransport.h    | 53 ++++++++-----------
 src/transport/TurbulentTransport.cpp          | 10 +++-
 2 files changed, 31 insertions(+), 32 deletions(-)

diff --git a/include/cantera/transport/TurbulentTransport.h b/include/cantera/transport/TurbulentTransport.h
index 69437875c..c3eaa4c5d 100644
--- a/include/cantera/transport/TurbulentTransport.h
+++ b/include/cantera/transport/TurbulentTransport.h
@@ -32,6 +32,27 @@ public:
 
     virtual void init(thermo_t* thermo, int mode=0, int log_level=0);
 
+    //! Returns the mixture thermal conductivity (W/m /K)
+    /*!
+     * The thermal conductivity is computed from the following mixture rule:
+     * \f[
+     *     \lambda = 0.5 \left( \sum_k X_k \lambda_k  + \frac{1}{\sum_k X_k/\lambda_k} \right)
+     * \f]
+     *
+     * It's used to compute the flux of energy due to a thermal gradient
+     *
+     * \f[
+     *     j_T =  - \lambda  \nabla T
+     * \f]
+     *
+     * The flux of energy has units of energy (kg m2 /s2) per second per area.
+     *
+     * The units of lambda are W / m K which is equivalent to kg m / s^3 K.
+     *
+     * @returns the mixture thermal conductivity, with units of W/m/K
+     */
+    virtual doublereal thermalConductivity();
+
     //! Returns the unity Lewis number approximation based diffusion
     //! coefficients [m^2/s].
     /*!
@@ -57,22 +78,8 @@ public:
     virtual void getMixDiffCoeffs(double* const d) {
         GasTransport::getMixDiffCoeffs(d);
 
-        doublereal t = m_thermo->temperature();
-        doublereal slope = (-(m_turbmodifier-1.0)/m_profwidth);
-        doublereal trans = slope*(t-m_rdlayerstart) + 1.0;
-
-        double Dm = thermalConductivity() / (m_thermo->density() * m_thermo->cp_mass());
-
-        if (t < m_rdlayerstart - m_profwidth) {
-            Dm *= m_turbmodifier;
-        } else if (t > m_rdlayerstart) {
-            Dm *= 1.0;
-        } else {
-            Dm *= trans;
-        }
-
         for (size_t k = 0; k < m_nsp; k++) {
-            d[k] = Dm;
+            d[k] += m_turbmodifier;
         }
     }
 
@@ -97,22 +104,8 @@ public:
     virtual void getMixDiffCoeffsMass(double* const d){
         GasTransport::getMixDiffCoeffsMass(d);
 
-        doublereal t = m_thermo->temperature();
-        doublereal slope = (-(m_turbmodifier-1.0)/m_profwidth);
-        doublereal trans = slope*(t-m_rdlayerstart) + 1.0;
-
-        double Dm = thermalConductivity() / (m_thermo->density() * m_thermo->cp_mass());
-
-        if (t < m_rdlayerstart - m_profwidth) {
-            Dm *= m_turbmodifier;
-        } else if (t > m_rdlayerstart) {
-            Dm *= 1.0;
-        } else {
-            Dm *= trans;
-        }
-
         for (size_t k = 0; k < m_nsp; k++) {
-            d[k] = Dm;
+            d[k] += m_turbmodifier;
         }
     }
 
diff --git a/src/transport/TurbulentTransport.cpp b/src/transport/TurbulentTransport.cpp
index 7b68925bc..1b066a891 100644
--- a/src/transport/TurbulentTransport.cpp
+++ b/src/transport/TurbulentTransport.cpp
@@ -29,12 +29,18 @@ void TurbulentTransport::init(ThermoPhase* thermo, int mode, int log_level)
     myfile.open ("turbulent-coefs.txt");
 
     myfile >> m_turbmodifier;
-    myfile >> m_profwidth;
-    myfile >> m_rdlayerstart;
 
     myfile.close();
 
     // m_reftemp = thermo->temperature();
     // cout << "check reference temperature" << m_reftemp << endl;
 }
+
+doublereal TurbulentTransport::thermalConductivity()
+{
+    MixTransport::thermalConductivity();
+    m_lambda += m_turbmodifier;
+    m_condmix_ok = true;
+    return m_lambda;
+}
 }