Analysis and reformulation of the $k$--$ω$ turbulence model for buoyancy-driven thermal convection

This study derives an analytical solution for the standard $k$--$\omega$ model in Rayleigh--Bénard convection to identify discrepancies in buoyancy treatment, leading to a reformulated model with two new algebraic functions that significantly improve predictions of mean temperature and turbulent heat flux across various buoyancy-driven flows.

The Gist

Imagine you are trying to predict how heat moves through a pot of water sitting on a stove. In the world of physics, this is called buoyancy-driven convection: hot fluid rises, cold fluid sinks, and they mix in a chaotic dance called turbulence.

For engineers designing things like nuclear reactors or building ventilation systems, they need a way to predict this heat movement without simulating every single swirling drop of water (which would take supercomputers years to calculate). Instead, they use a "shortcut" method called RANS (Reynolds-Averaged Navier-Stokes). Think of RANS as a weather forecast: it doesn't track every single raindrop, but it predicts the general pattern of the storm.

The most popular "forecasting tool" for this is a model called the k–ω model. However, for decades, this tool has had a blind spot. It works great for wind blowing over a wing (shear flow), but when it comes to heat rising from a hot floor (buoyancy), it often gets the numbers wrong. It's like a GPS that knows how to drive on a highway but gets completely lost in a city grid.

The Problem: The "Blind" GPS

The paper explains that the standard k–ω model doesn't know how to handle the "push" that heat gives to the fluid.

• The Old Way: Engineers tried to fix this by guessing. They added a "knob" (a mathematical constant) to the model, turning it up or down based on whether the air was stable or unstable. But there was no rulebook. One software turned the knob to 1, another to 0, and another to -2. It was a mess of guesswork, and the results were often inaccurate, especially for fluids that are very thick (high Prandtl number) or very thin (low Prandtl number).

The Solution: A New Map

The author, Da-Sol Joo, decided to stop guessing and start deriving.

1. The Laboratory: Instead of looking at a messy, real-world room, the author created a perfect, simplified "laboratory" in math: a flat, infinite layer of fluid heated from the bottom (Rayleigh-Bénard convection). In this perfect world, the fluid doesn't move sideways; it only moves up and down. This allowed the author to solve the equations on paper to see exactly how the model should behave.
2. The Discovery: The math revealed that the standard model was predicting the wrong relationship between heat, fluid thickness, and temperature. It was like a scale that always weighed heavy objects as if they were light.
3. The Fix: The author didn't throw the whole model away. Instead, they added two tiny, smart adjustments (algebraic functions) to the model's "brain":
   • Adjustment 1 (For thin fluids): A tweak that changes how the model handles the "dissipation" (how quickly turbulence dies out) when the fluid is thin.
   • Adjustment 2 (For thick fluids): A tweak that changes how heat diffuses right next to the walls when the fluid is thick.

Crucially, these adjustments are smart. They only turn on when buoyancy (heat rising) is present. If there is no heat, the model reverts to its original, standard form. It's like adding a special lens to a camera that only activates when you are taking a photo of a sunset; for regular photos, the camera works exactly as it always did.

The Results: A Better Forecast

The author tested this new "corrected" model against a wide variety of scenarios, not just the simple lab setup:

• Heated rooms: Where heat comes from inside the room (like a nuclear reactor core).
• Mixed flows: Where wind is blowing and heat is rising at the same time.
• Different shapes: Tall, narrow rooms vs. wide, short rooms.

The Outcome:

• The old model often missed the mark by 50% or more in predicting how much heat was transferred.
• The new, corrected model hit the target with high accuracy across all these different situations.
• It successfully predicted how heat moves in fluids that are very thick (like oil) and very thin (like liquid metals), areas where the old model failed miserably.

The Big Picture

The paper argues that we don't need to build a completely new, overly complex machine to solve this problem. The existing "GPS" (the k–ω model) was just missing a few specific instructions for heat. By deriving the correct instructions from first principles and adding them as simple, smart tweaks, the author created a tool that is:

• Accurate: It predicts heat transfer correctly.
• Simple: It doesn't require massive new computing power.
• Robust: It doesn't crash or give weird answers when the conditions change.

In short, the paper takes a broken compass, figures out exactly why it was spinning in circles, and adds a tiny magnet to make it point North again, allowing engineers to navigate the complex world of heat-driven turbulence with confidence.

Technical Summary

1. Problem Statement

Buoyancy-driven turbulence (e.g., natural convection, Rayleigh–Bénard convection) is ubiquitous in engineering and geophysical flows. While Reynolds-averaged Navier–Stokes (RANS) models, particularly two-equation models like $k$–$\varepsilon$ and $k$–$\omega$, are widely used due to their cost-effectiveness, they lack a universally accepted standard formulation for incorporating buoyancy effects.

• The Core Issue: There is no analytical guidance for how to model the buoyancy production term in the scale-determining equation (the $\varepsilon$ or $\omega$ equation).
• Current Limitations: Existing approaches are largely empirical. Some omit buoyancy terms entirely, while others add them analogously to shear production. These methods often fail to reproduce correct scaling laws for the Nusselt number ($Nu$) as a function of the Rayleigh ($Ra$) and Prandtl ($Pr$) numbers, particularly under unstable stratification.
• Specific Deficit: Standard models often predict incorrect asymptotic behaviors, such as $Nu \sim Pr^{-0.415}$ for high $Pr$, whereas experimental and DNS data suggest $Nu$ should be nearly independent of $Pr$($Nu \sim Pr^0$) in this regime.

2. Methodology

The study employs a rigorous analytical approach combined with numerical validation to derive and correct the standard Wilcox (2006) $k$–$\omega$ model.

A. Analytical Framework

• Canonical Test Case: The author derives an exact analytical solution for the standard $k$–$\omega$ model applied to statistically one-dimensional Rayleigh–Bénard convection in an infinite layer.
• Isolation of Mechanisms: In this configuration, the mean velocity is zero, and turbulent kinetic energy ($k$) is generated solely by buoyancy. This isolates the role of the buoyancy term in the $\omega$ equation, which is difficult to achieve in shear-driven flows.
• Derivation: The solution explicitly relates the Nusselt number ($Nu$) to $Ra$ and $Pr$. It analyzes the balance between buoyant production, dissipation, and diffusion in both the near-wall and bulk regions to determine the scaling exponents.

B. Reformulation and Correction

Based on the analytical discrepancies identified, the author proposes a reformulation involving two dimensionless algebraic functions:

1. Correction to the Buoyancy Coefficient ($C_{\omega b}$):

   • The standard constant $C_{\omega b}$ is replaced by a function of $Pr$ to correct the low-$Pr$ scaling.
   • For unstable stratification ($P_b > 0$), the coefficient is modified to:
     $$C_{\omega b}^+ = -0.9752 - 0.2988 Pr^{-5/16}$$
   • For stable stratification ($P_b < 0$), an analytical derivation based on homogeneous stably stratified shear flow yields a new constant:
     $$C_{\omega b}^- = -0.5385$$
     (This replaces the empirical $C_{\omega b} \approx -2$ commonly used, ensuring the model stabilizes at the correct flux Richardson number threshold).

2. Correction to Turbulent Thermal Diffusivity ($a_T$):

   • To correct the high-$Pr$ behavior, a near-wall modification function $\psi$ is introduced to the turbulent thermal diffusivity ($a_T = \nu_T / Pr_T + \psi$).
   • This function vanishes in the absence of buoyancy and enhances near-wall heat flux to recover the $Nu \sim Pr^0$ scaling for $Pr \gg 1$.
   • The function depends on the turbulent Reynolds number ($Re_T$), $Pr$, and the ratio of buoyant production to dissipation ($P_b/\varepsilon$).

C. Validation

The corrected model is validated against:

• Analytical Solutions: Comparing the derived scaling laws with 1D simulations.
• DNS and Experimental Data: Validating against a wide range of datasets including:
  • Rayleigh–Bénard convection (1D and 2D square cavities).
  • Internally heated convection (top-cooling and top-and-bottom cooling).
  • Unstably stratified Couette flow (mixed convection).
  • Side-heated natural convection in rectangular cavities (varying aspect ratios).

3. Key Contributions

1. Analytical Solution for $k$–$\omega$: The paper provides the first explicit analytical solution for the standard $k$–$\omega$ model in a buoyancy-dominated, statistically 1D Rayleigh–Bénard flow. This reveals the mathematical origin of the model's scaling errors.
2. Derivation of Scaling Laws: The analysis explicitly derives the model's predicted scalings ($Nu \sim Ra^{1/3}Pr^{1/3}$ for low $Pr$ and $Nu \sim Ra^{1/3}Pr^{-0.415}$ for high $Pr$) and contrasts them with the physically observed trends ($Nu \sim Pr^{1/8}$ and $Nu \sim Pr^0$).
3. Systematic Reformulation: Instead of ad-hoc empirical tuning, the paper derives two algebraic correction functions that:
   • Recover the correct $Nu$–$Pr$ scaling across the entire range ($10^{-3} \le Pr \le 10^3$).
   • Ensure full compatibility with the standard model (the corrections vanish when buoyancy is absent).
   • Provide a theoretically consistent value for stable stratification ($C_{\omega b} = -0.5385$).
4. Robustness Across Regimes: The corrections are shown to be effective not just in pure Rayleigh–Bénard convection, but also in mixed convection (Couette flow) and internally heated flows.

4. Results

• Scaling Accuracy: The corrected model successfully reproduces the experimental trends:
  • For $Pr \ll 1$: $Nu \sim Pr^{1/8}$ (previously predicted as $Pr^{1/3}$).
  • For $Pr \gg 1$: $Nu \sim Pr^0$ (previously predicted as $Pr^{-0.415}$).
• Quantitative Improvement:
  • In Rayleigh–Bénard convection, the corrected model reduces the error in $Nu$ predictions by 50–70% compared to the standard model (which used $C_{\omega b}=1$).
  • In internally heated convection, the model accurately predicts the maximum temperature and heat flux partitioning, whereas the standard model significantly overpredicts temperatures.
  • In stable stratification, the analytically derived $C_{\omega b} = -0.5385$ provides better agreement with DNS data for the suppression of turbulence compared to the standard $C_{\omega b} = -2$.
• Mixed Convection: In unstably stratified Couette flow, the model accurately captures the transition between shear-dominated and buoyancy-dominated regimes, predicting $Nu$ and friction Reynolds numbers within 10% of DNS data.
• Side-Heated Cavities: In shear-dominated side-heated flows where buoyancy plays a minor role, the corrected model yields results virtually identical to the standard model, demonstrating that the corrections do not degrade performance in non-buoyancy-dominated flows.

5. Significance

• Theoretical Clarity: The study bridges the gap between phenomenological scaling theories (like Grossmann–Lohse) and RANS modeling by showing how standard closure structures can be analytically tuned to match these scalings.
• Engineering Utility: It offers a practical, low-cost solution for engineers. By introducing only two algebraic functions, the model achieves high accuracy without the computational cost or complexity of higher-order closures (like Reynolds Stress Models).
• Paradigm Shift: The paper argues against the notion that standard two-equation models have reached a "glass ceiling." It demonstrates that with a systematic understanding of the model's intrinsic behavior, significant predictive improvements can be achieved through minimal, physically motivated modifications.
• Standardization: The proposed formulation provides a consistent, reproducible framework for buoyancy modeling in $k$–$\omega$ codes (e.g., OpenFOAM, ANSYS Fluent), addressing the current fragmentation where different solvers use inconsistent empirical coefficients.

In conclusion, this work provides a systematic, analytically grounded reformulation of the $k$–$\omega$ model that resolves long-standing discrepancies in buoyancy-driven turbulence predictions, offering a robust tool for a wide range of thermal convection problems.