Science and Technology of Nuclear Installations

Volume 2012 (2012), Article ID 939561, 10 pages

doi:10.1155/2012/939561

Research Article

CFD Analysis of a Void Distribution Benchmark of the NUPEC PSBT Tests: Model Calibration and Influence of Turbulence Modelling

E. Krepper and R. Rzehak

Institute of Fluid Dynamics, Helmholtz-Zentrum Dresden-Rossendorf, P.O. Box 510119, 01314 Dresden, Germany

Received 26 March 2012; Accepted 27 July 2012

Academic Editor: Diana Cuervo

Copyright © 2012 E. Krepper and R. Rzehak. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

The paper presents CFD calculations of the void distribution tests of the PSBT benchmark using ANSYS CFX-12.1. First, relevant aspects of the implemented wall boiling model are reviewed highlighting the uncertainties in several model parameters. It is then shown that the measured cross-sectionally averaged values can be reproduced well with a single set of calibrated model parameters for different test cases. For the reproduction of patterns of void distribution cross-sections, attention has to be focussed on the modelling of turbulence in the narrow channel. Only a turbulence model with the capability to resolve turbulent secondary flows is able to reproduce at least qualitatively the observed void distribution patterns.

1. Introduction

Based on NUPEC PWR Subchannel and Bundle Tests (PSBT), an international benchmark has been promoted by OECD and NRC and coordinated by Penn State University (PSU). The benchmark includes void distribution and departure from nucleate boiling exercises. In the first exercises the void fraction distribution was investigated in a steady state subchannel grade benchmark. In the present paper some of the single-channel steady state void fraction measurements are analysed to investigate the capabilities of present CFD modelling of wall boiling.

For engineering calculations, currently the most widely used CFD approach to model two-phase flows with significant volume fractions of both phases is the Eulerian two-fluid framework of interpenetrating continua (see, e.g., [1–3]). In this approach, balance equations for mass, momentum, and energy are written for each phase, that is, gas and liquid, separately and weighted by the so-called volume fraction which represents the ensemble averaged probability of occurrence for each phase at a certain point in time and space. Exchange terms between the phases appear as source/sink terms in the balance equations. These exchange terms consist of analytical or empirical correlations, expressing the interfacial forces, as well as the heat and mass fluxes, as functions of the average flow parameters. Since most of these correlations are highly problem-specific, their range of validity has to be carefully considered and the entire model has to be validated against experiments.

For the case of boiling flows, where heat is transferred into the fluid from a heated wall at such high rates that vapour is generated, additional source terms describing the physics of these processes at the heated wall have to be included. A CFD wall boiling model implemented in CFX following the lines of Kurul and Podowski [4, 5] was calibrated and validated by several authors, for example, Krepper et al. [6] against experimental results of Bartolomej and Chanturiya [7]. In these tests, subcooled flow boiling of water at high pressure flowing upwards in a vertical pipe heated from the outside was investigated and measurements of the axial development of void fraction, wall temperature and cross sectionally averaged liquid temperature were provided.

The aim of this work is to investigate the applicability of the CFX models to the PSBT tests. In Section 1 the PSBT tests for the void fraction distribution and the parameters of the selected tests are described briefly. In Section 2 the most important model details used for the simulations are described and critically reviewed. In Section 3 the sensitivity of the different parts of the models is investigated in view of the knowledge discussed in Section 2. The parameters are calibrated to only a single test. In the following Section 4 these calibrated models are applied to additional tests without any further change of model parameters. The cross sectionally averaged values of the selected tests are summarized for comparison to the calculations of other participants in the benchmark. In Section 5 the influence of the turbulence model on the void fraction distribution is investigated.

2. The OECD Benchmark Test

The PWR Subchannel and Bundle Tests (PSBT) were conducted by NUPEC (1987–1993) within an extensive experimental campaign aimed at verifying the reliability of fuel assemblies used for commercial nuclear power plants [8, 9]. Void fraction measurements and departure from nucleate boiling (DNB) tests were performed under PWR thermal-hydraulic conditions including steady states and transients such as power increase, flow reduction, depressurization, and temperature increase. The void fraction in each experiment was measured by gamma-ray transmission. These tests form the basis of an OECD benchmark for CFD and subchannel codes [10].

The subchannel test section, as shown in Figure 1, simulates a single subchannel of a PWR fuel assembly. The effective heated length is 1500 mm where the void measurement section is located near close to the top end at 1400 mm from the bottom of the heated section.

Figure 1: Test section for the subchannel void measurement [10].

For the analysis different tests in different pressure regions were selected (see Table 1).

Table 1: Selected tests.

Run
Pressure [MPa]
Mass flux [10⁶ kg m⁻² h⁻¹]
Inlet temperature [°C]
Power [kW]
Outlet quality

1.2211
15
11
295.4
90
−0.04

1.2223
15
11
319.6
70
0.04

1.2237
15
11
329.6
60
0.08

1.2422
15
5
284.1
60
0.02

1.2423
15
5
299.3
60
0.10

1.4324
10
5
238.9
60
0.00

1.4325
10
5
253.8
60
0.05

1.4326
10
5
268.8
60
0.11

1.4411
10
2
253.7
20
0.00

1.6222
5
5
204.2
50
0.02

3. The Models

3.1. Modelling of Boiling at a Heated Wall: The General Model Structure

In boiling, heat is transported from the hot wall to the fluid by several different mechanisms. On parts of the wall, where no bubbles reside, heat flows directly to the subcooled liquid in the same way as in single-phase flow. On parts of the wall where bubbles grow, heat is consumed by the generation of vapour which occurs at the so-called nucleation sites. Moreover, there is a liquid mixing mechanism due to the bubbles which leave the wall. As a consequence of the recirculation around the detaching bubbles, cold liquid from the bulk of the flow is brought into contact with the hot wall which leads to additional cooling. This mechanism is termed quenching. The total supplied heat flux is accordingly expressed as the sum of three contributions as

where , , and denote the heat flux components due to single-phase turbulent convection, quenching, and evaporation, respectively. The individual components in this heat flux partitioning are then modelled as functions of the wall temperature and other local flow parameters. Once this is accomplished, (1) can be solved iteratively for the local wall temperature , which satisfies the wall heat flux balance. Denoting the fraction of area influenced by the bubbles as , the heat flux components are expressed as discussed in the following.

The turbulent convection heat flux is calculated in the CFX model version in much the same way as for a pure liquid flow without boiling, but multiplied by the fraction of area unaffected by the bubbles, that is,

Here, is the heat transfer coefficient which is written using the temperature wall function known from Kader [11] as

where nondimensional variables (indicated by superscript “+”) and the friction velocity are defined as usual. Note that (2) and (3) may be evaluated at any location , provided it is used consistently whenever a variable depends on position. is represented in terms of the quenching heat transfer coefficient :

A grid-independent solution for is obtained by evaluating the non-dimensional temperature profile at a fixed value of . The evaporation heat flux is obtained via the evaporation mass flux at the wall:

where the generated vapour mass    is expressed in terms of the bubble diameter at detachment , bubble generation frequency , and nucleation site density    as

In terms of bubble detachment diameter and nucleation site density, the wall area fraction influenced by vapour bubbles is given by

Here, is the so-called bubble influence factor, for which a value of 2 is commonly used [4, 5]. Since corresponds to the case where the whole surface is under the influence of bubbles, as calculated by (7) has to be limited to values smaller than one. Moreover, it should be kept in mind that already as approaches 1 the assumptions of the model are not really satisfied anymore.

Correlations for the yet undetermined quantities used in the CFX wall boiling model are discussed in the following.

3.2. Bubble Size at Detachment

The bubble size at detachment depends on the liquid subcooling. Also the liquid properties which depend on the system pressure, the flow rate, and the heat flux have an influence.

An investigation of the bubble size at detachment was performed by Tolubinsky and Kostanchuk [12] for water at different pressures and subcoolings. The observed dependence on the liquid subcooling at atmospheric pressure can be fitted to a correlation:

A least squares fit gives parameter values are  m and  K under these conditions.

To match the tests of Bartolomej and Chanturiya [7] which were conducted at much higher pressures relevant for typical nuclear energy applications the values of and had to be adjusted to  mm and  K (e.g., [6]).

3.3. Nucleation Site Density

The situation concerning data on nucleation site density is much less clear. Most of the time, correlations are expressed in the form of power laws depending on the wall superheat as

However, a recent compilation [13] shows that vastly different parameter values are required to match different data sets. A likely reason for this fact is that nucleation site density is highly dependent on the microscale topography of the boiling surface, which in turn depends strongly on the processes that were used to finish the surface. These processes are very diverse and in most boiling experiments not specifically controlled or characterized.

For the tests of Bartolomej and Chanturiya [7] and  K were found to yield the best results in the model framework.

3.4. Further Model Aspects

The bubble detachment frequency    is given according to Cole [14] as a function of the detachment size . The quenching heat transfer coefficient is calculated as suggested by Mikic and Rohsenow [15].

In the bulk, vapour is assumed to be at saturation condition. Where the liquid is subcooled, that is, , vapour is condensing and calculated with the transfer rate:

Here, is the interfacial heat transfer coefficient, calculated according to Ranz and Marshall [16] and is the interfacial area. To close the phase transition model in the bulk bubbly flow with a mean bubble diameter , Kurul and Podowski [4] and also Anglart et al. [17] proposed to calculate the bubble diameter locally as a linear function of liquid subcooling :

Reference subcooling conditions for typical nuclear energy applications have been given as  mm at  K and  mm at   K in [17].

In a two-phase turbulent flow gaseous bubbles have an influence on the liquid turbulence. This effect is described by increasing the turbulent viscosity according to Sato et al. [18]. Conversely the liquid turbulence flow structures influence the gas distribution, which is described by a turbulent dispersion force according to Burns et al. [19].

For momentum exchange between the phases, finally, the Ishii and Zuber [20] drag law was used. Furthermore, a lift force according to Tomiyama et al. [21], a turbulent dispersion force according to Burns et al. [19], and a wall force according to Antal et al. [22] were included.

3.5. Evaluated Quantities

Quantities requested for the benchmark have been evaluated as follows. The mixed density was calculated according to

and the mixed enthalpy according to

with the temperature-dependent static enthalpies for liquid and gas . The quality was then calculated according to

with the enthalpies at saturation for liquid and gas .

4. Model Setup

4.1. General Setup

For most of the tests calculations were done for only 1/8 of the channel making use of the geometric symmetries. The mesh consisted of 720000 nodes; 1000 axial subdivisions were considered. To check the validity of this simplification for some tests furthermore 1/4 of the channel (850000 nodes) or the whole (1400000 nodes) were simulated with the same results. An inlet condition was set at the bottom. Upstream of the heated length an unheated part of 0.5 m length was simulated to model flow development. For the liquid and gas properties the water steam tables according to IAPWS IF97 were applied. In Section 3 the turbulence of the liquid phase was modelled by a Shear Stress Transport (SST) model [23]. In Section 4 two different Reynolds Stress Turbulence models were applied for comparison. The influence of the turbulence model is investigated in Section 5.

4.2. Calibration of Boiling Model Parameters

4.2.1. Bubble Size at Detachment

Earlier investigations [24] have shown that the bubble size at detachment has a sensitive influence on the calculated gas fraction since it determines the part of heat transferred by evaporation. To demonstrate the influence of this parameter different values (8) were investigated and the resulting cross sectionally averaged gas fraction for the test case 1.4324 is shown in Figure 2. Nevertheless, because of similar conditions as in the Bartolomej and Chanturiya [7] test cases,  mm was assumed for the further simulations.

Figure 2: Cross-sectionally averaged gas fraction for different values (8).

4.2.2. Nucleation Site Density

The nucleation site density has almost no influence on the liquid temperature, a small influence on the gas void fraction but a strong influence on the wall superheating . The common lack of information on the nucleation site density can be compensated by adjusting to match the measured wall temperature. Unfortunately, the wall temperatures are not available in the PSBT data. For the test case 1.4324 different values of the reference nucleation site density (9) were assumed to calculate the circumferentially averaged nucleation site density    (see Figure 3), bubble influence area (see Figure 4), wall temperature (see Figure 5), and the cross sectionally averaged gas volume fraction (see Figure 6). Figure 7 shows furthermore for   the distribution of the area fraction occupied by gas on the heated wall. For axial heights > 1 m (Figure 4) and for different areas (Figure 7) is found which does not seem reasonable. Therefore, was set to a smaller value of  m⁻² for all simulated cases. For all simulated cases compliance with the condition was checked.

Figure 3: Nucleation site density    (9).

Figure 4: Bubble influence area (7).

Figure 5: Wall temperature.

Figure 6: Cross sectionally averaged gas fraction.

Figure 7: Distribution of the area fraction occupied by gas (7) on the heated wall for (compare Figure 4). Only 1/8 of the channel is simulated considering the geometric symmetries.

5. Results

The calculations for the cross-sectional averages listed in Table 2 were performed applying a Shear Stress Transport (SST) model for the liquid phase turbulence. This model blends a model in the bulk and a model in the near wall region [23]. Furthermore, the corresponding Reynolds Stress model, the Baseline Reynolds Stress model (BSL-RSM) were investigated. This model has the same model structure as the SST model. In some cases the convergence of the BSL-RSM runs was questionable. The Reynolds Stress model (O-RSM) proved to be numerically more stable. The details of these approaches to turbulence modeling can be found in the ANSYS solver manual.

Table 2: Results.

Run
Pressure [MPa]
Turbulence model
Mixed density [kg m⁻³]
Void fraction [−]
Quality
[Pa m⁻¹]

1.2211
15
Exp
610
0.038
−0.040

BSL-RSM
588.6
0.090
−0.043
17481.9

SST
572.9
0.117
−0.037
17805.2

1.2223
15
Exp
456
0.311
0.040

O-RSM
471.2
0.266
0.073
13685.8

SST
494.5
0.212
0.059
15403.9

1.2237
15
Exp
390
0.440
0.080

O-RSM
431.9
0.327
0.110
12175.5

SST
447.2
0.291
0.101
15725.2

1.2422
15
Exp
522
0.182
0.020

BSL-RSM
520.8
0.187
0.015
8865.2

O-RSM
513.2
0.202
0.017
9363.2

SST
507.3
0.215
0.017
9468.4

1.2423
15
Exp
357
0.508
0.100

O-RSM
424.3
0.346
0.095
8927.7

SST
438
0.325
0.082
9401.9

1.4324
10
Exp
600
0.157
0.000

BSL-RSM
593.8
0.166
−0.002
5751.6

SST
578.9
0.190
0.000
6604.5

1.4325
10
Exp
478
0.335
0.050

BSL-RSM
473
0.336
0.048
9796.8

SST
470
0.343
0.046
10279.4

1.4326
10
Exp
353
0.531
0.110

BSL-RSM
363.3
0.503
0.099
9700

SST
361.9
0.506
0.098
10503.2

1.4411
10
Exp
599
0.152
0.000

BSL-RSM
614.6
0.124
0.000
7506.5

O-RSM
616.1
0.121
0.000
7541.9

SST
617.1
0.122
−0.004
7691

1.6222
5
Exp
549
0.306
0.020

O-RSM
578.9
0.263
0.016
7224.8

Exp: experiment.

BSL-RSM: Baseline Reynolds Stress.

O-RSM: Omega Reynolds Stress.

SST: Shear Stress Transport.

Calculating the cross sectionally averaged values, no significant differences between the different turbulence models could be discerned (examples for selected tests see Figure 8). In Section 5, however, a strong influence of turbulence modelling on the gas fraction distribution is shown.

Figure 8: Cross sectionally averaged values for runs 1.4324, 1.4325, and 1.4326; solid lines: Shear Stress Turbulence model (SST), dotted lines: Reynolds Stress Turbulence model (BSL).

The experiments do not provide values on bubble sizes or on temperatures of the heated wall. Exploiting this information gap for calibration of the correlations for the bubble size at detachment and the nucleation site density, at least cross-sectionally averaged values for the mixed density and gas fractions can be calculated with good agreement to experimental data.

6. Influence of the Turbulence Model on the Gas Fraction Distribution

Considering the cross-sectional distribution of the gas fraction, remarkable differences between the calculations using different turbulence models were found. Applying a turbulence model with a scalar turbulent viscosity, like the SST-model, the maximum of the gas was found near the heated wall in all cases. Figure 9(a) shows the calculated distribution for run 1.2422. In contrast, the measured gas distributions published in the benchmark specifications show significant amounts of gas also in the core of the channel. This is more similar to the results obtained by a Reynolds Stress model (see Figure 9(b)). The reason for these different patterns is probably that the SST model is not able to capture the secondary flows (see Figure 10).

(a)    Shear Stress Turbulence model (SST)

(b) Reynolds Stress Turbulence model (O-RST)

Figure 9: Cross sectionally averaged gas fraction at an axial distance of 1.4 m for run 1.2422.

(a)    Shear Stress Turbulence model (SST)

(b) Reynolds Stress Turbulence model (O-RST)

Figure 10: Velocity stream lines projected in a plane at an axial height of 1.4 m for run 1.2422.

Figure 11 shows the comparison of different runs with increasing inlet temperature (compare Table 1). With increasing inlet temperature the maximum gas fraction is shifted from the corner zones towards the core regions of the channel. The figure clearly shows the influence of the liquid flow field on the gas distribution similar to the measurements.

Figure 11: Tendencies of gas fraction distribution for the runs 1.4324, 1.4325, and 1.4326.

7. Conclusion

Current CFD implementations of a wall boiling model are able to calculate the gas fractions in good agreement with measurements. The model review in Section 3 has shown which information is necessary to perform meaningful calculations. With suitably calibrated correlations for the bubble size at detachment and the nucleation site density, at least cross sectionally averaged values for the mixed density and gas fractions can be calculated with good agreement to experimental data. In the corresponding experiments besides the measurement of gas fraction the determination of the temperature at the heated wall is urgent necessarily. Further research to reduce the number of model parameters that have to be determined by model calibration would be highly desirable. The distribution of gas within a cross section has been shown to depend strongly on the modelling of turbulence. Indication has been given that common two-equation models are too simplified for subchannel geometries while anisotropic Reynolds Stress models show a better potential. In order to make a stringent qualification of models, however, much more detailed CFD-grade data are urgently needed. While all of our calculations were performed for steady-state conditions, numerical convergence was a subject requiring special attention which could indicate that the real solution might contain transient elements.

Acknowledgment

The work was funded by the German Federal Ministry of Education and Research under the contract no. 02NUK010A.

References

D. A. Drew and S. L. Passman, Theory of Multicomponent Fluids, Springer, New York, NY, USA, 1998.

G. H. Yeoh and J. Y. Tu, Computational Techniques for Multiphase Flows—Basics and Applications, Butterworth-Heinemann, Elsevier Science and Technology, 2010.

M. Ishii and T. Hibiki, Thermo-Fluid Dynamics of Two-Phase Flow, Springer, NewYork, NY, USA, 2nd edition, 2011.

N. Kurul and M. Z. Podowski, “Multidimensional effects in forced convection subcooled boiling,” in Proceedings of the 9th International Heat Transfer Conference, Jerusalem, Israel, 1990.

N. Kurul and M. Podowski, “On the modeling of multidimensional effects in boiling channels ANS,” in Proceedings of 27th National Heat Transfer Conference, Minneapolis, Minn, USA, 1991.

E. Krepper, B. Končar, and Y. Egorov, “CFD modelling of subcooled boiling—concept, validation and application to fuel assembly design,” Nuclear Engineering and Design, vol. 237, no. 7, pp. 716–731, 2007.

G. G. Bartolomej and V. M. Chanturiya, “Experimental study of true void fraction when boiling subcooled water in vertical tubes,” Thermal Engineering, vol. 14, pp. 123–128, 1967, translated from Teploenergetika, vol. 14, no. 2, pp. 80–83, 1967

K. Hori, Y. Akiyama, K. Miyazaki, and etal, “In bundle void fraction measurement of PWR fuel assembly,” in Proceedings of ICONE 2, vol. 1, pp. 69–76, San Francisco, Calif, USA, March 1993.

K. Hori, Y. Akiyama, K. Miyazaki, T. Kurosu, and S. Sugiyama, “Void fraction in a single channel simulating one subchannel of a PWR fuel assembly,” in Two-Phase Flow Modeling and Experimentation, pp. 1013–1027, 1995.

JNES, “OECD/NEA Benchmark based on NUPEC PWR sub-channel and bundle tests (PSBT),” Assembly Specification and Benchmark Database (Volume 1), JNES/SAE-TH08-0019, April 2009

B. A. Kader, “Temperature and concentration profiles in fully turbulent boundary layers,” International Journal of Heat and Mass Transfer, vol. 24, no. 9, pp. 1541–1544, 1981.

V. I. Tolubinsky and D. M. Kostanchuk, “Vapour bubbles growth rate and heat transfer intensity at subcooled water boiling heat transfer,” in Proceedings of the 4th International Heat Transfer Conference, vol. 5, Paris, France, 1970, paper no. B-2.8.

N. I. Kolev, “Uniqueness of the elementary physics driving heterogeneous nucleate boiling and flashing,” Nuclear Engineering and Technology, vol. 38, no. 1, pp. 33–42, 2006.

R. Cole, “A photographic study of pool boiling in the region of the critical heat flux,” AIChE Journal, vol. 6, no. 4, pp. 533–542, 1960.

B. B. Mikic and W. M. Rohsenow, “A new correlation of pool-boiling data including the fact of heating surface characteristics,” ASME Journal of Heat Transfer, vol. 91, no. 2, pp. 245–250, 1969.

W. E. Ranz and W. R. Marshall, “Evaporation from drops,” Chemical Engineering Progress, vol. 48, no. 3, pp. 141–146, 1952.

H. Anglart, O. Nylund, N. Kurul, and M. Z. Podowski, “CFD prediction of flow and phase distribution in fuel assemblies with spacers,” Nuclear Engineering and Design, vol. 177, no. 1–3, pp. 215–228, 1997.

Y. Sato, M. Sadatomi, and K. Sekoguchi, “Momentum and heat transfer in two-phase bubble flow-I. Theory,” International Journal of Multiphase Flow, vol. 7, no. 2, pp. 167–177, 1981.

A. D. Burns, T. Frank, I. Hamill, and J. M. Shi, “The favre averaged drag model for turbulence dispersion in Eulerian multi-phase flows,” in Proceedings of the 5th International Conference on Multiphase Flow (ICMF'04), Yokohama, Japan, 2004.

M. Ishii and N. Zuber, “Drag coefficient and relative velocity in bubbly, droplet or particulate flows,” AIChE Journal, vol. 25, no. 5, pp. 843–855, 1979.

A. Tomiyama, H. Tamai, I. Zun, and S. Hosokawa, “Transverse migration of single bubbles in simple shear flows,” Chemical Engineering Science, vol. 57, no. 11, pp. 1849–1858, 2002.

S. P. Antal, R. T. Lahey, and J. E. Flaherty, “Analysis of phase distribution in fully developed laminar bubbly two-phase flow,” International Journal of Multiphase Flow, vol. 17, no. 5, pp. 635–652, 1991.

F. R. Menter, “Two-equation eddy-viscosity turbulence models for engineering applications,” AIAA Journal, vol. 32, no. 8, pp. 1598–1605, 1994.

E. Krepper and R. Rzehak, “CFD for subcooled flow boiling: simulation of DEBORA experiments,” Nuclear Engineering and Design, vol. 241, no. 9, pp. 3851–3866, 2011.