Boiling Model
Boiling: Theoretical Foundations
Overview
Professor, how does CFD for boiling work? Can you simulate water boiling in a kettle?
Of course you can. However, boiling phenomena are extremely complex, involving multiphase flow problems where bubble generation/detachment/condensation at the wall and liquid film behavior are intertwined. In CFD, we use wall heat flux partitioning models to describe the heat transfer associated with boiling.
What is a wall heat flux partitioning model?
The most representative one is the RPI model (Rensselaer Polytechnic Institute model), proposed by Kurul & Podowski (1990). It's a concept that decomposes the total heat flux from the wall into three components.
Governing Equations
Please show me the equation for the RPI model.
The total wall heat flux $q_w$ is partitioned into the following three components.
The meaning of each component is as follows.
- $q_{fc}$: Single-phase convective heat transfer (portion where liquid covers the wall)
- $q_{quench}$: Quenching heat flux (transient heat transfer when cold liquid contacts the wall after bubble detachment)
- $q_{evap}$: Evaporative heat flux (heat directly consumed for bubble generation)
What are the specific equations for each?
The single-phase convection component is based on the area fraction of the wall in contact with the liquid $(1 - A_b)$.
The quenching component is related to the bubble detachment frequency $f$ and waiting time.
The evaporation component is determined from the active nucleation site density $N_a$ on the wall, bubble departure diameter $d_w$, and detachment frequency $f$.
$A_b$ is the fraction of the wall area covered by bubbles, right?
Yes. $A_b = \min\left(1,\; K \frac{\pi d_w^2}{4} N_a\right)$, where $K$ is an empirical constant. The bubble departure diameter $d_w$ often uses the Tolubinsky & Kostanchuk (1970) model.
How do you determine the active nucleation site density $N_a$?
The correlation by Lemmert & Chawla (1977) is commonly used.
Here, $\Delta T_{sup} = T_w - T_{sat}$ is the wall superheat. $n$ is typically 1.805, and $C$ is a parameter determined from experiments.
Boiling Regimes
Are there different types of boiling?
Boiling regimes transition according to wall superheat. Taking pool boiling as an example, it is organized by the Nukiyama curve (boiling curve).
| Regime | Wall Superheat | Characteristics |
|---|---|---|
| Natural Convection | $\Delta T_{sup} < 5$ K | No bubbles, single-phase convection |
| Nucleate Boiling | 5โ30 K | Bubbles detach from wall, high heat transfer coefficient |
| Transition Boiling | 30โ100 K | Unstable, alternating liquid and vapor film formation |
| Film Boiling | $> 100$ K | Vapor film covers wall, low heat transfer coefficient |
The transition point from nucleate to film boiling is CHF (Critical Heat Flux: Critical Heat Flux), which is the most important parameter in nuclear reactor safety assessment.
So it's dangerous to exceed CHF.
Exceeding CHF causes a rapid rise in wall temperature leading to burnout. In nuclear engineering, DNBR (Departure from Nucleate Boiling Ratio) is managed as a safety margin during design.
Nukiyama CurveโThe Japanese Who Discovered Boiling's "Cliff"
In 1934, Shizuo Nukiyama of Tohoku University precisely measured the relationship between heat flux and wall temperature by varying the electrical power to a wire immersed in water. The result was an "S-shaped" curve (Nukiyama curve), where heat flux first increases to a maximum value (Critical Heat Flux, CHF), then wall temperature rises sharply before stabilizing again. This is because the moment CHF is exceeded, the wall becomes covered by a vapor film, drastically reducing the heat transfer coefficient. It represents a design limit that must never be exceeded in reactors or evaporators. This discovery became the starting point for boiling engineering and is still used worldwide 90 years later as a validation benchmark for CFD boiling models.
Computational Methods for Boiling
Details of Numerical Methods
When implementing a boiling model in CFD, what framework is used to solve it?
Boiling analysis is mainly based on the Euler-Euler two-fluid model, with the RPI model incorporated as a wall boundary condition. The continuity and momentum equations are solved for the liquid and vapor phases separately.
The transport equation for the vapor volume fraction $\alpha_v$ includes source terms $\dot{m}$ for evaporation and condensation.
The evaporation rate is determined from the RPI model's $q_{evap}$, right?
Exactly. The evaporative mass flux at the wall is as follows.
Condensation in the bulk is calculated from the interfacial heat transfer coefficient obtained via the Ranz-Marshall correlation, based on heat exchange with the subcooled liquid around the bubbles.
Bubble Force Models
How do bubbles generated by boiling move?
Modeling the forces acting on bubbles is important. The following interphase forces are considered.
| Force | Model | Role |
|---|---|---|
| Drag Force | Schiller-Naumann, Ishii-Zuber | Governs bubble velocity difference |
| Lift Force | Tomiyama | Lateral force due to velocity gradient |
| Wall Lubrication Force | Antal | Pulls bubbles away from the wall |
| Turbulent Dispersion Force | Lopez de Bertodano | Turbulent diffusion of bubbles |
| Virtual Mass Force | Auton | Acceleration effect |
Can Tomiyama's lift force change sign?
Good question. When the bubble diameter exceeds a critical value based on the Eรถtvรถs number $Eo$, the direction of the lift force reverses. Small bubbles move towards the wall, large bubbles move towards the pipe center. This is important physics determining wall-peaking and core-peaking void distributions.
Treatment of Wall Functions
Can normal wall functions be used on boiling surfaces?
They cannot. This is because bubble agitation significantly alters the near-wall flow structure compared to single-phase flow. Software like Fluent implements boiling-specific wall functions that perform wall temperature calculations consistent with the RPI model.
For wall meshing, it's more important that the first cell height is appropriate relative to the bubble departure diameter $d_w$ than constraints on $y^+$. Generally, it's recommended that the first cell height be at least $d_w$.
Time Step and Stability
Why does boiling analysis tend to diverge easily?
Because the volume expansion accompanying evaporation is rapid, generating large local volume source terms. The following countermeasures are effective.
- Gradually increase wall superheat (ramp-up)
- First obtain a single-phase steady solution, then activate boiling from there
- Set time step to $10^{-4}$ s or less
- Apply under-relaxation (0.3โ0.5) to the volume fraction equation
Dissecting the RPI Wall Boiling ModelโThe De Facto Standard for Nucleate Boiling CFD
The de facto standard for analyzing nucleate boiling in CFD is the RPI model. It decomposes the heat flux from the wall to the liquid into q_evap (evaporation) + q_quench (quenching) + q_conv (single-phase convection), calculating each from nucleation site density, bubble departure diameter, and detachment frequency. This model is implemented in ANSYS Fluent/CFX, but since the coefficient in the nucleation site density equation strongly depends on experimental conditions, caution is needed when extrapolating to new conditions. One researcher warns that "the RPI model is a stack of three correlations, so the uncertainty of each correlation is amplified multiplicatively," and errors in CHF prediction exceeding 30% are not uncommon.
Experience the theory firsthand with the interactive simulator for this field
All Simulators