<style> .reveal { font-size: 22px; } </style> ## CFD crash course ### Recognizing the bias in flow modeling <br/> ### --- Session iii - The Turbulence challenge--- --- ## When flow becomes unstable, ... Content - Introduction to turbulent flows - Common turbulence phenomena - Strategies handling turbulent flows - Turbulence models in CFD 1. Reynolds Averaging Navier-Stokes (RANS) method and the modeling concepts - The Reynolds stress representations - Boussinesq assumption and the eddy viscosity methods - High-Reynolds and Low-Reynolds turbulence model - The wall treatment 2. Eddy resolving transient simulation strategies - Large Eddy Simulation (LES) and the grid impact - RANS-LES models for seamless approaches 3. Steady and Transient approach --- #### --- Introduction to turbulence flows --- ### What is turbulence? <img src="https://www.gauss-centre.eu/fileadmin/research_projects/2022/CSE/Oberlack_pn73fu/Oberlack_Fig_1_Merged.jpg" width="23%" style="transform:rotate(270deg)"> <img src="https://res.cloudinary.com/engineering-com/image/upload/w_640,h_640,c_limit/turbuence-car-ansys_gcygpo.gif" width="40%"> <img src="https://www.researchgate.net/publication/333773336/figure/fig1/AS:769533098856449@1560482396718/A-sketch-of-a-turbulent-pipe-flow.png" width="45%" > <br/> Chaotic, 3-dimensional, instantaneous, irregular ... <br/>flow movements with numerous eddy structures ---- #### --- Introduction to turbulent flows --- ### Every thing starts from the instability Normalize the Navier-Stokes equation with: $$ \begin{align} u^* _ i = \frac{u_i}{U}\quad x^*_i = \frac{x_i}{ L},\quad p^* = \frac{p}{\rho U^2},\quad t^* = \frac{t}{T} \\ \frac{d u^*_i}{d t^*} =\frac{\partial u^*_i}{\partial t^*} + u^*_j\frac{\partial u^*_i}{\partial x^*_j} =-\frac{\partial p^*}{\partial x^*_i} + \frac{\nu}{UL} \frac{\partial ^2 u_i^*}{\partial {x_j^*}^2} \\ \rightarrow\rightarrow\rightarrow\frac{\nu}{UL} = \frac{1}{\text{Re}} \end{align} $$ - Reynolds number distinguishes laminar from turbulent states <img src="https://qph.cf2.quoracdn.net/main-qimg-c2598941cf4df97e4aebcd212b6b0d3e-lq" style="float:right" width=50%> <div style="text-align: left "> - &emsp;&emsp;for Re < $\text{Re}_{\text{cri}}$,<br/> &emsp;&emsp;&emsp;&emsp;&rarr; flow stabalizes, - &emsp;&emsp;for Re < $\text{Re}_{\text{cri}}$,<br/> &emsp;&emsp;&emsp;&emsp;&rarr; flow tends to unstable states </div> ---- #### --- Introduction to turbulent flows --- ### Life cylce of eddies <img src="https://i.imgur.com/AwLOY7V.png" style="float:right" width="35%"> <div style="text-align: left "> Evolution of turbulent strcutures <br/> <br/> - Instability flow structure leads to generation of large eddies <br/> <br/> - Large eddies pass energy to middle size eddies, further to small eddies. Meanwhile. <br/> <br/> - Smallest eddies dissipate. Energy transfer in the form of heat. Viscous effect dominates. <br/> <br/> <br/> &rarr; A trubulent field obtains __*all scales of eddy sizes*__, represented in the energy spectrum &rarr; The more turbulent (higher Reynolds number), the longer the evolution track, leading to higher wave number ($\kappa = 1/\lambda$). </div> ---- #### --- Introduction to turbulent flows --- ### The challenge in the CFD framework How fine should be the mesh? <img src="https://i.imgur.com/dgqlTIA.png" width="70%"> If we attempt to capture the entire flow structure, i.e. - Resolving till the smallest lengthscale: $\large \eta = (\frac{\nu}{\epsilon})^{1/4}$ - The higher Reynolds number, the smaller the smallest lengthscale. - Number of cell (N) ~ $Re^{9/4}$ &rarr; Completely resovling flow structure is not affortable. --- #### --- Common turbulence phenomena --- ### Statistical operating point of view <img src="https://i.imgur.com/Bc2nrbY.png" style="float:right" width= "20%"> Probability Density Function (PDF): $\mathbf{f} (V; x, t)$ - Mean: $\small\overline{U}(x,t) = \int\limits^{\infty} _{-\infty} V \mathbf{f}(V;x,t) dV = \lim\limits_{N \to \infty} \frac{1}{N}\sum\limits_{\xi=1}^N V^{(\xi)}$ $\small\qquad\qquad\Leftrightarrow\qquad U_i = \overline{U_i} + u_i^{\prime}$ <br/> - Variance: $\small\text{var}(U)= \int\limits^{\infty} _{-\infty} (V-\overline{U})^2 \mathbf{f}(V;x,t) dV$ $\small\qquad\qquad\qquad= \lim\limits_{N \to \infty} \frac{1}{N}\sum\limits_{\xi=1}^N \left[(V^{(\xi)} - \overline{U})(V^{(\xi)} - \overline{U})\right] = \overline{u^{\prime}_iu^{\prime}_i}$ <br/> - Covariance: $$ \small \begin{align} \text{cov}(U_i,U_j) = &\int\limits^{\infty} _{-\infty}\int\limits^{\infty} _{-\infty} (V_i-\overline{U}_i)(V_j-\overline{U}_j) \mathbf{f}_{ij}(V_i,V_j;x,t) dV_idV_j\\ = &\lim_{N \to \infty} \frac{1}{N}\sum\limits_{\xi=1}^N \left[(V_i^{(\xi)} - \overline{U}_i)(V_j^{(\xi)} - \overline{U}_j)\right] = \overline{u^{\prime}_iu^{\prime}_j} \end{align} $$ ---- #### --- Introduction to turbulent flows --- ### Statistical operating point of view - Correlated behavior $\mathbb{R}_{ij}(x,r_k,t, \Delta t) = \overline{u^\prime_i(x, t)u^\prime_j(x+r_k,t+\Delta t)}$ - In homogeneous, isotropic state and at the same time: - $\Delta t = 0$, $\quad j = i$ <img src="https://i.imgur.com/qfpo1Xc.png" width="80%"> $\large\kappa = \frac{1}{l}\qquad$ &rarr; the larger the lengthscale, the smaller the wave number ---- #### --- Common turbulence phenomena --- ### Energy spectrum and the decaying process [Decay of homogeneous isotropic turbulence](https://vimeo.com/307065762) <img src="https://www.researchgate.net/publication/243404571/figure/fig3/AS:669311685361675@1536587749335/Decay-of-turbulent-kinetic-energy-DNS-and-LES-results-at-low-Re.png" width="40%"> - Decay rate of HIT is time dependent $\frac{\partial \overline{u_i^\prime u_i^\prime}(t)}{\partial t} = a\ \overline{u_i^\prime u_i^\prime}(t_0) \ (t-t_0)^{-n}$ $\large n= 1.2 \sim 1.7$ ---- #### --- Common turbulence phenomena --- ### Velocity profile on turbulent pipe flow <img src="https://i.imgur.com/zjXbipL.png" width="50%"> Thinking on: - The velocity profile for laminar flow (green line) yield $U(r) = U_{max}(1- (r/R)^2)$ as solution of Navier-Stokes Equation. - What is the equation, for which the mean velocity profile can be solved? ---- #### --- Common turbulence phenomena --- ### Turbulent boundary layer and Velocity profile No analytic solution for flow profile based on slip/non-slip condition anymore. <img src="https://www.mdpi.com/fluids/fluids-06-00369/article_deploy/html/images/fluids-06-00369-g002.png" width="30%"> <img src="https://upload.wikimedia.org/wikipedia/commons/thumb/4/4c/Law_of_the_wall_%28English%29.svg/500px-Law_of_the_wall_%28English%29.svg.png" width="30%"> $y^+ = y \ u_\tau /\nu\ ;\qquad$ $u_\tau = \sqrt{\frac{\tau_w}{\rho} }\ ;\qquad$ $\tau_w = \mu \frac{\partial \overline{U}}{\partial y}\ ;$ $\qquad u^+ = \frac{\overline{U}}{u_\tau}$ Three zones in near wall region categorized: - viscous sublayer: $y^+ < 5\ :\qquad$ $\qquad u^+ = y^+$ - buffer layer: $5< y^+ <13$ - __*log-law*__ region: $y^+ > 30\ : \qquad$ $\qquad u^+ = \frac{1}{\kappa} \text{ln}\ y^+ + C^+$ --- #### --- Strategies handling turbulent flows --- ### Reynolds Averaging Navier-Stokes method <img src="https://i.imgur.com/hlI0Vp2.png" width="70%"> The averaging principles: - $\phi_i = \overline{\phi}_i \ + \phi_i^{\prime}$ - $\overline{\overline{\phi}_i} = \overline{\phi}_i \quad; \quad \overline{\phi_i^{\prime}} =0$ Thinking on: - How about __$\overline{\phi_i \phi_j}$, $\overline{\overline{\phi}_i \overline{\phi}_j}$ and $\overline{\overline{\phi}_i \phi_j^\prime}$ ?__ - Is $\overline{\phi_i \phi_j} = \overline{\overline{\phi}_i \overline{\phi}_j}$? - What are the physical properties in the momentum transport to be averaged? ---- #### --- Strategies handling turbulent flows --- ### Reynolds Averaged Navier-Stokes (RANS) method Let $\quad U_i = \overline{U}_i+ u^\prime_i\qquad P = \overline{P}+P^\prime$ <br/> Averaging the Navier Stokes Equation in incompressible flow yields: $$ \begin{align} \overline{\frac{\partial U_j}{\partial x_j}} = \frac{\partial \overline{U}_j}{\partial x_j} & = 0\ \\ \overline{\frac{D \ U_i}{Dt}} = \overline{\frac{\partial U_i}{\partial t}} + \overline{\frac{\partial}{\partial x_j} U_i U_j }& = \overline{-\frac{1}{\rho}\frac{\partial P}{\partial x_i}} + \overline{\frac{\partial}{\partial x_j}\left( \nu \frac{\partial U_i}{\partial x_j}\right)} \\ = \frac{\partial}{\partial x_j} ( \overline{U}_i \overline{U}_j & + \overline{u_i^\prime u_j^\prime} ) = -\frac{1}{\rho}\frac{\partial \overline{P}}{\partial x_i} + \frac{\partial}{\partial x_j}\left( \nu \frac{\partial \overline{U}_i}{\partial x_j}\right) \end{align} $$ Thinking on: - Where is the time derivative term ($\Large\frac{\partial \overline{U}_i}{\partial t}$)? - Are the equations ready to solve? - Is there anything unknown in the momentum transport? ---- #### --- Strategies handling turbulent flows --- ### Reynolds stress in mathematical intepretation __*Reynolds stress ($\overline{u^\prime_i u^\prime_j }$)*__ --- the stumbling block in the transport equation - The covariance __*at the same point*__ with varying directions. $i,j = 1, 2, 3$ - $\overline{u^\prime_i u^\prime_j } = \lim_{N \to \infty} \frac{1}{N}\sum\limits_{\xi=1}^N \left[(V_i^{(\xi)} - \overline{U}_i)(V_j^{(\xi)} - \overline{U}_j)\right]$ - A tensor 2nd order, symmetric, totally 6 elements - $\overline{u_i^\prime u_j^\prime} = \begin{pmatrix} \overline{u^\prime u^\prime} & \overline{u^\prime v^\prime} & \overline{u^\prime w^\prime} \\ \overline{v^\prime u^\prime} & \overline{v^\prime v^\prime} & \overline{v^\prime w^\prime }\\ \overline{w^\prime u^\prime} & \overline{w^\prime v^\prime} & \overline{w^\prime w^\prime} \end{pmatrix}$ - Trace of Reynolds stress $\left(tr(\overline{u^\prime_i u^\prime_j })\right)$ = 2 xTurbulent kinetic energy ($k$) - $k = \frac{1}{2} (\overline{u^\prime u^\prime}+ \overline{v^\prime v^\prime} + \overline{w^\prime w^\prime}) = \frac{1}{2} \overline{u_i^\prime u_i^\prime}$ ---- #### --- Strategies handling turbulent flows --- ### Movement of continuum and the Reynolds stress tensor <img src="https://i.imgur.com/U02jTBa.png" style="float:right" width ="25%"> Movement of continuum are described in - $\large \overline{S}_{ij} = \frac{1}{2}(\frac{\partial \overline{U}_i}{\partial x_j}+\frac{\partial \overline{U}_j}{\partial x_i})$ : rate of strain (deformation) - $\large \overline{\Omega}_{ij} = \frac{1}{2}(\frac{\partial \overline{U}_i}{\partial x_j}-\frac{\partial \overline{U}_j}{\partial x_i})$ : rate of rotation <br/> <br/> Cayley-Hamilton's Theorem for every tensor 2nd order: $A_{ij} ^3 - \sigma_1 A_{ij} ^2 + \sigma_2 A_{ij} -\sigma_3 I_{ij} =0$ $\sigma_1$, $\sigma_2$ and $\sigma_3$ are the invariance of the $A_{ij}$. - higher order tensor can be downgraded into combination of $A_{ij}^3$, $A_{ij}^2$, $A_{ij}$ and $I_{ij}$ ---- #### --- Strategies handling turbulent flows --- ### Exact solution via flow pattern description Applying rate of strain and rotation to describe Reynolds stress, $\overline{u_i^\prime u_j^\prime} = \mathfrak{F}\left(k, \overline{S}_{ij}, \overline{S}_{ij}^2, \overline{S}_{ij}^3, \overline{\Omega}_{ij}, \overline{\Omega}_{ij}^2, \overline{\Omega}_{ij}^3,.(\text{their combinations and invariants}).. \right)$ $$ \begin{align} \frac{\overline{u_i^\prime u_j^\prime}}{k}-\frac{2}{3}\delta_{ij} = & \beta_1 \overline{S}_{ij}\\ +& \beta_2 (\overline{S}_{ij}^2 -\frac{1}{3}\text{II}_S \mathbb{I})+\beta_3 (\overline{\Omega}_{ij}^2 -\frac{1}{3}\text{II}_{\Omega} \mathbb{I}) + \beta_4 (\overline{S}_{ij}\overline{\Omega}_{ij}-\overline{\Omega}_{ij}\overline{S}_{ij}) \\ +&\beta_5 ( \overline{S}_{ij} ^2\overline{\Omega}_{ij} - \overline{\Omega}_{ij} ^2\overline{S}_{ij}) + \beta_6 ( \overline{S}_{ij} ^2\overline{\Omega}_{ij} + \overline{\Omega}_{ij} ^2\overline{S}_{ij} - \frac{2}{3} \text{IV} \mathbb{I})\\ +&\beta_7 ( \overline{S}_{ij} ^2\overline{\Omega}_{ij} ^2 + \overline{\Omega}_{ij} ^2\overline{S}_{ij} ^2 - \frac{2}{3} V \mathbb{I})+ \beta_8 (\overline{S}_{ij} \overline{\Omega}_{ij} \overline{S}_{ij} ^2 - \overline{S}_{ij} ^2\overline{\Omega}_{ij} \overline{S}_{ij})\\ +&\beta_9 ( \overline{\Omega}_{ij} \overline{S}_{ij} \overline{\Omega}_{ij} ^2 - \overline{\Omega}_{ij} ^2\overline{S}_{ij} \overline{\Omega}_{ij}) + \beta_{10} (\overline{\Omega}_{ij} \overline{S}_{ij} ^2 \overline{\Omega}_{ij} ^2 - \overline{\Omega}_{ij} ^2 \overline{S}_{ij} ^2 \overline{\Omega}_{ij}) \end{align} $$ $\small\text{II}_{S} = tr\{\overline{S}_{ij}^2\},\; II_{\Omega} = tr\{\overline{\Omega}_{ij}^2\},\;\text{III}_{S} = tr\{\overline{S}_{ij}^3\},\;\text{IV}=tr\{{\overline{S}_{ij}} \overline{\Omega}^2_{ij}\},\; \text{V} = tr\{\overline{\Omega}_{ij}^2\overline{S}_{ij}^2\}$ - Which terms are known, and which terms are unknown? - $k$, $\overline{S}_{ij}$, $\overline{\Omega}_{ij}$, $\beta_{1-10}$ ---- #### --- Strategies handling turbulent flows --- ### Exact transport of Reynolds stress tensor $$ \begin{align} \frac{D\ \overline{u^\prime_i u^\prime_j}}{Dt} = \frac{\partial\ \overline{u^\prime_i u^\prime_j}}{\partial t} + \underbrace{\overline{U}_k \frac{\partial \overline{u^\prime_i u^\prime_j}}{\partial x_k}}_{\text{Convection}} &= \underbrace{- \left(\overline{u^\prime_i u^\prime_k}\frac{\partial \overline{U}_j}{\partial x_k} + \overline{u^\prime_j u^\prime_k}\frac{\partial \overline{U}_i}{\partial x_k} \right)}_{\text{Production}} - \underbrace{2\nu\overline{\frac{\partial u^\prime_i}{\partial x_k}\frac{\partial u^\prime_j}{\partial x_k}}}_{\text{Dissipation}}\\ +\underbrace{\overline{\frac{p^\prime}{\rho} \left( \frac{\partial u_i^\prime}{\partial x_j}+\frac{\partial u_j^\prime}{\partial x_i}\right)}}_{\text{Pressure-Strain correlation}} &+\underbrace{\frac{\partial}{\partial x_k}\left(\nu\frac{\partial \overline{u^\prime_i u^\prime_j}}{\partial x_k} - \overline{u^\prime_i u^\prime_j u^\prime_k} - \overline{\frac{p^\prime}{\rho}\left( u_i^\prime\delta_{jk} + u_j^\prime\delta_{ik}\right)}\right)}_{\text{Diffusion}} \end{align} $$ <br/> - Which terms are known, and which terms are unknown? - $\overline{u^\prime_i u^\prime_j}$, $\quad\overline{u^\prime_i u^\prime_k}\frac{\partial \overline{U}_j}{\partial x_k}$, $\quad 2\nu\overline{\frac{\partial u^\prime_i}{\partial x_k}\frac{\partial u^\prime_j}{\partial x_k}}$, $\quad\overline{\frac{p^\prime}{\rho} \left( \frac{\partial u_i^\prime}{\partial x_j}+\frac{\partial u_j^\prime}{\partial x_i}\right)}$ - $\nu\frac{\partial \overline{u^\prime_i u^\prime_j}}{\partial x_k}$, $\quad\overline{u^\prime_i u^\prime_j u^\prime_k}$, $\quad\overline{\frac{p^\prime}{\rho}\left( u_i^\prime\delta_{jk} + u_j^\prime\delta_{ik}\right)}$ --- #### --- Turbulence modeling in CFD --- ### RANS (Reynolds Averaged Navier-Stokes) Models <img src="https://i.imgur.com/nSYiV1d.png" width="65%"> <br/> - Showing averaged physical characteristics ($\overline{\Phi}$) - Smooth representation in all the properties - No turbulent strtures revealed - Mostly applied in all turbulence relevant engineering assignements - Aerodynamics, Heat transfer, multiphase - Computational cost in control - convergence can sometimes be challenging ---- #### --- Turbulence modeling in CFD --- ### RANS: RSM ---The Reynolds Stress Models - The unkown terms in Reynolds stress transport are - $\quad 2\nu\overline{\frac{\partial u^\prime_i}{\partial x_k}\frac{\partial u^\prime_j}{\partial x_k}}$, $\quad\overline{\frac{p^\prime}{\rho} \left( \frac{\partial u_i^\prime}{\partial x_j}+\frac{\partial u_j^\prime}{\partial x_i}\right)}$, $\quad\overline{u^\prime_i u^\prime_j u^\prime_k}$, $\quad\overline{\frac{p^\prime}{\rho}\left( u_i^\prime\delta_{jk} + u_j^\prime\delta_{ik}\right)}$ - Solving the transport of Reynolds stress $(\large\frac{D \overline{u^\prime_i u^\prime_j}}{Dt})$, 6 elements with additional transport of dissipation $\large(\frac{D\epsilon}{Dt})$ - Capable in capturing __*rotation, curvature, anisotropic turbulent characteristics.*__ - Models (e.g. LRR, LG, SSG, ...) varies only from different modeling principles of unknown terms. - Stable numerics in complex system not easy to achieve. ---- #### --- Turbulence modeling in CFD --- ### RANS: EVM ---The Eddy Viscosity Models <div style="text-align: left "> $\quad \qquad\large\frac{\overline{u_i^\prime u_j^\prime}}{k}-\frac{2}{3}\delta_{ij} = \beta_1 \overline{S}_{ij}$ </div> $+ \beta_2 (\overline{S}_{ij}^2 -\frac{1}{3}\text{II}_S \mathbb{I})+\beta_3 (\overline{\Omega}_{ij}^2 -\frac{1}{3}\text{II}_{\Omega} \mathbb{I}) + \beta_4 (\overline{S}_{ij}\overline{\Omega}_{ij}-\overline{\Omega}_{ij}\overline{S}_{ij}) + ...$ <br/> __Boussinesq hypothesis__: <br/>Taking the __first order__ flow movement as dependency for approximation - $\large\overline{u_i^\prime u_j^\prime}=k\beta_1 \overline{S}_{ij}+\frac{2}{3}k\delta_{ij};\qquad k\beta_1 = -2\nu_t$ - $k$ and $\beta_1$ are the unknown variables, which require __additional transports__ to approach. Thinking on: - What is the dimension of $\beta_1$, $k$ and turbulent viscosity $\nu_t$? - Why is there no $\overline{\Omega}_{ij}$ contribution? ---- #### --- Turbulence modeling in CFD --- ### RANS: EVM --- A dimensional approach of Eddy viscosity $\large \overline{u_i^\prime u_j^\prime}=-2\nu_t \overline{S}_{ij}+\frac{2}{3}k\delta_{ij};$ <br/> $\large \nu_t \ [m^2/s] = (L\ [m])^2/(\tau\ [s])$ $L: \text{trubulent length scale;}\quad \tau: \text{turbulent time scale}$ Transport equations ($\frac{Dk}{Dt}$ and $\frac{D\Phi}{Dt}$ ) utilized to acquired the implemented properties <br/> Thinking on: - in the dependency of turbulent kinetic energy $k$ and dissipation $\epsilon$, <br/> what is the representation of eddy viscosity $\nu_t$? - Introducing an attribute called frequency $\omega \ [1/s]$, <br/> how is the representation of eddy viscosity $\nu_t$? ---- #### --- Turbulence modeling in CFD --- ### RANS: EVM --- the k- somthing high-Re eddy viscosity models $\overline{u_i^\prime u_j^\prime}=-2\nu_t \overline{S}_{ij}+\frac{2}{3}k\delta_{ij};$ - Transport of $k$ is always necessary <br/> | Model | 2nd transport variable ($\Phi$)| dimension | $\nu_t$| | -------- | :--------: | -------- | -------- | | $k-\epsilon$ | $\epsilon$ | $[m^2/s^3]$ | $C_\mu k^2/\epsilon$ | | $k-\omega$ | $\omega$ | $[1/s]$ | $C_\mu k/\omega$ | | $k-L$ | $L$ | $[m]$ | $C_\mu k^{1/2}L$ | | $k-kL$ | $kL$ | $[m^3/s^2]$ | $C_\mu kL/k^{1/2}$ | <br/> Features of all the eddy viscosity models: - The turbulent anisotropy $(\frac{\overline{u^\prime_i u^\prime_j}}{k}-\frac{2}{3}\delta_{ij})$ imitates the flow shear stress ($\overline{S}_{ij}$) - Generally not capable capturing curved flow type - Since no rate of rotation ($\overline{\Omega}_{ij}$) contribution obtained - The wall information of High-Reynolds models goes only till log-law region. - No detail information beyond buffer region available ---- #### --- Turbulence modeling in CFD --- ### RANS: The $k$-$\epsilon$ model example $$ \begin{align} \frac{\partial ku_j}{\partial x_j} &= P_k -\epsilon + \frac{\partial}{\partial x_j}\left[ \frac{\nu_t}{\sigma_k} \frac{\partial k}{\partial x_j}\right] \\ \frac{\partial \epsilon u_j}{\partial x_j} &= c_{\epsilon 1} \frac{\epsilon}{k}P_k -c_{\epsilon 2}\frac{\epsilon^2}{k} + \frac{\partial}{\partial x_j}\left[ \frac{\nu_t}{\sigma_k} \frac{\partial \epsilon}{\partial x_j}\right] \\ \\ & -\overline{u^\prime_i u^\prime_j} = 2\nu_t \overline{S}_{ij}-\frac{2}{3}k\delta_{ij} \qquad \nu_t = c_\mu \frac{k^2}{\epsilon} \end{align} $$ <img src="https://i.imgur.com/NeIWghm.png" width ="50%" style="float:right"> $c_\mu = 0.09$ : calibrated via shear stress $\overline{u^\prime_i u^\prime_j}$ to $k$ behavior. $c_{\epsilon 1} = 1.44$ : calibrated by equilibrium of $P_k$ and $\epsilon$ in log region $c_{\epsilon 2} = 1.92$ : calibrated via decay of homogeneous isotropic turbulence - No near wall characteristics obtained --> high- Reynolds description ---- #### --- Turbulence modeling in CFD --- ### RANS: yplus and the Wall functions $\large y^+ = y u_\tau /\nu$ <img src="https://i.imgur.com/3FL39zJ.png" width="40%" > - All the High-Re model should be equipped with wall function as boundary condition - $y^+$ of the first cell to wall determines the wall bounded flow type - Standard wall function modify the $\nu_t$, $P_k$ and $\epsilon$ according to the physics in the log-law region. i.e. - finer mesh ($y^+ < 13$) with high-Re models will show counter-effects <img src="https://i.imgur.com/mpBF4uG.png" width="50%"> ---- #### --- Turbulence modeling in CFD --- ### RANS: High and Low-Reynolds models Low-Reynolds models: addtional information for neal wall approach - Models are capable to integrate till the wall, no wall function necessary standard high-Re k-epsilon model: $$ \begin{align} \frac{\partial \epsilon u_j}{\partial x_j} = &c_{\epsilon 1} \frac{\epsilon}{k}P_k -c_{\epsilon 2}\frac{\epsilon^2}{k} + \frac{\partial}{\partial x_j}\left[ \frac{\nu_t}{\sigma_k} \frac{\partial \epsilon}{\partial x_j}\right]\\ & \nu_t = C_\mu k^2/\epsilon \end{align} $$ [Low-Re k-epsilon models](https://www.cfd-online.com/Wiki/Low-Re_k-epsilon_models): $$ \begin{align} \frac{\partial \epsilon u_j}{\partial x_j} = &c_{\epsilon 1} f_1\frac{\epsilon}{k}P_k -c_{\epsilon 2}f_2\frac{\epsilon^2}{k} + \frac{\partial}{\partial x_j}\left[ \frac{\nu_t}{\sigma_k} \frac{\partial \epsilon}{\partial x_j}\right]+ E\\ & \nu_t = C_\mu f_\mu k^2/\epsilon \end{align} $$ $f_\mu$, $f_1$, $f_2$ and $E$ are the __near wall ($y^+<13$) damping approach__, <br/>serve as __seamless transition__ to standard k-epsilon representation in outer region ---- #### --- Turbulence modeling in CFD --- ### RANS: High and Low-Reynolds models | | High-Re models | Low-Re models | | -------- | -------- | -------- | | Concepts | integrates till <br/>log-law region | integrates till <br/> viscous sublayer | | Handling | classical transport of $\Phi$ | additional terms regarding<br/> neall wall calibration taken | | Examples | k-Epsilon | k-Omega SST | | model related <br/> Mesh dependency | $y^+ > 15$ | no restriction | | Boundary condition <br/> for wall| standard wall function | viscous treatment or <br/> $y^+$ correlated hybrid treatment | ---- #### --- Turbulence modeling in CFD --- ### RANS: setups of inlet turbulence condition <img src="https://doc.cfd.direct/notes/cfd-general-principles/img/index2975x.png" width="58%" style="float:right"> [Turbulent indensity $I$ [%]](https://www.cfd-online.com/Wiki/Turbulence_intensity) $$ \large \begin{align} I\equiv & u^\prime /\overline{U} \\ u^\prime = &\sqrt{2k /3}\\ k = &\frac{3}{2}(I\ \overline{U})^2 \end{align} $$ [Turbulent length scale $L$](https://www.cfd-online.com/Wiki/Turbulence_length_scale) $$ \large \begin{align} L &\equiv C_\mu^{3/4} k^{3/2}/\epsilon \\ &\equiv C_\mu^{-1/4} \sqrt{k}/\omega \quad C_\mu = 0.09 \\ \epsilon &= C_\mu^{3/4} k^{3/2}/L \ ;\quad \omega = C_\mu^{-1/4} \sqrt{k}/L \end{align} $$ <div style="text-align: left "> - $L$ can be selected as dimension of the domain/characteristic object For free stream inflow hitting on immersed body, effects of inlet turbulence properties are generally insignificant. </div> ---- #### --- Turbulence modeling in CFD --- ### RANS: setups of inlet turbulence boundary layer <img src="https://i.imgur.com/nrrlWR6.jpg" width= "60%"> For __pipe/duct/boundary layer__ type flows targeting on boundary layer <br/> on immersed body, __consistent characteristics__ are important. <img src="https://i.imgur.com/5qacjsa.png" width="60%"> - using idea of precursor run to generate a fully developed turbulent profile based on the applied model. - Boundary layer flow profile is determined by applied transport equations, but not a reckless log-law curve ---- #### --- Turbulence models in CFD --- ### RANS: Energy equation related handling Transport equations for compressible flows $$ \begin{align} \frac{\partial \overline{\rho}}{\partial t} + \frac{\partial \overline{\rho}\overline{U}_j}{\partial x_j} & = 0\ \\ \frac{\partial \overline{\rho}\overline{U}_i}{\partial t} + \frac{\partial}{\partial x_j} ( \overline{\rho}\overline{U}_i \overline{U}_j + \overline{\rho}\overline{u_i^\prime u_j^\prime} ) &= -\frac{\partial \overline{P}}{\partial x_i} + \frac{\partial}{\partial x_j}\overline{\tau}_{ij} + \overline{\rho} g_i\\ \frac{\partial K}{\partial t} + \frac{\partial}{\partial x_j} \left(K \overline{U}_j + \overline{u^\prime_i h^\prime}\right) &=-\overline{P}\frac{\partial \overline{U}_i}{\partial x_i} -\overline{\tau}_{ij}\frac{\partial \overline{U}_j}{\partial x_i} + \frac{\partial}{\partial x_j}\left(\frac{c_p\mu}{Pr} \frac{\partial \overline{T}}{\partial x_j} \right)\\ K = \frac{1}{2}\overline{\rho} \overline{U}_i^2 + \overline{\rho} \overline{h} \quad&\qquad \overline{\tau}_{ij} = \mu\left[\left(\frac{\partial \overline{U}_j}{\partial x_i} + \frac{\partial \overline{U}_i}{\partial x_j}\right)-\frac{2}{3}\frac{\partial \overline{U}_k}{\partial x_k}\delta_{ij}\right] \end{align} $$ - Prandtl number: $Pr = \large\frac{c_p \mu}{\kappa}$ behavior of viscous $\mu$ to thermal diffusion $\kappa$. - Turbulence contribution of mechanical aspect in energy equation neglected - $\overline{u^\prime_i u^\prime_j}$ and $\overline{u^\prime_i h^\prime_j}$ are the unknowns in the main transport equations - $\overline{u^\prime_i h^\prime_j} = c_p\overline{u^\prime_i T^\prime_j}$, $\quad c_p$: specific heat capacity. - $\overline{u^\prime_i T^\prime_j}$ : __turbulent heat fluxs__ ---- #### --- Turbulence modeling in CFD --- ### RANS: Handling of turbulent heat flux -- $\overline{u^\prime_i T^\prime_j}$ Analogy to the Reynolds stress handling with the introduction of <br/>___$Pr_t$ --- turbulent Prandtl number___<br/> $\rho \overline{u_i^\prime T^\prime} = \Large\frac{\rho\nu_t}{Pr_t}\frac{\partial T}{\partial x_j}$ - Linear behavior to temperature gradient as the Bousinesq hypothesis for Reynolds stress - $Pr_t$ used to set as a constant value $\sim 0.9$, as single constant for turbulent heat flux approach - $\mu_t$ followed by Eddy Viscosity Modeling approach, e.g. $\mu_t = C_\mu k^2/\epsilon$ - applying the transport of turbulent characterisitcs ($k$, $\epsilon$, $\omega$) for turbulent heat flux representation ---- #### --- Turbulence modeling in CFD --- ### RANS: Wallfunction for turbulent heat flux Utilizing the same strategy in thermal wall treatment as to the momentum wall boundary conditions || [__Thermal Boundary Layer__](https://doc.cfd.direct/notes/cfd-general-principles/thermal-boundary-layers) | | Velocity Boundary Layer | |-------| :--------: | -------- | :--------: | || <img src="https://doc.cfd.direct/notes/cfd-general-principles/img/index3351x.png" style="float:right"> | | <img src="https://www.cfd-online.com/W/images/6/6a/Img_lawOfTheWall_whiteBG.png" width="50%"> | ||$T_\tau =- \Large\frac{q_w}{\kappa}\frac{\nu}{\text{Pr}_t}\frac{1}{u_\tau}$||$u_\tau = \sqrt{\tau_w/\rho}$| |$y^+ = y u_\tau/\nu$|$T^+ = (\overline{T}-T_w)/{T_\tau}$||$u^+ = \overline{U}/u_\tau$| |$y^+<5$|$T^+ = \text{Pr}_t\ y^+$||$u^+ = y^+$| |$y^+>30$|$T^+ = (\text{Pr}_t/\kappa)\ \text{ln}(y^+)+ B_t$||$u^+ =(1/\kappa) \ \text{ln}(y^+)+B$| |wall function updates|$\alpha_t$ $\quad(\kappa_t)$||$\nu_t$, $P_k$ and $\epsilon$| ---- #### --- Turbulence modeling in CFD --- ### RANS: Summary - I - RANS: Raynolds Averaged Navier-Stokes - physical properties are averaged and sorted to __$\overline{\text{mean}}$__ and __$\text{fluctuation}^\prime$__ parts. - Averaging transport equation system (continuety, momentum and energy), the unknown terms remain to Reynolds stress $(\overline{u_i^\prime u_j^\prime})$ and turbulent heat fluxs $(\overline{u_i^\prime T^\prime})$. - Averaging concept in CFD can be barely achieved in time, but not space/phase. - Modeling of the Reynolds stress classified to - RSM: direct transport of Reynolds stress $\frac{D\ \overline{u_i^\prime u_j^\prime} }{Dt}$ - EVM: linear behavior of turbulent anisotropy to main shear stress $\overline{S}_{ij}$ - and higher order $\overline{S}_{ij}$, $\overline{\Omega}_{ij}$ combinations - All the type of turbulence models can be High-Re and Low-Re, depending on the near wall treatment obtained in the model - Smooth profiles of physical properties $(\overline{\Phi})$ expected as simulation results ---- #### --- Turbulence modeling in CFD --- ### RANS: Summary - II - Boundary conditions of applied turbulence model depend on the modeling (high-Re/low-Re) strategy - $y^+$ is the characteristics to identify the mesh resolution on the first cell. - For high-Re models, wall function is necessary. The first to-wall cell size is restricted to be $y^+ > 13$ - (hybrid/commerically standard) wall function adpates the $\nu_t$, $\alpha_t$, $P_k$, $\epsilon$, ... on the fisrt cell - based on the evaluated $y^+$ value. - switch automatically between viscous sublayer and log-law region - RANS is generally recommanded only to apply in __steady state simulations__. - Generally poor performance in curved flow of EVMs, due to the lacking of rotational information - Delayed response to separated flow. Validation needed. --- #### --- Turbulence modeling in CFD --- ### Eddy resolving simulations (LES) <img src="https://www.cttc.upc.edu/sites/default/files/public/windenergy/airfoll/1.jpg" width="50%"> - Resolving phyiscal properties in scales instead of averaging - multi-scale structures show up, i.e. - diverse size of 3D eddies as velocity structures - unsteady effects even in statistically steady states. - Higher accuracy achieved by __resolving flow structure__ instead of completely modeling. - Reduction of the modeling portion lead to increase of resolving ---- #### --- Turbulence modeling in CFD --- ### LES: Filtering instead of averaging Large eddies captured, small get through the filter. <img src="https://i.imgur.com/a7UDxXw.png" width="50%"> $\widetilde{\Phi} (x_i,t) = \int \Phi(x_i, t) G(r,x_i,t) dr$ - The filtering process $G(r, x_i, t)$ is steered by a defined length scale $r$. - It can be $k^{3/2}/\epsilon$, $\sqrt{k/\omega^2}$ or grid size $\Delta$. - The filtering departs physical properties into: $\Phi(x_i, t) = \widetilde{\Phi} (x_i,t) + \phi^*$ - $\widetilde{\widetilde{\Phi}} (x_i,t) = \widetilde{\Phi} (x_i,t)$, $\quad\widetilde{\phi^*} = \phi^*$ ---- #### --- Turbulence modeling in CFD --- ### Eddy resolving simulations The filtered governing equations in incompressible flows, newtonian fluid: $$ \begin{align} & \frac{\partial \widetilde{U_i}}{\partial x_i} = 0\\ \frac{\partial \widetilde{U_i}}{\partial t} + \frac{\partial}{\partial x_j}& \left( \widetilde{U_i} \widetilde{U_j}+ \widetilde{\tau}_{ij}\right) = -\frac{1}{\rho}\frac{\partial \widetilde{P}}{\partial x_i}+\frac{\partial}{\partial x_j}\left(\nu\frac{\partial \widetilde{U_i}}{\partial x_j}\right)\\ \\ \widetilde{\tau}_{ij} = \widetilde{U_iU_j} -\widetilde{U}_i\widetilde{U}_j& \quad\text{}\\ = \widetilde{\widetilde{U}_i\widetilde{U}_j} - \widetilde{U}_i\widetilde{U}_j & \quad \text{: interaction between large structures} \\+ \widetilde{u_i^*\widetilde{U}_j}+ \widetilde{u_j^*\widetilde{U}_i}& \quad \text{: interaction between large structures and small eddies}\\+ \widetilde{u_i^*u_j^*} & \quad \text{: interaction between small eddies} \end{align} $$ - __Same structure to RANS representation__, but obtaining different physical meanings ---- #### --- Turbulence modeling in CFD --- ### Working with Eddy resolving simulations <img src="https://i.ytimg.com/vi/qEtcCjln-0Q/maxresdefault.jpg" width= "40%"><img src="https://www.gauss-centre.eu/fileadmin/research_projects/2016/Env_Energy/weihing_lutz_WEAloads_Fig01.jpg" width="41%"> - Simulations have to be carried out in __transient state__ of __3-dimensional mesh__, - turbulent structures are intrinstically 3D and unsteady. - computational cost is much higher than steady simulations - Instandt simulation results should be __realistic like__, - __eddy structures__ vary in every moment. - If not, __creating generic eddies__ (inflow generator) and/or __refining the mesh__ to create instability - Final flow characteristics are the __time/phase averaging__ of the simulation results, denoting - certain amount (>10) of physically run-though-time necessary - instant results can only serve as reference, no general physical contribution. ---- #### --- Turbulence modeling in CFD --- ### The Scale resolving turbulence models Common strategies are: - __Smagorinsky model (LES)__: filtering via cell size ($\Delta$) - flaw in near wall treatment, no modeling specified. Usually coupled with van Drist wall function - flaw in large grid size. Mesh study necessary - __Dynamic SGS eddy viscosity model (LES)__: slightly adjusting the filtering grid size. - still not capable for overrated cell size. - __WALE (Wall Adapted Local Eddy viscosity) (LES)__: imporved Smagorinsky with wall approach - __OneEquation LES (Deardorf)__: compensate the disadvantages in grid size restriction - Solving $k$ transport equation in adapting large cell size disadvantage. Mesh study not-necessary - __DES (Detached Eddy Simulation)__/__PANS (Partially Averaged Navier Stokes)__/__VLES (Very Large Eddy Simulation)__: involving RANS models structure to estimate the sub-scale properties. Switching deterministic length scale in $\Delta$ and turbulent length scale $L$. - __SAS (Scale Adpating Simulation__): URANS like simulation strategy, adapting filtered length scale from turbulent length scale $L$ to resolved flow length scale. --- #### --- Turbulence modeling in CFD --- ### URANS and Scale resolving simulations URANS: unsteady RANS, using RANS model for unsteady flows | | Unsteady RANS Simulation | Scale resolving Simulation | Laminar approach | | :-------- | :--------: | :--------: | :-------: | | Continuity Eqn. | $\large\frac{\partial \overline{U_i}}{\partial x_i} = 0$ | $\frac{\partial \widetilde{U_i}}{\partial x_i} = 0$ | $\large\frac{\partial {U_i}}{\partial x_i} = 0$| |Momentum Eqn. |$\frac{\partial \overline{U}_i}{\partial t} + \frac{\partial}{\partial x_j} ( \overline{U}_i \overline{U}_j + \overline{u_i^\prime u_j^\prime} ) \\= -\frac{1}{\rho}\frac{\partial \overline{P}}{\partial x_i} + \frac{\partial}{\partial x_j}\left( \nu \frac{\partial \overline{U}_i}{\partial x_j}\right)$ |$\frac{\partial \widetilde{U_i}}{\partial t} + \frac{\partial}{\partial x_j} \left( \widetilde{U_i} \widetilde{U_j}+ \widetilde{\tau}_{ij}\right) \\= -\frac{1}{\rho}\frac{\partial \widetilde{P}}{\partial x_i}+\frac{\partial}{\partial x_j}\left(\nu\frac{\partial \widetilde{U_i}}{\partial x_j}\right)$| $\frac{\partial {U}_i}{\partial t} + \frac{\partial}{\partial x_j} ( {U}_i {U}_j ) \\= -\frac{1}{\rho}\frac{\partial {P}}{\partial x_i} + \frac{\partial}{\partial x_j}\left( \nu \frac{\partial {U}_i}{\partial x_j}\right)$| <br/> Thinking on: In the aspect of CFD operation and matrix solving, <br/> what is the difference between $\widetilde{U}$, $\overline{U}$ and $U$ or $\widetilde{P}$, $\overline{P}$ and $P$? ---- #### --- Turbulence modeling in CFD --- ### URANS and Scale resolving simulations __No difference__ in Continuity and Momentum equations <br/> regarding URANS and Scale resolving methods. i.e., <img src="https://www.computationalfluiddynamics.com.au/wp-content/uploads/2013/09/Q-criterion_RANS-SRS.jpg" width="65%"> - All the transient simulation are under filtering concept. __Averaging does not really exist.__ - The __modeling length scale__ ($L (= k^{3/2}/\epsilon)$, $\Delta$, ) determines the filtering grid. - The smaller the length scale (filtering grid), the more structure will be resolved. - Large modeling length scales can suppress the coherent flow movement. - leads to __partially steady state__ results in certain region, result agreement to reality can be worse than steady state RANS. - URANS generally not recommanded for statistically steady state flows. i.e. no moving domain, no transient boundary condition. ---- #### --- Turbulence modeling in CFD --- ### Scale resolving simulations: Summary - By suppressing the modeling term ($\tau_{SGS}$), momentum equation __resolves more flow patterns__. __Eddy structure__ as feature of instant simulation results. - __3D domain, trasient setup__ necessary for simulation - Averaging of simulation results necessary to persue statistic flow patterns - Calculation time has to be considered to acquire sufficient data samples (> 10 x flow-through-time) - high computational costs - Instant result shall be different size of eddies resolved - Mostly applied in Aeroacoustic and combustion application. - Small scale structure as main focus. (noise, flame generation) - Unsteady RANS generally not recommanded for statistically steady state flow types. - no limitation on modeling length leads to oversupressing of flow patterns - LES relevant models are more appropriate capturing correct manner. --- #### --- Turbulence modeling in CFD --- ### Scenario: dynamic stall prediction <img src="https://upload.wikimedia.org/wikipedia/commons/thumb/f/f2/1915ca_abger_fluegel_%28cropped_and_mirrored%29.jpg/1024px-1915ca_abger_fluegel_%28cropped_and_mirrored%29.jpg" width="20%"> <img src="https://www.flexcompute.com/uploads/20220402-fig1-cl-vs-aoa.png" width="18%"> __[Stall](https://www.youtube.com/watch?v=SiOiVHUEYao)__ - Flow of the airfoil suction side separates. - The flow separation creates recirculation bubbles which evolute in vortex shapes. - The Airfoil lost its lift due to of non-smooth streamline profile. __Challenges__ - No defined separation point due to curved surface - separation point fluttering, leads to von Karman vortex street __Thinking on the flow structure resolution of simulation results:__ - What kind of Lift-AoA curve would I expect when EVM models is applied for airfoil? - For scale resolving strategy, what shall I take care in creating simulation case? ---- #### --- Turbulence modeling in CFD --- ### Scenario: dynamic stall prediction Validation of turbulence model for profile simulation: - pressure coefficient distribution along 2D airfoile <img src="https://www.mdpi.com/ijtpp/ijtpp-07-00026/article_deploy/html/images/ijtpp-07-00026-g005.png" width="50%"> - lift and drag coefficient in relation to angle of attack. <img src="https://media.springernature.com/lw685/springer-static/image/chp%3A10.1007%2F978-3-030-31307-4_12/MediaObjects/460925_1_En_12_Fig6_HTML.png" width="50%">