## Mathematical Model
### Eddy Viscosity and Heat Diffusivity Model  
In order to model the eddy viscosity \(\
u_{t}\) and heat diffusivity \(\alpha_{t}\), we need to introduce different time scales. For the dynamical case, let \(\tau_{lu}\) be the dynamical turbulence time scale and the eddy viscosity \(\
u_{t}\) be defined as [8; 16; 17; 18; 27]  
\[\
u_{t}=C_{\mu}k\tau_{lu}\,,\qquad\tau_{lu}=f_{1\mu}A_{1\mu}+f_{2\mu}A_{2\mu}\,, \tag{17}\]  
where  
\[f_{1\mu}A_{1\mu}=(1-e^{-R_{d}/14})^{2}\,\tau_{u}\,,\qquad f_{2\mu}A_{2\mu}=f_ {1\mu}e^{-2.5\times 10^{-5}R_{t}^{2}}\,\tau_{u}\frac{5}{R_{t}^{3/4}}\,. \tag{18}\]  
The far wall characteristic dynamical time scale \(\tau_{u}\) is defined as \(k/\varepsilon=1/C_{\mu}\epsilon^{\Omega}\), \(f_{1\mu}A_{1\mu}\) and \(f_{2\mu}A_{2\mu}\) are the correction terms for the bulk and the near wall behavior, respectively. As the wall distance increases, \(f_{2\mu}\) decreases, and it is significantly different from zero only in the near wall region. This function forces the correct near wall behavior of \(\
u_{t}\) when the wall distance \(\delta\) tends to zero; namely, \(\
u_{t}\) has to be proportional to \(\delta^{3}\). For the thermal case, let \(\tau_{\theta}\) be the dynamical turbulence time scale and the eddy heat diffusivity \(\alpha_{t}\) be defined by [8; 16; 17; 18; 19; 20; 27]  
\[\alpha_{t}=C_{\theta}k\tau_{\theta}\,,\qquad\tau_{\theta}=f_{1\theta}B_{1 \theta}+f_{2\theta}B_{2\theta}\,, \tag{19}\]  
where  
\[f_{1\theta}B_{1\theta} =\left(1-e^{-\sqrt{P_{t}}R_{d}/19}\right)\,\left(1-e^{-R_{d}/14} \right)\,\,\frac{\tau_{u}}{P_{t}\tau_{t,\infty}}\,, \tag{20}\] \[f_{2\theta}B_{2\theta} =\tau_{u}\left(f_{2\theta}\frac{2R}{C_{\gamma}+R}+f_{2\theta} \sqrt{\frac{2R}{P_{\gamma}}}\frac{1.3}{\sqrt{P_{t}}R_{t}^{3/4}}\right)\,,\] (21) \[f_{2\theta} =f_{1\theta}e^{-(R_{t}/500)^{2}}\,,\qquad f_{2\theta}=f_{1 \theta}e^{-(R_{t}/200)^{2}}\,. \tag{22}\]  
As in the dynamical case the time scale \(\tau_{l\theta}\) is calculated as the sum of a near and a far wall term. \(B_{1\theta}\) is the bulk term and depends only on the dynamical time scale \(\tau_{u}\). It is the most important term in the region far from the wall, while \(B_{2\theta}\) dominates in the near and medium regions. In the near and medium region, another characteristic thermal time scale should be introduced. Let \(\tau_{\theta}\) be the thermal time scale defined by \(k_{\theta}/\varepsilon_{\theta}=1/C_{\mu}e^{\tau_{\theta}}\). Let \(R\) be the thermal to dynamical time scale ratio \(\tau_{\theta}/\tau_{u}=k_{\theta}\varepsilon/(k\varepsilon_{\theta})=e^{\Omega -\Omega_{\theta}}\). \(B_{2\theta}\), which consists of two contributions, and is the sum of a near wall term proportional to the time scale \(\tau_{u}\sqrt{R}=\sqrt{\tau_{u}\tau_{\theta}}\) and a term proportional to the mixed time scale \(\tau_{m}\) calculated as the harmonic mean between the time scales \(\tau_{u}\) and \(\tau_{\theta}\); namely, \(\tau_{m}=\tau_{u}R/(C_{\gamma}+R)\). This setting allows the model to reproduce the near wall behavior of the eddy thermal diffusivity with \(\alpha_{t}\) proportional to the cubic power of the distance from the wall \(\delta\) in the presence of no thermal fluctuations on wall surfaces, and \(\alpha_{t}\) proportional to the square power of the distance from the wall \(\delta\), when thermal fluctuations are considered [8; 16; 17; 27]. In this study we shall also consider the Kays model which has been used in many RANS simulations [5; 13]. This model computes \(\alpha_{t}\) as a function of the eddy viscosity \(\
u_{t}\). We consider \(\
u_{t}\) as in Equation (17) with the appropriate time scales and model functions already defined. \(\alpha_{t}\) is computed through the definition of a variable turbulent Prandtl number as  
\[Pr_{t}=0.85+\frac{0.7\,\
u}{Pr\,\
u_{t}}\,,\qquad\qquad\alpha_{t}=\frac{\
u_{ t}}{Pr_{t}}\,, \tag{23}\]  
where \(\
u_{t}/\
u\) is the eddy viscosity ratio. This model has some advantages because it is very simple, not computationally expensive and allows one to compute the \(Pr_{t}\) as a local function.