Nucleation Site Number Density

The nucleation site number density determines the number of locations on the heated surface where bubbles form, per unit area.

This is the leading factor determining the evaporation rate in a mechanistic model of subcooled boiling (see Eqn. (2110)). The diagram below illustrates this parameter:



A nucleation site is a gas-filled cavity in the wall. In practice, there is a range of cavity sizes. A larger wall superheat allows bubble growth in more cavities of smaller size.

Three standard submodels have been implemented for the nucleation site number density:

  • Lemmert Chawla – A correlation against wall superheat, and is often used together with the Tolubinsky Kostanchuk model for bubble departure diameter.
  • Hibiki Ishii – A mechanistic model for active nucleation site number density as a function of critical cavity size and contact angle. It is designed for use with the Kocamustafaogullari model for bubble departure diameter, and is applicable up to 198 bar.
  • Li – A mechanistic model based on a parametric analysis of existing experimental data. It is a function of three variables: wall superheat, pressure, and contact angle. The model is applicable for pressures ranging from 0.101 MPa–19.8 MPa.

The Lemmert Chawla model is the default option, as its fixed power for superheat dependence means that it is easier to obtain initial solutions. The more recent Hibiki Ishii model has a wider range of applicability, but it can be difficult to reach a converged solution in high flux cases. The Li model has a wider range applicability than the Hibiki Ishii model, particularly in accurately predicting results for small and moderate wall superheat values.

Under the 45 bar, 0.57 MW/m2 conditions of the first test in Bartolomei and Chanturiya [432], these two methods (each used with its recommended departure diameter model) can produce comparable results for axial void and wall superheat. However, the departure diameter that is associated with the Hibiki Ishii and Kocamustafaogullari solution is much smaller.

At higher pressures, and higher heat fluxes, the Hibiki Ishii model has an increasingly stronger dependence on wall superheat. This model often predicts a much smaller wall superheat than the Lemmert Chawla model for the same heat flux. At higher levels of wall superheat and when no limiter is applied to the nucleation site density, the Li model predicts smaller and more realistic values than the Hibiki-Ishii model.

Lemmert Chawla

The original Lemmert Chawla model [528][499] for nucleation site number density is:

1. EQUATION_DISPLAY
n [ m - 2 ] = ( m Δ T s u p ) p
(2126)

where:

  • n is the nucleation site density
  • m is a calibration constant, with a default value of 185.0/K
  • p is a superheat exponent, with a default value of 1.805
  • Δ T s u p is the wall superheat, which is calculated as:
2. EQUATION_DISPLAY
Δ T s u p = min  ( max  ( T w a l l - T s a t , 0 ) , Δ T m a x )
(2127)

where Δ T m a x is the maximum superheat applied to the Lemmert Chawla model.

Simcenter STAR-CCM+ uses a modified form of the Lemmert Chawla equations. This modification allows that the superheat exponent does not have to be constant, but can itself increase significantly with wall superheat.

3. EQUATION_DISPLAY
n n 0 " = ( Δ T s u p Δ T 0 ) A + B ( Δ T / Δ T 0 )
(2128)

where:

  • A and B are calibration constants
  • n 0 " = m p is the reference nucleation site number density. The default value is 12366.44 /m2 (which is 185.01.805).
  • Δ T 0 is the reference wall superheat

You can adjust the Lemmert Chawla Method so that it is suitable for use with the Kocamustafaogullari bubble departure diameters by comparing it with the Hibiki-Ishii Method, at a given operating pressure, and for a useful range of wall superheats. For example, the table below shows a fit of the Lemmert Chawla method to the properties of water assuming a contact angle of 0.722 radians:

P

(bar)

Superheat range (K) n 0

(m-2)

Δ T 0

(K)

A B Fit

Accuracy

45 1 - 10 138,100 1 1.348 0.35 7%
70 1 - 7 302,300 1 1.566 0.58 10%
150 1 - 3 5,665,000 1 3.095 1.7 0.1%

Using a fitted Lemmert Chawla method instead of using the Hibiki Ishii method gives similar results. The main advantage of such a fit is that it indicates the sensitivity of either method to wall superheat as an explicit superheat exponent. For example, using the calibrations above, the table below shows that, even at 45 bar, the nucleation site number density can vary as strongly as n = f ( Δ T 5 ) :

P (bar) Superheat range (K) Exponent range
45 1 - 10 1.7 - 4.8
70 1 - 7 2.2 - 5.6
150 1 - 3 4.8 - 8.2
Hibiki Ishii

For high wall superheats, the Hibiki Ishii method can produce extreme values of wall heat fluxes. To improve robustness, a limit is applied to the amount of superheat that can be applied:

4. EQUATION_DISPLAY
Δ T s u p = min  ( max  ( T w a l l - T s a t , 0 ) , Δ T m a x )
(2129)

The default value for the maximum superheat Δ T m a x is 25K. Do not allow this limit to affect the converged solution. If necessary, adjust the limit during the calculation.

Also, use system reference pressure in place of local liquid pressure to improve robustness.

5. EQUATION_DISPLAY
P l P R E F
(2130)

Following Hibiki and Ishii [473], the temperature of the vapor next to the wall is identified with the wall temperature.

6. EQUATION_DISPLAY
T g T s a t + Δ T s u p
(2131)

For the contact angle between the working fluid and the surface during boiling, θ , Hibiki and Ishii used a value at room temperature. They suggest that a temperature dependence of contact angle is implicitly included in the calibration of the pressure effects function f ( ρ + ) .

The nucleation site number density is then modeled as:

7. EQUATION_DISPLAY
n " = n " ¯ ( 1 exp [ θ 2 8 μ 2 ] ) ( exp [ f ( ρ + ) λ R c ] 1 )
(2132)

where:

  • n ¯ is the average cavity density, with a default value of 4.74E5 sites/m-2
  • θ is the wall contact angle, with a default value of 0.722 radians (41.37o), which is suitable for water systems
  • μ is the wall contact angle scale, with a standard value of 0.722 radians (41.37o)
  • λ is the cavity length scale, with default value 2.50E-6 m
  • f ( ρ + ) is the logarithmic, non-dimensional density function representing the effect of pressure
  • R c is the critical cavity radius

The critical cavity radius is:

8. EQUATION_DISPLAY
R c = 2 σ [ 1 + ( ρ g / ρ l ) ] / P l exp [ h lg ( T g T s a t ) / ( R T g T s a t ) ] 1
(2133)

where:

  • σ is the surface tension
  • ρ g is the vapor phase density
  • ρ l is the liquid phase density
  • P l is the liquid pressure
  • h l g is the latent heat
  • R is the specific gas constant for a vapor of given molecular weight.

The logarithmic, non-dimensional density function is:

9. EQUATION_DISPLAY
f ( ρ + ) = C 0 + C 1 ρ + + C 2 ρ +2 + C 3 ρ +3
(2134)

where:

10. EQUATION_DISPLAY
ρ + = log 10 ( ρ * )
(2135)

and:

11. EQUATION_DISPLAY
ρ * = ρ l - ρ g ρ g
(2136)

The standard values for the density function constants are:

12. EQUATION_DISPLAY
C 0 = - 0.01064 C 1 = 0.48246 C 2 = - 0.22712 C 3 = 0.05468
(2137)
Li

The Li model ([502] for nucleation site number density is based on a parametric analysis of existing experimental data. It is a function of three variables: wall superheat, pressure, and contact angle. The model is applicable for pressures ranging from 0.101 MPa–19.8 MPa. The superheat exponent of the Li model is not constant, but can increase significantly with wall superheat.

The Li model for nucleation site number density is:

13. EQUATION_DISPLAY
n = n 0 ( 1 cos θ ) exp [ f ( P ) ] Δ T s u p A Δ T s u p + B
(2138)

where:

  • n is the nucleation site density.
  • n 0 is the average cavity density, with a default value of 1000 sites/m2.
  • Δ T s u p is the wall superheat.
  • [ f ( P ) ] is a non-dimensional function representing the effect of pressure, given by:
    14. EQUATION_DISPLAY
    f ( P ) = 26.006 3.678 exp ( 2 P ) 21.907 exp ( P 24.065 )
    (2139)
  • P is the liquid pressure.
The term 1 cos θ accounts for the contact angle effects, defined as:
15. EQUATION_DISPLAY
( 1 cos θ ) = ( 1 cos θ 0 ) ( T c T s a t T c T 0 ) γ
(2140)
where:
  • θ 0 is the wall contact angle, with a default value of 0.722 radians (41.37o), which is suitable for water systems.
  • T c is the contact angle temperature with a default value of 647.15 K.
  • T 0 is the zero contact angle temperature with a default value of 287.15 K.

A and B are functions based on experimental data, which can be defined as a linear correlation on pressure:

16. EQUATION_DISPLAY
A = 0.0002 P 2 + 0.0108 P + 0.0119
(2141)
and:
17. EQUATION_DISPLAY
B = 0.122 P + 1.998
(2142)

Nucleation Site Number Density Limiting

Not all nucleation sites can be active at the same time, because some of the sites can be inside areas that are covered by existing bubbles that are still attached to the wall during the growth time of the bubbles. To account for this effect, Gilman [467] proposes a correction to the nucleation site number density calculation:

18. EQUATION_DISPLAY
n e f f = n exp [ f t g n π d w 2 4 ]
(2143)

where n e f f is the effective nucleation site number density after correction, n is the uncorrected nucleation site number density, f is the bubble departure frequency, t g is the bubble growth time, and d w is the bubble departure diameter.

However, in this model, n can be very large at high pressure such that the exponential term in Eqn. (2143) becomes zero, which is unphysical. To keep the calculated nucleation site density within physical limits, a further correction is required.

The area occupied by nucleating bubbles while they are growing and attached to the wall during the growth time is:

19. EQUATION_DISPLAY
A = t g π d w 2 4 n
(2144)

where n is the number of nucleating bubbles on the wall and t g is the bubble growth time.

The fraction of area and time that is occupied by nucleating bubbles in one bubble departure cycle is:

20. EQUATION_DISPLAY
f A T = t g ( t g + t w ) π d w 2 4 n A = f t g π d w 2 4 n
(2145)
where ( t g + t w ) is the time of one departure cycle and t w is the wait time between successive nucleating bubbles.

The bubble departure frequency is given by:

21. EQUATION_DISPLAY
f = 1 t g + t w
(2146)

The maximum n can be obtained from:

22. EQUATION_DISPLAY
f A T = 1 = f t g π d w 2 4 n max
(2147)
23. EQUATION_DISPLAY
n max = 1 f t g A p
(2148)

where A p = π d w 2 4 is the projected area of a bubble at departure. Note that f A T is the factor inside exponential function of Eqn. (2143). Instead of the exponential function, the following function is used to keep n within n max :

24. EQUATION_DISPLAY
n e f f = 1 1 n + 1 n max = n n max n + n max
(2149)

When n is much smaller than the maximum, Eqn. (2149) gives similar values to Gilman’s model. However, as n approaches the maximum, Eqn. (2149) returns the maximum value, hence keeping n within the physical limit.