- Email: [email protected]
A subchannel based annular flow dryout model
Annals of Nuclear Energy 94 (2016) 313–324
Contents lists available at ScienceDirect
Annals of Nuclear Energy journal homepage: www.elsevier.com/locate/anucene
A subchannel based annular flow dryout model Najmeddine Hammouda a,⇑, Zhong Cheng b, Yanfei F. Rao b a b
Canadian Nuclear Safety Commission (CNSC), Ottawa, Ont. K1P 5S9, Canada Canadian Nuclear Laboratories (CNL), Chalk River, Ont. K0J 1J0, Canada
a r t i c l e
i n f o
Article history: Received 3 February 2016 Received in revised form 31 March 2016 Accepted 1 April 2016 Available online 8 April 2016 Keywords: Subchannel thermalhydraulics ASSERT-PV Safety Annular flow Dryout power CHF
a b s t r a c t This paper assesses a popular tube-based mechanistic critical heat flux model (Hewitt and Govan’s annular flow model (based on the model of Whalley et al.), and modifies and implements the model for bundle geometries. It describes the results of the ASSERT subchannel code predictions using the modified model, as applied to a single tube and the 28-element, 37-element and 43-element (CANFLEX) CANDU bundles. A quantitative comparison between the model predictions and experimental data indicates good agreement for a wide range of flow conditions. The comparison has resulted in an overall average error of 0.15% and an overall root-mean-square error of 5.46% with tube data representing annular film dryout type critical heat flux, and in an overall average error of 0.9% and an overall RMS error of 9.9% with Stern Laboratories’ CANDU-bundle data. Crown Copyright Ó 2016 Published by Elsevier Ltd. All rights reserved.
1. Introduction The objective of this study is to assess a popular tube-based annular flow (or liquid film) dryout model (Hewitt and Govan, 1990), and to modify and to implement the model for bundle geometries using the advanced solution of subchannel equation in reactor thermal-hydraulics (ASSERT) subchannel code (Carver et al., 1990; Hammouda et al., 2001; Rao et al., 2014a). ASSERT is a computer code that has been developed to model the subchannel flow and phase distribution in horizontal Canada uranium deuterium (CANDUÒ1) pressurised heavy water reactor fuel channels (Carver et al., 1990; Hammouda et al., 2001; Rao et al., 2014a). The code has been designed to be general enough to accommodate other geometries and orientations. These include single subchannels of different shapes, and multiple subchannels of CANDU PHWR, PWR and BWR designs, in both vertical and horizontal orientations. As well, the code can accommodate a range of fluids, including single- and two-phase heavy water, light water, various Freons, and two-phase air–water. ASSERT (Hammouda et al., 2001; Rao et al., 2014a; Rowe et al., 1993) is based on ASSERT-IV (Rowe et al., 1988), which, in turn originated from the COBRA-IV computer program (Wheeler and et al., 1976; Stewart and et al., 1977); ASSERT has been enhanced to meet the specific requirements for the
⇑ Corresponding author. E-mail addresses: [email protected] (N. Hammouda), zhong.cheng@cnl. ca (Z. Cheng), [email protected] (Y.F. Rao). 1 CANDUÒ is a registered trademark of Atomic Energy of Canada Limited (AECL). http://dx.doi.org/10.1016/j.anucene.2016.04.002 0306-4549/Crown Copyright Ó 2016 Published by Elsevier Ltd. All rights reserved.
thermalhydraulic analysis of two-phase flow in the horizontally oriented CANDU fuel. Currently, empirical equations are the most widely used critical heat flux (CHF) prediction methods, because of a lack of sufficient understanding of the phenomena governing CHF. A large number of CHF correlations are available in the literature; and some of the most popular for vertical flow in tubes are those of Bowring (1972), Griffith et al. (1977) and Biasi et al. (1967). Most CHF correlations were derived for high-pressure steam-water flow in tubes, and have relatively narrow ranges of validity. The number of physically based (or mechanistic) CHF models is very small. In general, mechanistic CHF models are developed for two main flow regimes: (i) the annular flow regime at highquality flows, where the mechanism for CHF is dryout of the liquid film; and (ii) the bubbly flow regime at subcooled or low-quality flows, where the most widely suggested mechanism for CHF is bubble crowding and coalescence near the wall to form a vapour film. Two popular models of CHF are the film dryout model of Whalley et al. (1974) for a rod-centred subchannel formulation, and the departure from nucleate boiling (DNB) model of Weisman and Pei (1983). ASSERT has a built-in capability to systematically determine CHF occurrence, based on its prediction of local subchannel flow and rod heat-flux conditions. It includes a number of options for computing CHF in the form of empirical correlations. Of these options, the CHF table look-up method of Groeneveld et al. (2005) is the most used, because it is applicable over a wide range of conditions. However, the CHF look-up table was derived on the
314
N. Hammouda et al. / Annals of Nuclear Energy 94 (2016) 313–324
Nomenclature A C d D E F Fw ~ g G h hf hg hfg hin, Dhin k P Ph Pw q00 q000 t Vg ~ V Vr X Xlo x
flow area (m2) mass concentration of droplets (kg m3) hydraulic diameter (m) deposition mass flux (kg m2 s1) entrainment mass flux (kg m2 s1) mass-flow rate [kg/s] frictional force in Table 1 gravitational acceleration vector in Table 1 mass flux (kg m2 s1) mixture specific enthalpy in Table 1 specific enthalpy of liquid phase in Table 1 specific enthalpy of vapour phase in Table 1 latent heat of vaporisation [J/kg] specific enthalpy at inlet, subcooling at inlet deposition mass-transfer coefficient (m s1) pressure (Pa) heated perimeter (m) wetted perimeter (m) heat flux (W m2) volumetric heat flux (W m3) time (s) vapour velocity perpendicular to liquid film interface (m s1) velocity vector in Table 1 lateral relative velocity Table 1 quality [–] fraction of liquid entrained at the start of annular flow [–] axial location [m]
basis of 8-mm tube CHF data, and some uncertainties remain in its application to subchannel geometries. Application of the table look-up method for bundle subchannel geometries required use of phenomenological correction models (Rao et al., 2014a). The CHF look-up table is a method of directly representing experimental data obtained for upward flow in 8-mm vertical tubes. Therefore, it may be viewed as a limit to what can be achieved using an empirical approach for CHF in vertical tubes. The parametric trends in CHF of the look-up table are those observed for steam-water flow inside vertical tubes and, for the most part, they have not been ascertained for flow in rod-bundle subchannels. For fluids other than light water, the CHF look-up table requires the use of fluid-to-fluid modelling techniques, which may add uncertainties to the table predictions. However, mechanistic CHF models can be very useful tools to assess parametric trends in subchannel geometries. Other advantages of physically based models are that they can be directly applied to other fluids, can readily be used to extrapolate outside the original data base, and can be used for various geometries and varying heat-flux shapes. Also, in regions where the CHF look-up table values are uncertain, because of scarcity of data, mechanistic CHF models can be used as possible alternatives. It should be noted that this is the first time the annular flow dryout model in ASSERT-PV is published together with an assessment of the model against tube and CANDU-bundle CHF experiments. Recent development and validation of ASSERT-PV CHF models have been focused on the CHF table look-up method plus the associated phenomenological correction models, as summarised in Rao et al. (2014a,b) for ASSERT-PV 3.2. The ASSERT annular flow dryout model was implemented into an earlier version 3.0 (Hammouda et al., 2001) and has stayed unchanged in version 3.2. The Canadian nuclear industry is planning to put more concentrated effort in developing and improving the annular flow dryout model for more accurate subchannel thermalhydraulic analyses.
Greek symbols: void fraction Cs deposition mass-flux suppression due to evaporation (kg m2 s1) Cvap vapour mass-flux generation source function (kg m2 s1) l viscosity (N s m2) q mass density (kg m3) r surface tension (N m1)
a
Subscripts: exp experimental f saturated liquid or liquid phase g saturated vapour or vapour phase l subcooled liquid LE liquid entrained LF liquid film i interface if interface to liquid ig interface to vapour pre predicted sat saturation v vapour w wall
2. Annular flow dryout model The majority of the available annular flow dryout models are based on that of Whalley et al. (1974). The difference from Whalley et al.’s model is in the constitutive relationships representing the mechanisms of entrainment and deposition. The basic Whalley et al. model is a three-field model representing two-phase interactions between vapour, liquid film and entrained droplets. The one-dimensional equations for conservation of mass in annular flow can be derived by considering mass balances written separately for each field,
dðql aLF Þ dGLF Pw Ph þ ¼ ðD EÞ ðCs þ Cv ap Þ dt dx A A
ð1Þ
dðql aLE Þ dGLE Pw Ph þ ¼ ðE DÞ þ Cs dt dx A A
ð2Þ
dðqv aÞ dGv Ph þ ¼ Cv ap dt dx A
ð3Þ
The droplet deposition rate is calculated from
D ¼ kC
ð4Þ
where C is the concentration of droplets (kg m3) in the vapour core, and k (m s1) is the deposition mass transfer coefficient. The concentration, C, is defined as
C¼F
F LE þ qF v
LE
ql
ð5Þ
v
In Eq. (5), FLE and Fv are the droplet and vapour flow rates, respectively.
N. Hammouda et al. / Annals of Nuclear Energy 94 (2016) 313–324
315
2.1. Deposition and entrainment models
2.3. Model of deposition suppression
The development of constitutive relations for entrainment and deposition is a critical step in calculating dryout in annular flow through a mechanistic approach. In fact, it is these relations that distinguish one mechanistic dryout model from another, since the basic conservation equations are the same for any dryout model in the annular flow regime. Hewitt and Govan (1990) derived improved models for deposition and entrainment in annular flow. The new entrainment and deposition correlations successfully predicted a wide range of equilibrium (where the rate of droplet entrainment is equal and opposite to the rate of droplet deposition) and non-equilibrium data. Specifically, the correlations were used to predict equilibrium entrainment and the variation of entrainment in developing flows of air–water and steam-water data (Govan et al., 1988). The distinguishing feature of Hewitt and Govan (1990) model is in its modelling of hydrodynamic non-equilibrium, which includes the fact that a limiting film flow rate (for very thin films) exists, below which no entrainment can take place. The correlation for droplet deposition includes the experimentally observed effect of a decreasing deposition mass-transfer coefficient with increasing droplet concentration. Subsequent studies by Vaitekunas et al. (1996) and Hoyer (1998) showed that Hewitt and Govan’s model of deposition and entrainment gives satisfactory results, when compared to tube data of various fluids. The correlation for the deposition rate coefficient is as follows:
Film evaporation may inhibit droplet deposition. The vapour flowing away from the evaporating film may prevent liquid droplets in the core from reaching the liquid film surface. Hoyer (1998) gives the following correlation to account for the suppressed deposition of droplets due to evaporation:
k¼
8 > <
0:18
qffiffiffiffiffiffi
r qv d
; qC 6 0:3
v
; qC > 0:3
qffiffiffiffiffiffi 0:65 > : 0:083 r C q d q v
v
ð6Þ
v
The entrainment correlation is
8 h i0:316 < E 5:75 105 ðGLF GLFC Þ2 rdqq2l ¼ v Gv : 0
; GLF > GLFC ; GLF 6 GLFC
ð7Þ
where GLFC, the critical film mass flux for the onset of entrainment, is given by
GLFC d
ll
rffiffiffiffiffiffi l ql ¼ exp 5:8504 þ 0:4249 v
ll
qv
ð8Þ
2.2. Model of additional entrainment The majority of entrainment and deposition correlations were derived from data for fully developed annular flow, which mainly describe entrainment and deposition in the absence of heat transfer. However, it has been experimentally shown by Milashenko et al. (1989) that the applied heat flux can have a significant effect on entrainment in annular flow through the boiling process. The mechanism of heat-flux-induced entrainment is the main cause of additional entrainment in diabatic annular flows. The additional entrainment is generated by the boiling process, and is due to the release of bubbles from the liquid film (Hoyer, 1998). The total entrainment (dynamic entrainment Ed, and bubble entrainment, Eq) is given by the correlation of Milashenko et al. (1989) as
Aw Cq ¼ 1:75F LF q00 106
qg qf
!1:3 ð9Þ
where Cq ¼ Ed þ Eq , and Aw is the heated area of the channel. Eq. (9) is applicable in the range of pressures from 5 to 10 MPa, mass fluxes from 1 to 3 Mg m2 s1 and heat fluxes from 0.15 to 4 MW m2. Hoyer (1998) recommends using the maximum value of entrainment obtained from either Eqs. (7) or (9).
Cs ¼ k q C
ð10Þ
where the suppression coefficient kq is given by
kq ¼
qv
0:065ql
Vg
ð11Þ
and where
V g ¼ Cv ap
A
qv Ai
ð12Þ
where Vg is the vapour velocity perpendicular to the liquid film area Ai (Hoyer, 1998). The numerical implementation of the annular film dryout model in ASSERT is described in the next section. 3. Modification and implementation in ASSERT-PV 3.1. Thermal-hydraulic model and solution scheme of ASSERT ASSERT uses the subchannel approach to model single- and two-phase flow in rod bundles. Subchannels are defined by the coolant flow areas between rods, which are bounded by the rods and imaginary planes linking adjacent rod centre lines. Subchannels are divided axially into control volumes, which communicate axially in the same subchannel, and laterally across fictitious boundaries with control volumes in neighbouring subchannels. The code is based on a five-equation, advanced drift-flux thermal-hydraulic model that accounts for the thermal and mechanical non-equilibrium and transverse flow found in CANDU-type fuel channels. The thermal-hydraulic conservation equations used in ASSERT are derived from the two-fluid equations (Ishii, 1975). Heat transfer, friction, and the exchange of mass, momentum and energy between subchannels affect the distribution of flow and enthalpy among individual subchannels. In addition to the conservation equations, ASSERT uses a number of closure relationships to model wall-to-fluid heat and momentum transfer and inter-subchannel mass, momentum and energy exchanges due to turbulent mixing, void drift and buoyancy drift (Carver et al., 1990; Hammouda et al., 2001). The governing equations are given in Table 1. ASSERT uses a staggered grid numerical solution scheme based on a pressure–velocity algorithm (Hammouda et al., 2001). 3.2. The annular flow dryout model: application to rod-bundle geometry For the complex geometry of rod bundles, subchannel analysis is required to calculate CHF, and to appropriately account for thermal non-equilibrium and intersubchannel mixing effects. Several attempts have been made to investigate the occurrence of CHF using subchannel analysis tools in the annular flow regime using the film dryout criterion (Whalley, 1977; Saito et al., 1983; Lim and Weisman, 1988). Whalley (1977) was among the first to apply subchannel analysis to predict film dryout in rod bundles. In his analysis, the bundle was divided into subchannels centred on the rods. On each rod, mass balances on the liquid film, liquid droplets and vapour were applied, accounting for inter subchannel mixing.
316
N. Hammouda et al. / Annals of Nuclear Energy 94 (2016) 313–324
Table 1 Governing equations in ASSERT.a Mixture mass: @q VÞ ¼ 0; þ r ðq~ @t
where q ¼ ðaqÞg þ ðaqÞf ¼ ag qþg þ af qþf VÞ ¼ ðaq~ VÞg þ ðaq~ VÞf ¼ q~ V ðq~ Mixture momentum h i ðaqÞg ðaqÞf ~ ~ V~ Vþ V r V r þ rP ¼ F w þ q~ þ r q~ g q
@ðq~ VÞ @t
Mixture energy: ~ q @h @t þ qV rh þ r
hðaqÞ ðaqÞ g f q
i 00 00 ðhg hf Þ~ V r ¼ q000 w r ½ðaq Þg þ ðaq Þf
Phasic energy (transportive form) Liquid: @h 00 000 VÞf rhf ¼ q000 ðaqÞf @tf þ ðaq~ wf r ðaq Þf þ qif Vapour: ðaqÞg a
@hg @t
00 000 VÞg rhg ¼ q000 þ ðaq~ wg r ðaq Þg þ qig
Here, + = variables to be defined by state relationships and * = variables to be defined by constitutive relationships.
The constitutive relations described in Sections 2.1–2.3 were derived from data obtained for tube geometry. In the analysis presented here, it is assumed that the mechanisms of hydrodynamic entrainment and deposition found in a tube are applicable also to rod-bundle subchannels. This assumption is based on the general notion that any one subchannel within a rod bundle is geometrically equivalent to a round tube, so that similar thermalhydraulic behaviour may be expected (Lim and Weisman, 1988). Since each subchannel of a rod bundle is bounded by a number of rods, it is expected that the liquid film flow is not necessarily the same for all rod sectors facing the same subchannel. The liquid film flow rate in each rod segment depends on the average conditions in the subchannel (void fraction, pressure, liquid droplet concentration, etc.), and on the heat flux of that rod segment. Hydrodynamic entrainment and deposition rates depend on the average subchannel flow conditions and the liquid film flow rate. On the other hand, heatflux-induced entrainment and deposition suppression explicitly depend on the rod heat flux. Therefore, within the same subchannel, mainly the rod having the highest heat flux or the highest radial heat-flux factor determines CHF. In this work, the conventional subchannel approach, in contrast to a rod-centred approach (Whalley, 1977), is used to solve Eqs. (1)–(3). The liquid film is calculated from Eq. (1) using rod heatflux-based correlations (Eq. (9) and Eqs. (10)–(12)) and subchannel-based correlations (Eqs. (5)–(7)). Eqs. (2) and (3) are not needed, since ASSERT calculates the average subchannel flows (mixture and phasic), void fraction and heat flux of each rod. ASSERT calculates mixture flow rates by solving mixture equations of momentum and mass in each subchannel, and provides phasic flow rates using the mixture flow and a drift flux model. ASSERT can provide only the total phasic flow of liquid, since a multi-field two-fluid model is necessary to compute droplet flow and liquid film flow separately. Given the total liquid flow, predicted by ASSERT, and the liquid film flow, predicted from Eq. (1), the entrained droplet flow rate can easily be determined. As a result, Eq. (2) is not needed. Similarly, since the phasic vapour flow is already computed by ASSERT model, Eq. (3) is not used.
The subscripts n, i and j denote rod n, subchannel i, and axial node j. For the purpose of CHF computation, only the perimeters of heated rod segments are considered in the calculation. Therefore, the rod segment wetted perimeter Pw can be replaced by the heated perimeter Ph in Eq. (13) to obtain
dF LF n ¼ ðD E Csn Cv apn Þi;j Phnj dx i;j
ð14Þ
The rod segment heated perimeter is
Phn ¼ p dn wi;n
ð15Þ
where dn is fuel rod n diameter, and wi,n is the fraction of rod n facing subchannel i. The vapour mass-flux generation source Cvap is expressed in terms of the rod segment heat flux as
ðCv apn Þi;j ¼
q00wn hfg
ð16Þ i;j
Inserting Eqs. (15) and (16) into Eq. (14), we obtain
" 00 # q dF LF n ¼ ðD Csn Þi;j Ei;j wn pdn wi;n dx i;j hfg i;j
ð17Þ
Eq. (17) is a first-order differential equation that can be integrated along each rod segment for all subchannels until the first occurrence of a zero liquid-flow rate, FLF, is obtained. For every rod segment facing a subchannel, Eq. (17) is solved as an initial value problem in the x-direction, using a fourth-order Runge–Kutta method. The initial values of all the variables in Eq. (17) are specified at axial location j at the start of the annular flow regime, as described below. Starting at j, the Runge–Kutta method advances the solution of Eq. (17) over an interval Dxj. Using the Runge–Kutta formula, Eq. (17) can be expressed as
dF LF n dx
i;j
¼
ðF LF F LF j
Þ
j1 n;i
Dxj
¼ 16 f ððF LF n Þi;j1 ; xj1 Þþ ð18Þ
1 f ððF LF n Þi;j1=2 ; xj1=2 Þ 3 1 þ 13 f ððF LF n Þi;j1=2 ; xj1=2 Þ þ 6 f ððF LF n Þi;j ; xj Þ
3.3. Numerical Scheme
where In steady state, the transient term in the mass balance of the liquid film given by Eq. (1) can be dropped to obtain
dGLF n Pwn P hn ¼ ðD EÞi;j ðCs þ Cv apn Þi;j dx i;j A i;j A i;j n
ð13Þ
ðF LF n Þi;j1=2 ¼ ðF LF n Þi;j1 þ
Dxj f ðF LF n Þi;j1 ; xj1 2
ð19Þ
ðF LF n Þi;j1=2 ¼ ðF LF n Þi;j1 þ
Dxj f ðF LF n Þi;j1=2 ; xj1=2 2
ð20Þ
N. Hammouda et al. / Annals of Nuclear Energy 94 (2016) 313–324
ðF LF n Þi;j ¼ ðF LF n Þi;j1 þ Dxj f ðF LF n Þi;j1=2 ; xj1=2
ð21Þ
The function f represents the RHS of Eq. (17), and is evaluated at the appropriate axial level. For example, at axial level xj-1, the function f for ðF LF n Þi;j1 is
"
00 # q f ððF LF n Þi;j1 ; xj1 Þ ¼ ðD Csn Þi;j1 Ei;j1 wn pdn wi;n hfg i;j1
ð22Þ The
intermediate
ðF LF n Þi;j1=2 ; andðF LF n Þi;j
function
values
of
f
at
ðF LF n Þi;j1=2 ;
must be computed in the order given by
Eqs. (19)–(21). The Runge–Kutta method requires the function f to be calculated four times: once at j 1, twice at j 1/2 and once at j. The implicit dependence of D, E and ðCsn Þ on the liquid film flow ðF LF n Þ is obvious from Eqs. (4)–(22). The computation of these parameters at intermediate values of the liquid film flow requires the use of averages and donor assignments. For example, at axial level j 1/2, the phasic vapour flow, Fv, is taken to be the average of the vapour flow at axial level j 1 and j. For a given subchannel i, the deposition flux Di,j is calculated using the concentration Ci,j given by Eq. (5). The concentration Ci,j, and hence the deposition Di,j is computed in two different ways. The first method computes F LEi;j of Eq. (5) from a mass balance of the total phasic flow of liquid F li;j , obtained from ASSERT-PV’s predictions, and the total liquid film flow, obtained from Eq. (17), from all rod segments facing subchannel i. The second method computes F LEi;j using a mass balance of the local phasic liquid flow, and the liquid film flow for a given rod segment. The latter approach gives the concentration Ci,j as pdn wi;n Pwi;j
C i;j ¼ pdn wi;n Pw i;j
F li;j ðF LF n Þi;j
F l ðF LF n Þi;j i;j
ql
i;j
ð23Þ
F v i;j
þq
317
In contrast, in vertical annular flow there is a wealth of information on droplet concentration distribution, as well as on entrainment and deposition rates. As a consequence, several attempts have been made to extend an annular flow model from vertically upward flow to horizontal flow, where separate descriptions of liquid droplet deposition and entrainment were allowed. Of these attempts, it is shown that a variable deposition flux could lead to significant improvement in annular horizontal flow models (Paras and Karabelas, 1991). In modifying and applying Hewitt and Govan (1990) model in ASSERT, it is implicitly assumed that gravity effects on the liquid film of a rod segment are at least partially accounted for through the gravity terms in the momentum equations and buoyancy force in the drift flux equations of ASSERT. The non-uniformity of droplet concentration distribution in the subchannel, however, is not considered in the current film dryout model of ASSERT. Thus, it is expected that in horizontal flow, disagreement between the dryout model predictions and experiment data will be most significant at low flow rates, where the influence of gravity on droplet distribution may be important. Under the influence of gravity, droplets emitted from the liquid film surface follow a parabolic trajectory, and subsequently may redeposit on the liquid film surface at some downstream location. Also, observations in horizontal annular flow seem to indicate that gravity enhances deposition near the bottom of a tube (James et al., 1987). Other than gravity effect, mechanisms for liquid entrainment and deposition that can influence the liquid film thickness distribution in horizontal flow were proposed by Butterworth (1973). These are: 1. Spreading of the liquid film by interfacial wave actions. 2. Spreading by circumferential shear forces due to secondary gas flow. 3. Spreading by surface tension forces.
v i;j
where P wi;j is the total subchannel wetted perimeter. Both methods for computing Di,j are included in the current coding, and the user is given the choice to select either method for CHF calculations. The difference in the calculated deposition rate Di,j by both methods is small. 3.4. Mixing, cross flow, channel orientation and liquid film spreading In all calculations, the average subchannel properties, already determined by ASSERT, are used as input to the subchannel computations for each axial segment. Therefore, mixing and crossflow effects are already considered for the various subchannels after computations are completed. However, turbulent mixing and void diffusion are accounted for in the framework of the drift flux model of ASSERT, which does not compute mixing effects based on specific inter subchannel droplet mixing, as is usually done in multifield two-fluid models. In horizontal flow, especially when the flow rate is low, gravitational force tends to destroy the axial symmetry of the liquid film in annular flow, resulting in a thicker liquid film at the rod bottom than at the top. This gravity-imposed asymmetry affects not only the liquid distribution in the film but also the entrained-droplets distribution. Therefore, in horizontal flow, it is expected that the circumferential distributions of entrainment and deposition rates in subchannels are generally non-uniform. The non-uniformity of liquid film, entrainment and deposition cause major difficulties in developing reliable theoretical models and conducting reliable experiments. In fact, data on local deposition and entrainment rates in horizontal flows are very limited.
These mechanisms are assumed to act against the liquid drainage force of gravity, and they transfer liquid towards the top of a tube or a rod. None of these mechanisms is considered in the current annular flow model of ASSERT, because they require extensive multidimensional theoretical analyses on the liquid film that are beyond the capability of the current subchannel-based analysis. Many workers (Laurinat et al., 1985; Lin et al., 1985; Fukano and Ousaka, 1989), however, investigated the effects of these mechanisms on the circumferential distribution of film thickness in horizontal flow. Laurinat et al. (1985) and Lin et al. (1985) neglected the effects of surface tension in their analysis, and Fukano and Ousaka (Fukano and Ousaka, 1989) concluded that shear forces due to secondary flow are the most important factor to affect the liquid film distribution in a horizontal annular flow. 3.5. Boundary conditions To initiate calculation, the initial deposition and entrainment rates must be specified at the start of annular flow. Unfortunately, accurate specification of initial conditions is not possible, because a determination of initial conditions at the start of annular flow is very difficult to obtain, especially for diabatic flows. This is one of the major difficulties of applying an annular flow dryout model. The boundary conditions are usually arbitrarily estimated, by assuming that a given fraction of the liquid phase is in the liquid film at a given quality. Normally, calculations are started as early as possible (typically at a quality of 1%), in order to allow the mechanisms of entrainment and deposition to correct themselves for any errors in the initial conditions. Fortunately, for long enough channels, the predictions are relatively insensitive to boundary
318
N. Hammouda et al. / Annals of Nuclear Energy 94 (2016) 313–324
conditions. However, predictions can be significantly affected by the estimation of boundary conditions for short heated lengths and high mass fluxes (Govan et al., 1988). Whalley et al. (1974) started the calculations at the point where the quality is 1%, and with an assumed entrained liquid mass fraction of 99%. Although these initial conditions do not appear to be realistic, especially because annular flow may start at qualities much higher than 1%, they were found to give reasonable predictions (Hoyer, 1998). Various other boundary conditions can be assumed. For instance, Govan et al. (1988) found that the predictions of dryout are more sensitive to initial conditions, and using Hewitt and Govan’s (Hewitt and Govan, 1990) entrainment and deposition correlations, obtained good predictions when the fraction of liquid entrained at the start of annular flow Xlo was specified by
X lo ¼
0:99
if G=Gcrit < 0:2
0:90 if G=Gcrit P 0:2
ð24Þ
where
0:85 x rql 0:5 Gcrit ¼ 1:42 De De
ð25Þ
In Eq. (24), G is the total mass flux, and De is the hydraulic equivalent diameter. In the current model, calculation of the initial entrainment and deposition rates is left to the user to specify. The default calculations of the initial conditions, however, are set to
8 > < X lo at > : X 8 > < X lo at > : X
9 ¼ 0:99 > = if G=Gcrit < 0:2 > ; ¼ 0:10 9 ¼ 0:90 > = if G=Gcrit P 0:2 > ; ¼ 0:0
ð26Þ
At a value of 0.90 for Xlo, we follow the approach of Hoyer (1998), where calculations are started in the bubble region (typically at qualities near 0%) for high flow rates. 4. Assessment against tube data This exercise comprises comparisons of predicted CHF against selected tube data of Bennett et al. (1967) and Becker (1965). These data sets are considered to be reliable and are frequently used for CHF model validation. The data cover the following conditions: Tube diameter: 0.00922–0.01262 m Heated length: 1–5.6 m Pressure: 7–14 MPa Mass flux: 0.19–5.63 Mg m2 s1 Inlet subcooling: 33.2–1150 kJ kg1 Orientation: vertical Heat-flux profile: uniform To check the overall performance of the annular dryout model, the measured and calculated dryout heat fluxes are compared in Fig. 1 for Becker and Bennett et al.’s data. The average and RMS errors between measurements and calculations are also shown in the figure. The predictions show satisfactory variation with a small systematic deviation. Note that the model was iterated to predict CHF at the exit of the tube. This yielded a different CHF and hence a different critical quality from that of the experiment. As such, dryout heat-flux comparisons at the critical qualities of the experiments cannot be made.
Figs. 2–6 show the dependency of the prediction error on the inlet subcooling, pressure, mass flux, diameter, and heated length. Fig. 2 shows that prediction errors are almost independent of inlet subcooling, with the exception of some overprediction at very low subcoolings. The largest relative errors occurred at the lower subcoolings (approximately less than 200 kJ kg1). Fig. 3 shows that the prediction errors appear to be independent of pressure. There is, however, a slight overprediction bias at all pressures, except for the 8 MPa cases, where there is a slight underprediction bias. The largest relative errors are observed for the highest pressure of 14 MPa. Fig. 4 shows that the prediction errors are almost independent of mass flux. There is a slight over-prediction bias for mass fluxes below 3 Mg m2 s1. This is mostly caused by the higher scatter of the data at these conditions, and the shorter test sections from which these data points were obtained. A detailed explanation is provided at the end of this section. Fig. 5 shows that while there is a systematic overprediction bias of the dryout heat flux for the 9.22 mm, 9.4 mm and 10 mm diameter test sections, the dryout heat flux is on average underpredicted for the 12.62 mm test section. However, the limited number of data sets (only two are used here), as well as the small number of data points tested for a given test-section diameter, prohibit any general conclusions. Fig. 6 shows that there is no specific trend in the dryout heat flux errors with heated length. As discussed in Section 3.5, the main difficulty in applying the film dryout model is in determining the initial conditions at the start of annular flow. The parameters that are needed at the start of the annular flow are the initial deposition rate, the initial entrainment rate and the flow mass quality. These are neither measured nor easy to determine, and an ad hoc approach is usually required. For instance, in advance of model verification or validation, a series of calculations are usually performed by changing the mass quality and the initial entrained liquid mass fraction at the start of annular flow until predictions give the best agreement with the data (Govan et al., 1988; Sugawara, 1990). Figs. 7 and 8 show the impact of initial conditions on the predicted dryout for a long tube (5.563 m) and a short tube (1.829 m). For this comparison, three sets of initial conditions were used: 1. the default initial conditions as given by Eq. (26), 2. Xlo = 0.99 at X = 0.01, and 3. Xlo and X adjusted, for each case, to give best agreement with the data. Fig. 7 shows that for the long test section the dryout predictions are insensitive to the initial conditions at low mass-flow rates. As the mass flow increases, however, the predictions become more sensitive to initial conditions. Fig. 8 shows that, for the short test section, dryout predictions are very sensitive to initial conditions at all mass flow rates. In general, especially for short test sections, it is expected that the initial conditions must depend on mass flux, inlet subcooling and pressure. Therefore, it is desirable to develop a relationship that is able to reflect this dependency for a wide range of inlet conditions, to eliminate the need for random guessing. This should be further investigated in future work. Furthermore, in the short test sections, used in the Bennett et al. (1967) and Becker (1965) experiments, not only the phenomena of droplet deposition and entrainment do not have sufficient length to reach equilibrium conditions but also the CHF is not entirely a liquid film dryout type. In such short test sections, to reach CHF conditions, the test section local heat flux, at relatively higher subcooling, must be significantly higher than that in the longer test sections at equivalent conditions. This usually produces favourable conditions for departure from nucleate boiling (DNB) type CHF or
N. Hammouda et al. / Annals of Nuclear Energy 94 (2016) 313–324
Fig. 1. Comparison of predictions vs measured dryout heat flux (Becker and Bennett data).
Fig. 2. Prediction accuracy against inlet subcooling (Becker and Bennett data).
Fig. 3. Prediction accuracy against pressure (Becker and Bennett data).
319
320
N. Hammouda et al. / Annals of Nuclear Energy 94 (2016) 313–324
Fig. 4. Prediction accuracy against mass flux (Becker and Bennett data).
Fig. 5. Prediction accuracy against diameter (Becker and Bennett data).
transition type CHF (CHF that is not quite DNB type or film dryout type). DNB is the CHF point that marks the end of the nucleate boiling regime, where the most widely suggested mechanism for CHF is bubble crowding and coalescence near the wall to form a vapour film. This type of CHF is usually encountered at pressurised light water reactor (PWR) conditions, and requires fundamentally different modelling concepts and assumptions compared to the liquid film dryout of annular flow. The work described herein covers the liquid film dryout in annular flow, as this is the type of CHF encountered in CANDU heavy water pressurised reactors for the postulated design basis events, such as loss of flow and small loss of coolant accidents. Hoyer (1998) suggested using Eqs. (9) and (10) to account for heat-flux-induced entrainment and deposition inhibition through film evaporation. The deposition-inhibition correlation, however, is empirical, and is expected to be valid only for the conditions examined by Hoyer for flow in vertical tubes. The heat-fluxinduced entrainment correlation is valid only for the conditions specified in Section 2.2. Fig. 9 shows the impact of applying Eqs. (9) and (10) on the dryout heat flux. Fig. 9 shows that the
prediction accuracy deteriorates when using these equations, and hence they are not recommended for the ASSERT-PV code.
5. Assessment against bundle data The bundle experimental data sets selected for model prediction comparisons are the Stern Laboratories (SL) 28-element Fortman, 2010, 37-element (Fortman et al., 1997; Fortman, 2011) and 43-element (Dimmick et al., 1999; Leung et al., 2001) bundles. These particular datasets were chosen because (a) the bundle geometries were typical of CANDU designs, (b) the tests included axially non-uniform heat flux profiles, and (c) the data were obtained using the latest experimental techniques, and hence the data are relatively reliable. The test sections of the three bundles used for the comparison are shown in Fig. 10 (Rao et al., 2014b). This section compares the annular flow dryout model’s prediction of the onset of dryout power with the experimental measurements in the SL 28-element, 37-element and 43-element bundles. The data cover the following conditions (Rao et al., 2014b):
N. Hammouda et al. / Annals of Nuclear Energy 94 (2016) 313–324
321
Fig. 6. Prediction accuracy against heated length.
Fig. 7. Impact of initial conditions on the predictions using long test section (Becker and Bennett data).
Rod diameter: 15.3 mm (28e), 13.1 mm (37e), 11.5 and 13.5 mm (43e) Channel length: 12 bundles: 5.76 m (heated), 5.94 m (total) Exit pressure: 5–11 MPa Channel flow rate: 7–29 kg s1 Inlet temperature: 200–290 °C Applied power: 3.25–12.6 MW Orientation: horizontal Axial heat-flux profile: Cosine (28e), downstream skewed (37e, 43e) The detailed experiment data and ASSERT modelling for the bundle tests are described in reference (Rao et al., 2014b). The measured and calculated dryout powers are compared in Fig. 11 for the three types of SL bundles. It is clear that the model predicts most of the data within the ±20% error band. Given the scatter in the experimental data, the model’s results for dryout power are in reasonable agreement with the experiment. The CHF look-up table results are also in excellent agreement with the measurements. ASSERT predicts the CHF using the LW-T-
2005 CHF table (Groeneveld et al., 2005) and its associated correction factors for diameter, quality, gap size, flow-orientation and enhancement effects (Rao et al., 2014b). The annular flow dryout model, however, does not fully account for gravity effect, especially at low flow rates, or CHF enhancement due to appendages. The possible impact of these on dryout power predictions will be discussed in Section 6. Fig. 12–14 show the dependency of the prediction error on the inlet temperature, exit pressure and inlet mass flow rate. Fig. 12 shows a relatively large scattering in prediction errors at all inlet temperatures. Fig. 13 shows that the largest relative errors are observed for the higher pressure of 9–11 MPa. The prediction errors appear to be independent of exit pressure and inlet temperature. Fig. 14 shows that while the predictions are very satisfactory at high mass flow rates, the model clearly overpredicts the data at low mass flow rates. This is partly a consequence of not accounting for gravity effects on droplet concentration distribution in the subchannel at low flow rates, as discussed in section 3.4. This is an area that deserves further investigation in the future, to improve the dryout model’s performance at low-flow rates in horizontal channels.
322
N. Hammouda et al. / Annals of Nuclear Energy 94 (2016) 313–324
Fig. 8. Impact of initial conditions on the predictions using short test section (Becker and Bennett data).
Fig. 9. Effect of heat-flux controlled deposition and entrainment rates on the predictions (Becker and Bennett data).
Fig. 10. Stern Lab test section, axial view (top) and cross-section view (28-element, bottom-left, 37-element, bottom-middle, 43-element, bottom-right) (Dimmick et al., 1999).
N. Hammouda et al. / Annals of Nuclear Energy 94 (2016) 313–324
323
Fig. 11. Comparison predictions vs measured dryout power (Stern Lab data). Fig. 13. Prediction accuracy against exit pressure (Stern Lab data).
Fig. 12. Prediction accuracy against inlet temperature (Stern Lab data).
6. Conclusions and recommendations 6.1. Conclusions 1. An annular flow dryout model, based on that of Whalley et al. (1974) and using the entrainment and deposition correlations of Hewitt and Govan, 1990, was modified for rod bundle geometries and implemented in ASSERT-PV. 2. The ASSERT-PV predictions of dryout power using the mechanistic CHF model have been compared with data from tube and from CANDU-bundle geometries. The overall average and RMS errors of the predictions for the two sets of tube data of water are 0.15% and 5.46%, respectively. The overall average and RMS errors of the predictions for the SL 28-element, 37element and CANFLEX bundle data are 0.9% and 9.9%, respectively.
Fig. 14. Prediction accuracy against flow rate (Stern Lab data).
3. The model was found to be relatively sensitive to the initial conditions, especially at high mass fluxes and short heated length. 4. The model was found to be sensitive to flow orientation (e.g., horizontal flow), due to the increased importance of gravity on liquid film thickness, entrainment rate, deposition rate, and droplet concentration at lower mass-flow rates. 6.2. Recommendations From the analysis presented in this paper, it is recommended that both the annular flow dryout modes and the CHF look-up table be considered as primary methods for predicting CHF in tubes and bundles. More specifically, the CHF look-up table should remain the default method of predicting dryout-type CHF until the liquid film dryout model predictions are improved for low mass-flow rates in horizontal channels.
324
N. Hammouda et al. / Annals of Nuclear Energy 94 (2016) 313–324
The effects of appendages on the liquid film and droplet dynamics in annular flow have not been included in the present dryout model formulation. They are expected to have a significant impact on the CHF, and should be investigated thoroughly. It is well established experimentally that dryout occurs just upstream of the spacer in annular flow for BWR grid-type spacers, and the same occurrence has been confirmed for CANDU-type appendages, such as spacer pads, bearing pads and end plates. It should be noted that work is currently under way to qualify a modified annular flow dryout model that includes the effect of bundle appendages on the balance of the droplet entrainment and deposition in a new version of ASSERT-PV (v3.2.1). This study shows that while the annular flow dryout model predictions are satisfactory at high mass fluxes, the model clearly overpredicts the data at low mass fluxes. This is partly a consequence of not accounting for gravity effects on droplet concentration distribution in the subchannel at low flow rates, as discussed in Section 3.4. This is an area that deserves further investigation in the future, to improve the dryout power prediction at low flow rates in horizontal channels. As discussed in Section 3.5, one of the difficulties in applying the annular flow dryout model is in determining the initial conditions at the start of annular flow, and an ad hoc approach is usually required. In general, it is expected that the initial conditions must depend on mass flux, inlet subcooling and pressure. Therefore, it is desirable that a relationship be developed that can reflect this dependency for a wide range of inlet conditions, to minimise the uncertainty associated with the initial conditions. References Becker, K.M., 1965. An Analytical and Experimental Study of Burnout Conditions in Vertical Round Ducts Report AE-178. Aktiebolaget Atomenergi, Stockholm, Sweden. Bennett, A.W., Hewitt, G.F., Kearsey, H.A., Keeys, R.K.F., 1967. Heat Transfer to Steam-Water Flowing in Uniformly Heated Tubes in which the Critical Heat Flux has been Exceeded UKAEA Report, AERE- R5373. Biasi, L., Clerici, G.C., Garribba, S., Sala, R., Tozzi, A., 1967. Studies on burnout. Part 3: a new correlation for round ducts and uniform heating and its comparison with world data. Energia Nucl. 14, 530–536. Bowring, R.W., 1972. A Simple but Accurate Round Tube Uniform Heat Flux Dryout Correlation Over the Range 0.7–17 MN/m2 UKAEA Report AEEW R789. Butterworth, D., 1973. An Analysis of Film Flow for Horizontal Flow and Condensation in a Horizontal Tube Report AERE-R7575, pp. 1–18. Carver, M.B., Kiteley, J.C., Tahir, A., Banas, A.O., Rowe, D.S., 1990. Simulation of flow and phase distribution in vertical and horizontal bundles using the ASSERT subchannel code. Nucl. Eng. Des. 122, 413–424. Dimmick, G.R., Inch, W.W.R., Jun, J.S., Suk, H.C., Hadaller, G.I., Fortman, R.A., Hayes, R.C., 1999. Full scale water CHF testing of the CANFLEX bundle. In: The 6th International Conference on CANDU Fuel. Canadian Nuclear Society. Fortman, R.A., 2010. PDO and CHF Tests On CANDU 28-Element Fuel Stern Laboratories Report (internal report) SL-195. Fortman, R.A., 2011. Critical Heat Flux and Post-Dryout Tests Using the Reference 37-Element Fuel Simulation in Water Stern Laboratories Report (internal report) SL-212. Fortman, R.A., Hadaller, G.I., Hayes, R.C., Stern, F., 1997. Heat transfer studies with CANDU fuel simulators. In: The 5th International Conference on Nuclear Engineering (ICONE5), Nice, France.
Fukano, T., Ousaka, A., 1989. Prediction of the circumferential distribution of film thickness in horizontal and near-horizontal gas-liquid annular flows. Int. J. Multiphase Flow 15 (3), 403–419. Govan, A.H., Hewitt, G.F., Owen, D.G., Bott, T.R., 1988. An improved CHF modelling code. In: 2nd UK National Heat Transfer Conf., Glasgow. Griffith, P., Pearson, J.F., Lepkowski, R.J., 1977. Critical heat flux during a loss-ofcooling accident. Nucl. Saf. 18. Groeneveld, D.C., Shan, J.Q., Vasic, A.Z., Leung, L.K.H., Durmayaz, A., Yang, J., Cheng, S.C., Tanase, A., 2005. The 2005 CHF look-up table. In: 11th International Topical Meeting on Nuclear Reactor Thermal-Hydraulics (NURETH-11), Avignon, France, 2005 October 2–6, Paper No. 166. Hammouda, N., Carlucci, L.N., Waddington, G.M., 2001. ASSERT-PV IST V3R0: Software Theory Manual AECL Report, RC-2545, FFC-FCT-322. Hewitt, G.F., Govan, A.H., 1990. Phenomenological modelling of non-equilibrium flows with phase change. Int. J. Heat Mass Transfer 33 (2), 229–242. Hoyer, N., 1998. Calculation of dryout and post-dryout heat transfer for tube geometry. Int. J. Multiphase Flow 24 (2), 319–334. Ishii, M., 1975. Thermo-Fluid Dynamic Theory of Two-Phase Flow. Eyrolles, Paris. James, P.W., Wilkes, N.S., Conkie, W., Burns, A., 1987. Developments in the modelling of horizontal annular two-phase flow. Int. J. Multiphase Flow 13 (2), 173–198. Laurinat, J.E., Hanratty, T.J., Jepson, W.P., 1985. Film thickness distribution for gasliquid annular flow in horizontal pipes. Physicochem. Hydrodyn. 6, 179–195. Leung, L.K.H., Jun, J.S., Bullock, D.E., Inch, W.W.R., Suk, H.C., 2001. Dryout power in a CANFLEX bundle string with raised bearing pads. In: The 7th International Conference on CANDU Fuel. Canadian Nuclear Society. Lim, J.C., Weisman, J., 1988. A phenomenologically based prediction of rod-bundle dryout. Nucl. Eng. Des. 105, 363–371. Lin, T.F., Jones, O.C., Lahey, R.T., Block, R.C., Murase, M., 1985. Film thickness measurements and modeling in horizontal annular flows. Physicochem. Hydrodyn. 6, 197–206. Milashenko, V.I., Nigmatulin, B.I., Petukhov, V.V., Tribunkin, N.I., 1989. Burnout and distribution of liquid in evaporative channels of various lengths. Int. J. Multiphase Flow 15, 393–402. Paras, S.V., Karabelas, A.J., 1991. Droplet entrainment and deposition in horizontal annular flow. Int. J. Multiphase Flow 17 (4), 455–468. Rao, Y.F., Cheng, Z., Waddington, G.M., Nava-Dominguez, A., 2014a. ASSERT-PV 3.2: advanced subchannel thermalhydraulics code for CANDU fuel bundles. Nucl. Eng. Des. 275, 69–79. Rao, Y.F., Cheng, Z., Waddington, G.M., 2014b. Assessment of subchannel code ASSERT-PV for prediction of critical heat flux in CANDU bundles. Nucl. Eng. Des. 276, 216–227. Rowe, D.S., Carver, M.B., Kiteley, J.C., Tahir, A., Banas, A.O., 1988. ASSERT-IV Theory Manual AECL Report, CRNL-4221. Rowe, D.S., Carver, M.B., Banas, A.O., Kiteley, J.C., 1993. ASSERT-PV Theoretical Description and Numerical Solution ARD-TD-391, COG-92-371, RC-925. Saito, T., Ishizuka, T., Kimura, J., Kagawa, T., 1983. Application of multi-fluid model in dryout prediction. In: Merilo, M. (Ed.), Thermal Hydraulics of Nuclear Reactors, vol. I. American Nuclear Society, LaGrange Park, IL. Stewart, C.W. et al., 1977. COBRA-IV: The Model and the Method Battelle Pacific Northwest Laboratories Report, BNWL-2214. Sugawara, S., 1990. Droplet deposition and entrainment modeling based on the three-fluid model. Nucl. Eng. Des. 122, 67–84. Vaitekunas, D.A., Tain, R.M., Cheng, S.C., 1996. Assessment of Critical Heat Flux Models AECL Report COG-95-448 (also UO-MCG-TH-95002). Weisman, J., Pei, B.S., 1983. Prediction of critical heat flux in flow boiling at low qualities. Int. J. Heat Mass Transfer 26, 1463–1477. Whalley, P.B., 1977. The calculation of dryout in a rod bundle. Int. J. Mutiphase Flow 3, 501–515. Whalley, P.B., Hutchinson, P., Hewitt, G.F., 1974. The calculation of critical heat flux in forced convection boiling. In: Fifth Int. Heat Transfer Conference, Tokyo, Paper B6.11. Wheeler, C.L. et al., 1976. COBRA-IV-I: An Interim Version of COBRA for ThermalHydraulic Analysis of Rod Bundle Nuclear Fuel Elements and Cores Battelle Pacific Northwest Laboratories Report, BNWL-1962.
Contents lists available at ScienceDirect
Annals of Nuclear Energy journal homepage: www.elsevier.com/locate/anucene
A subchannel based annular flow dryout model Najmeddine Hammouda a,⇑, Zhong Cheng b, Yanfei F. Rao b a b
Canadian Nuclear Safety Commission (CNSC), Ottawa, Ont. K1P 5S9, Canada Canadian Nuclear Laboratories (CNL), Chalk River, Ont. K0J 1J0, Canada
a r t i c l e
i n f o
Article history: Received 3 February 2016 Received in revised form 31 March 2016 Accepted 1 April 2016 Available online 8 April 2016 Keywords: Subchannel thermalhydraulics ASSERT-PV Safety Annular flow Dryout power CHF
a b s t r a c t This paper assesses a popular tube-based mechanistic critical heat flux model (Hewitt and Govan’s annular flow model (based on the model of Whalley et al.), and modifies and implements the model for bundle geometries. It describes the results of the ASSERT subchannel code predictions using the modified model, as applied to a single tube and the 28-element, 37-element and 43-element (CANFLEX) CANDU bundles. A quantitative comparison between the model predictions and experimental data indicates good agreement for a wide range of flow conditions. The comparison has resulted in an overall average error of 0.15% and an overall root-mean-square error of 5.46% with tube data representing annular film dryout type critical heat flux, and in an overall average error of 0.9% and an overall RMS error of 9.9% with Stern Laboratories’ CANDU-bundle data. Crown Copyright Ó 2016 Published by Elsevier Ltd. All rights reserved.
1. Introduction The objective of this study is to assess a popular tube-based annular flow (or liquid film) dryout model (Hewitt and Govan, 1990), and to modify and to implement the model for bundle geometries using the advanced solution of subchannel equation in reactor thermal-hydraulics (ASSERT) subchannel code (Carver et al., 1990; Hammouda et al., 2001; Rao et al., 2014a). ASSERT is a computer code that has been developed to model the subchannel flow and phase distribution in horizontal Canada uranium deuterium (CANDUÒ1) pressurised heavy water reactor fuel channels (Carver et al., 1990; Hammouda et al., 2001; Rao et al., 2014a). The code has been designed to be general enough to accommodate other geometries and orientations. These include single subchannels of different shapes, and multiple subchannels of CANDU PHWR, PWR and BWR designs, in both vertical and horizontal orientations. As well, the code can accommodate a range of fluids, including single- and two-phase heavy water, light water, various Freons, and two-phase air–water. ASSERT (Hammouda et al., 2001; Rao et al., 2014a; Rowe et al., 1993) is based on ASSERT-IV (Rowe et al., 1988), which, in turn originated from the COBRA-IV computer program (Wheeler and et al., 1976; Stewart and et al., 1977); ASSERT has been enhanced to meet the specific requirements for the
⇑ Corresponding author. E-mail addresses: [email protected] (N. Hammouda), zhong.cheng@cnl. ca (Z. Cheng), [email protected] (Y.F. Rao). 1 CANDUÒ is a registered trademark of Atomic Energy of Canada Limited (AECL). http://dx.doi.org/10.1016/j.anucene.2016.04.002 0306-4549/Crown Copyright Ó 2016 Published by Elsevier Ltd. All rights reserved.
thermalhydraulic analysis of two-phase flow in the horizontally oriented CANDU fuel. Currently, empirical equations are the most widely used critical heat flux (CHF) prediction methods, because of a lack of sufficient understanding of the phenomena governing CHF. A large number of CHF correlations are available in the literature; and some of the most popular for vertical flow in tubes are those of Bowring (1972), Griffith et al. (1977) and Biasi et al. (1967). Most CHF correlations were derived for high-pressure steam-water flow in tubes, and have relatively narrow ranges of validity. The number of physically based (or mechanistic) CHF models is very small. In general, mechanistic CHF models are developed for two main flow regimes: (i) the annular flow regime at highquality flows, where the mechanism for CHF is dryout of the liquid film; and (ii) the bubbly flow regime at subcooled or low-quality flows, where the most widely suggested mechanism for CHF is bubble crowding and coalescence near the wall to form a vapour film. Two popular models of CHF are the film dryout model of Whalley et al. (1974) for a rod-centred subchannel formulation, and the departure from nucleate boiling (DNB) model of Weisman and Pei (1983). ASSERT has a built-in capability to systematically determine CHF occurrence, based on its prediction of local subchannel flow and rod heat-flux conditions. It includes a number of options for computing CHF in the form of empirical correlations. Of these options, the CHF table look-up method of Groeneveld et al. (2005) is the most used, because it is applicable over a wide range of conditions. However, the CHF look-up table was derived on the
314
N. Hammouda et al. / Annals of Nuclear Energy 94 (2016) 313–324
Nomenclature A C d D E F Fw ~ g G h hf hg hfg hin, Dhin k P Ph Pw q00 q000 t Vg ~ V Vr X Xlo x
flow area (m2) mass concentration of droplets (kg m3) hydraulic diameter (m) deposition mass flux (kg m2 s1) entrainment mass flux (kg m2 s1) mass-flow rate [kg/s] frictional force in Table 1 gravitational acceleration vector in Table 1 mass flux (kg m2 s1) mixture specific enthalpy in Table 1 specific enthalpy of liquid phase in Table 1 specific enthalpy of vapour phase in Table 1 latent heat of vaporisation [J/kg] specific enthalpy at inlet, subcooling at inlet deposition mass-transfer coefficient (m s1) pressure (Pa) heated perimeter (m) wetted perimeter (m) heat flux (W m2) volumetric heat flux (W m3) time (s) vapour velocity perpendicular to liquid film interface (m s1) velocity vector in Table 1 lateral relative velocity Table 1 quality [–] fraction of liquid entrained at the start of annular flow [–] axial location [m]
basis of 8-mm tube CHF data, and some uncertainties remain in its application to subchannel geometries. Application of the table look-up method for bundle subchannel geometries required use of phenomenological correction models (Rao et al., 2014a). The CHF look-up table is a method of directly representing experimental data obtained for upward flow in 8-mm vertical tubes. Therefore, it may be viewed as a limit to what can be achieved using an empirical approach for CHF in vertical tubes. The parametric trends in CHF of the look-up table are those observed for steam-water flow inside vertical tubes and, for the most part, they have not been ascertained for flow in rod-bundle subchannels. For fluids other than light water, the CHF look-up table requires the use of fluid-to-fluid modelling techniques, which may add uncertainties to the table predictions. However, mechanistic CHF models can be very useful tools to assess parametric trends in subchannel geometries. Other advantages of physically based models are that they can be directly applied to other fluids, can readily be used to extrapolate outside the original data base, and can be used for various geometries and varying heat-flux shapes. Also, in regions where the CHF look-up table values are uncertain, because of scarcity of data, mechanistic CHF models can be used as possible alternatives. It should be noted that this is the first time the annular flow dryout model in ASSERT-PV is published together with an assessment of the model against tube and CANDU-bundle CHF experiments. Recent development and validation of ASSERT-PV CHF models have been focused on the CHF table look-up method plus the associated phenomenological correction models, as summarised in Rao et al. (2014a,b) for ASSERT-PV 3.2. The ASSERT annular flow dryout model was implemented into an earlier version 3.0 (Hammouda et al., 2001) and has stayed unchanged in version 3.2. The Canadian nuclear industry is planning to put more concentrated effort in developing and improving the annular flow dryout model for more accurate subchannel thermalhydraulic analyses.
Greek symbols: void fraction Cs deposition mass-flux suppression due to evaporation (kg m2 s1) Cvap vapour mass-flux generation source function (kg m2 s1) l viscosity (N s m2) q mass density (kg m3) r surface tension (N m1)
a
Subscripts: exp experimental f saturated liquid or liquid phase g saturated vapour or vapour phase l subcooled liquid LE liquid entrained LF liquid film i interface if interface to liquid ig interface to vapour pre predicted sat saturation v vapour w wall
2. Annular flow dryout model The majority of the available annular flow dryout models are based on that of Whalley et al. (1974). The difference from Whalley et al.’s model is in the constitutive relationships representing the mechanisms of entrainment and deposition. The basic Whalley et al. model is a three-field model representing two-phase interactions between vapour, liquid film and entrained droplets. The one-dimensional equations for conservation of mass in annular flow can be derived by considering mass balances written separately for each field,
dðql aLF Þ dGLF Pw Ph þ ¼ ðD EÞ ðCs þ Cv ap Þ dt dx A A
ð1Þ
dðql aLE Þ dGLE Pw Ph þ ¼ ðE DÞ þ Cs dt dx A A
ð2Þ
dðqv aÞ dGv Ph þ ¼ Cv ap dt dx A
ð3Þ
The droplet deposition rate is calculated from
D ¼ kC
ð4Þ
where C is the concentration of droplets (kg m3) in the vapour core, and k (m s1) is the deposition mass transfer coefficient. The concentration, C, is defined as
C¼F
F LE þ qF v
LE
ql
ð5Þ
v
In Eq. (5), FLE and Fv are the droplet and vapour flow rates, respectively.
N. Hammouda et al. / Annals of Nuclear Energy 94 (2016) 313–324
315
2.1. Deposition and entrainment models
2.3. Model of deposition suppression
The development of constitutive relations for entrainment and deposition is a critical step in calculating dryout in annular flow through a mechanistic approach. In fact, it is these relations that distinguish one mechanistic dryout model from another, since the basic conservation equations are the same for any dryout model in the annular flow regime. Hewitt and Govan (1990) derived improved models for deposition and entrainment in annular flow. The new entrainment and deposition correlations successfully predicted a wide range of equilibrium (where the rate of droplet entrainment is equal and opposite to the rate of droplet deposition) and non-equilibrium data. Specifically, the correlations were used to predict equilibrium entrainment and the variation of entrainment in developing flows of air–water and steam-water data (Govan et al., 1988). The distinguishing feature of Hewitt and Govan (1990) model is in its modelling of hydrodynamic non-equilibrium, which includes the fact that a limiting film flow rate (for very thin films) exists, below which no entrainment can take place. The correlation for droplet deposition includes the experimentally observed effect of a decreasing deposition mass-transfer coefficient with increasing droplet concentration. Subsequent studies by Vaitekunas et al. (1996) and Hoyer (1998) showed that Hewitt and Govan’s model of deposition and entrainment gives satisfactory results, when compared to tube data of various fluids. The correlation for the deposition rate coefficient is as follows:
Film evaporation may inhibit droplet deposition. The vapour flowing away from the evaporating film may prevent liquid droplets in the core from reaching the liquid film surface. Hoyer (1998) gives the following correlation to account for the suppressed deposition of droplets due to evaporation:
k¼
8 > <
0:18
qffiffiffiffiffiffi
r qv d
; qC 6 0:3
v
; qC > 0:3
qffiffiffiffiffiffi 0:65 > : 0:083 r C q d q v
v
ð6Þ
v
The entrainment correlation is
8 h i0:316 < E 5:75 105 ðGLF GLFC Þ2 rdqq2l ¼ v Gv : 0
; GLF > GLFC ; GLF 6 GLFC
ð7Þ
where GLFC, the critical film mass flux for the onset of entrainment, is given by
GLFC d
ll
rffiffiffiffiffiffi l ql ¼ exp 5:8504 þ 0:4249 v
ll
qv
ð8Þ
2.2. Model of additional entrainment The majority of entrainment and deposition correlations were derived from data for fully developed annular flow, which mainly describe entrainment and deposition in the absence of heat transfer. However, it has been experimentally shown by Milashenko et al. (1989) that the applied heat flux can have a significant effect on entrainment in annular flow through the boiling process. The mechanism of heat-flux-induced entrainment is the main cause of additional entrainment in diabatic annular flows. The additional entrainment is generated by the boiling process, and is due to the release of bubbles from the liquid film (Hoyer, 1998). The total entrainment (dynamic entrainment Ed, and bubble entrainment, Eq) is given by the correlation of Milashenko et al. (1989) as
Aw Cq ¼ 1:75F LF q00 106
qg qf
!1:3 ð9Þ
where Cq ¼ Ed þ Eq , and Aw is the heated area of the channel. Eq. (9) is applicable in the range of pressures from 5 to 10 MPa, mass fluxes from 1 to 3 Mg m2 s1 and heat fluxes from 0.15 to 4 MW m2. Hoyer (1998) recommends using the maximum value of entrainment obtained from either Eqs. (7) or (9).
Cs ¼ k q C
ð10Þ
where the suppression coefficient kq is given by
kq ¼
qv
0:065ql
Vg
ð11Þ
and where
V g ¼ Cv ap
A
qv Ai
ð12Þ
where Vg is the vapour velocity perpendicular to the liquid film area Ai (Hoyer, 1998). The numerical implementation of the annular film dryout model in ASSERT is described in the next section. 3. Modification and implementation in ASSERT-PV 3.1. Thermal-hydraulic model and solution scheme of ASSERT ASSERT uses the subchannel approach to model single- and two-phase flow in rod bundles. Subchannels are defined by the coolant flow areas between rods, which are bounded by the rods and imaginary planes linking adjacent rod centre lines. Subchannels are divided axially into control volumes, which communicate axially in the same subchannel, and laterally across fictitious boundaries with control volumes in neighbouring subchannels. The code is based on a five-equation, advanced drift-flux thermal-hydraulic model that accounts for the thermal and mechanical non-equilibrium and transverse flow found in CANDU-type fuel channels. The thermal-hydraulic conservation equations used in ASSERT are derived from the two-fluid equations (Ishii, 1975). Heat transfer, friction, and the exchange of mass, momentum and energy between subchannels affect the distribution of flow and enthalpy among individual subchannels. In addition to the conservation equations, ASSERT uses a number of closure relationships to model wall-to-fluid heat and momentum transfer and inter-subchannel mass, momentum and energy exchanges due to turbulent mixing, void drift and buoyancy drift (Carver et al., 1990; Hammouda et al., 2001). The governing equations are given in Table 1. ASSERT uses a staggered grid numerical solution scheme based on a pressure–velocity algorithm (Hammouda et al., 2001). 3.2. The annular flow dryout model: application to rod-bundle geometry For the complex geometry of rod bundles, subchannel analysis is required to calculate CHF, and to appropriately account for thermal non-equilibrium and intersubchannel mixing effects. Several attempts have been made to investigate the occurrence of CHF using subchannel analysis tools in the annular flow regime using the film dryout criterion (Whalley, 1977; Saito et al., 1983; Lim and Weisman, 1988). Whalley (1977) was among the first to apply subchannel analysis to predict film dryout in rod bundles. In his analysis, the bundle was divided into subchannels centred on the rods. On each rod, mass balances on the liquid film, liquid droplets and vapour were applied, accounting for inter subchannel mixing.
316
N. Hammouda et al. / Annals of Nuclear Energy 94 (2016) 313–324
Table 1 Governing equations in ASSERT.a Mixture mass: @q VÞ ¼ 0; þ r ðq~ @t
where q ¼ ðaqÞg þ ðaqÞf ¼ ag qþg þ af qþf VÞ ¼ ðaq~ VÞg þ ðaq~ VÞf ¼ q~ V ðq~ Mixture momentum h i ðaqÞg ðaqÞf ~ ~ V~ Vþ V r V r þ rP ¼ F w þ q~ þ r q~ g q
@ðq~ VÞ @t
Mixture energy: ~ q @h @t þ qV rh þ r
hðaqÞ ðaqÞ g f q
i 00 00 ðhg hf Þ~ V r ¼ q000 w r ½ðaq Þg þ ðaq Þf
Phasic energy (transportive form) Liquid: @h 00 000 VÞf rhf ¼ q000 ðaqÞf @tf þ ðaq~ wf r ðaq Þf þ qif Vapour: ðaqÞg a
@hg @t
00 000 VÞg rhg ¼ q000 þ ðaq~ wg r ðaq Þg þ qig
Here, + = variables to be defined by state relationships and * = variables to be defined by constitutive relationships.
The constitutive relations described in Sections 2.1–2.3 were derived from data obtained for tube geometry. In the analysis presented here, it is assumed that the mechanisms of hydrodynamic entrainment and deposition found in a tube are applicable also to rod-bundle subchannels. This assumption is based on the general notion that any one subchannel within a rod bundle is geometrically equivalent to a round tube, so that similar thermalhydraulic behaviour may be expected (Lim and Weisman, 1988). Since each subchannel of a rod bundle is bounded by a number of rods, it is expected that the liquid film flow is not necessarily the same for all rod sectors facing the same subchannel. The liquid film flow rate in each rod segment depends on the average conditions in the subchannel (void fraction, pressure, liquid droplet concentration, etc.), and on the heat flux of that rod segment. Hydrodynamic entrainment and deposition rates depend on the average subchannel flow conditions and the liquid film flow rate. On the other hand, heatflux-induced entrainment and deposition suppression explicitly depend on the rod heat flux. Therefore, within the same subchannel, mainly the rod having the highest heat flux or the highest radial heat-flux factor determines CHF. In this work, the conventional subchannel approach, in contrast to a rod-centred approach (Whalley, 1977), is used to solve Eqs. (1)–(3). The liquid film is calculated from Eq. (1) using rod heatflux-based correlations (Eq. (9) and Eqs. (10)–(12)) and subchannel-based correlations (Eqs. (5)–(7)). Eqs. (2) and (3) are not needed, since ASSERT calculates the average subchannel flows (mixture and phasic), void fraction and heat flux of each rod. ASSERT calculates mixture flow rates by solving mixture equations of momentum and mass in each subchannel, and provides phasic flow rates using the mixture flow and a drift flux model. ASSERT can provide only the total phasic flow of liquid, since a multi-field two-fluid model is necessary to compute droplet flow and liquid film flow separately. Given the total liquid flow, predicted by ASSERT, and the liquid film flow, predicted from Eq. (1), the entrained droplet flow rate can easily be determined. As a result, Eq. (2) is not needed. Similarly, since the phasic vapour flow is already computed by ASSERT model, Eq. (3) is not used.
The subscripts n, i and j denote rod n, subchannel i, and axial node j. For the purpose of CHF computation, only the perimeters of heated rod segments are considered in the calculation. Therefore, the rod segment wetted perimeter Pw can be replaced by the heated perimeter Ph in Eq. (13) to obtain
dF LF n ¼ ðD E Csn Cv apn Þi;j Phnj dx i;j
ð14Þ
The rod segment heated perimeter is
Phn ¼ p dn wi;n
ð15Þ
where dn is fuel rod n diameter, and wi,n is the fraction of rod n facing subchannel i. The vapour mass-flux generation source Cvap is expressed in terms of the rod segment heat flux as
ðCv apn Þi;j ¼
q00wn hfg
ð16Þ i;j
Inserting Eqs. (15) and (16) into Eq. (14), we obtain
" 00 # q dF LF n ¼ ðD Csn Þi;j Ei;j wn pdn wi;n dx i;j hfg i;j
ð17Þ
Eq. (17) is a first-order differential equation that can be integrated along each rod segment for all subchannels until the first occurrence of a zero liquid-flow rate, FLF, is obtained. For every rod segment facing a subchannel, Eq. (17) is solved as an initial value problem in the x-direction, using a fourth-order Runge–Kutta method. The initial values of all the variables in Eq. (17) are specified at axial location j at the start of the annular flow regime, as described below. Starting at j, the Runge–Kutta method advances the solution of Eq. (17) over an interval Dxj. Using the Runge–Kutta formula, Eq. (17) can be expressed as
dF LF n dx
i;j
¼
ðF LF F LF j
Þ
j1 n;i
Dxj
¼ 16 f ððF LF n Þi;j1 ; xj1 Þþ ð18Þ
1 f ððF LF n Þi;j1=2 ; xj1=2 Þ 3 1 þ 13 f ððF LF n Þi;j1=2 ; xj1=2 Þ þ 6 f ððF LF n Þi;j ; xj Þ
3.3. Numerical Scheme
where In steady state, the transient term in the mass balance of the liquid film given by Eq. (1) can be dropped to obtain
dGLF n Pwn P hn ¼ ðD EÞi;j ðCs þ Cv apn Þi;j dx i;j A i;j A i;j n
ð13Þ
ðF LF n Þi;j1=2 ¼ ðF LF n Þi;j1 þ
Dxj f ðF LF n Þi;j1 ; xj1 2
ð19Þ
ðF LF n Þi;j1=2 ¼ ðF LF n Þi;j1 þ
Dxj f ðF LF n Þi;j1=2 ; xj1=2 2
ð20Þ
N. Hammouda et al. / Annals of Nuclear Energy 94 (2016) 313–324
ðF LF n Þi;j ¼ ðF LF n Þi;j1 þ Dxj f ðF LF n Þi;j1=2 ; xj1=2
ð21Þ
The function f represents the RHS of Eq. (17), and is evaluated at the appropriate axial level. For example, at axial level xj-1, the function f for ðF LF n Þi;j1 is
"
00 # q f ððF LF n Þi;j1 ; xj1 Þ ¼ ðD Csn Þi;j1 Ei;j1 wn pdn wi;n hfg i;j1
ð22Þ The
intermediate
ðF LF n Þi;j1=2 ; andðF LF n Þi;j
function
values
of
f
at
ðF LF n Þi;j1=2 ;
must be computed in the order given by
Eqs. (19)–(21). The Runge–Kutta method requires the function f to be calculated four times: once at j 1, twice at j 1/2 and once at j. The implicit dependence of D, E and ðCsn Þ on the liquid film flow ðF LF n Þ is obvious from Eqs. (4)–(22). The computation of these parameters at intermediate values of the liquid film flow requires the use of averages and donor assignments. For example, at axial level j 1/2, the phasic vapour flow, Fv, is taken to be the average of the vapour flow at axial level j 1 and j. For a given subchannel i, the deposition flux Di,j is calculated using the concentration Ci,j given by Eq. (5). The concentration Ci,j, and hence the deposition Di,j is computed in two different ways. The first method computes F LEi;j of Eq. (5) from a mass balance of the total phasic flow of liquid F li;j , obtained from ASSERT-PV’s predictions, and the total liquid film flow, obtained from Eq. (17), from all rod segments facing subchannel i. The second method computes F LEi;j using a mass balance of the local phasic liquid flow, and the liquid film flow for a given rod segment. The latter approach gives the concentration Ci,j as pdn wi;n Pwi;j
C i;j ¼ pdn wi;n Pw i;j
F li;j ðF LF n Þi;j
F l ðF LF n Þi;j i;j
ql
i;j
ð23Þ
F v i;j
þq
317
In contrast, in vertical annular flow there is a wealth of information on droplet concentration distribution, as well as on entrainment and deposition rates. As a consequence, several attempts have been made to extend an annular flow model from vertically upward flow to horizontal flow, where separate descriptions of liquid droplet deposition and entrainment were allowed. Of these attempts, it is shown that a variable deposition flux could lead to significant improvement in annular horizontal flow models (Paras and Karabelas, 1991). In modifying and applying Hewitt and Govan (1990) model in ASSERT, it is implicitly assumed that gravity effects on the liquid film of a rod segment are at least partially accounted for through the gravity terms in the momentum equations and buoyancy force in the drift flux equations of ASSERT. The non-uniformity of droplet concentration distribution in the subchannel, however, is not considered in the current film dryout model of ASSERT. Thus, it is expected that in horizontal flow, disagreement between the dryout model predictions and experiment data will be most significant at low flow rates, where the influence of gravity on droplet distribution may be important. Under the influence of gravity, droplets emitted from the liquid film surface follow a parabolic trajectory, and subsequently may redeposit on the liquid film surface at some downstream location. Also, observations in horizontal annular flow seem to indicate that gravity enhances deposition near the bottom of a tube (James et al., 1987). Other than gravity effect, mechanisms for liquid entrainment and deposition that can influence the liquid film thickness distribution in horizontal flow were proposed by Butterworth (1973). These are: 1. Spreading of the liquid film by interfacial wave actions. 2. Spreading by circumferential shear forces due to secondary gas flow. 3. Spreading by surface tension forces.
v i;j
where P wi;j is the total subchannel wetted perimeter. Both methods for computing Di,j are included in the current coding, and the user is given the choice to select either method for CHF calculations. The difference in the calculated deposition rate Di,j by both methods is small. 3.4. Mixing, cross flow, channel orientation and liquid film spreading In all calculations, the average subchannel properties, already determined by ASSERT, are used as input to the subchannel computations for each axial segment. Therefore, mixing and crossflow effects are already considered for the various subchannels after computations are completed. However, turbulent mixing and void diffusion are accounted for in the framework of the drift flux model of ASSERT, which does not compute mixing effects based on specific inter subchannel droplet mixing, as is usually done in multifield two-fluid models. In horizontal flow, especially when the flow rate is low, gravitational force tends to destroy the axial symmetry of the liquid film in annular flow, resulting in a thicker liquid film at the rod bottom than at the top. This gravity-imposed asymmetry affects not only the liquid distribution in the film but also the entrained-droplets distribution. Therefore, in horizontal flow, it is expected that the circumferential distributions of entrainment and deposition rates in subchannels are generally non-uniform. The non-uniformity of liquid film, entrainment and deposition cause major difficulties in developing reliable theoretical models and conducting reliable experiments. In fact, data on local deposition and entrainment rates in horizontal flows are very limited.
These mechanisms are assumed to act against the liquid drainage force of gravity, and they transfer liquid towards the top of a tube or a rod. None of these mechanisms is considered in the current annular flow model of ASSERT, because they require extensive multidimensional theoretical analyses on the liquid film that are beyond the capability of the current subchannel-based analysis. Many workers (Laurinat et al., 1985; Lin et al., 1985; Fukano and Ousaka, 1989), however, investigated the effects of these mechanisms on the circumferential distribution of film thickness in horizontal flow. Laurinat et al. (1985) and Lin et al. (1985) neglected the effects of surface tension in their analysis, and Fukano and Ousaka (Fukano and Ousaka, 1989) concluded that shear forces due to secondary flow are the most important factor to affect the liquid film distribution in a horizontal annular flow. 3.5. Boundary conditions To initiate calculation, the initial deposition and entrainment rates must be specified at the start of annular flow. Unfortunately, accurate specification of initial conditions is not possible, because a determination of initial conditions at the start of annular flow is very difficult to obtain, especially for diabatic flows. This is one of the major difficulties of applying an annular flow dryout model. The boundary conditions are usually arbitrarily estimated, by assuming that a given fraction of the liquid phase is in the liquid film at a given quality. Normally, calculations are started as early as possible (typically at a quality of 1%), in order to allow the mechanisms of entrainment and deposition to correct themselves for any errors in the initial conditions. Fortunately, for long enough channels, the predictions are relatively insensitive to boundary
318
N. Hammouda et al. / Annals of Nuclear Energy 94 (2016) 313–324
conditions. However, predictions can be significantly affected by the estimation of boundary conditions for short heated lengths and high mass fluxes (Govan et al., 1988). Whalley et al. (1974) started the calculations at the point where the quality is 1%, and with an assumed entrained liquid mass fraction of 99%. Although these initial conditions do not appear to be realistic, especially because annular flow may start at qualities much higher than 1%, they were found to give reasonable predictions (Hoyer, 1998). Various other boundary conditions can be assumed. For instance, Govan et al. (1988) found that the predictions of dryout are more sensitive to initial conditions, and using Hewitt and Govan’s (Hewitt and Govan, 1990) entrainment and deposition correlations, obtained good predictions when the fraction of liquid entrained at the start of annular flow Xlo was specified by
X lo ¼
0:99
if G=Gcrit < 0:2
0:90 if G=Gcrit P 0:2
ð24Þ
where
0:85 x rql 0:5 Gcrit ¼ 1:42 De De
ð25Þ
In Eq. (24), G is the total mass flux, and De is the hydraulic equivalent diameter. In the current model, calculation of the initial entrainment and deposition rates is left to the user to specify. The default calculations of the initial conditions, however, are set to
8 > < X lo at > : X 8 > < X lo at > : X
9 ¼ 0:99 > = if G=Gcrit < 0:2 > ; ¼ 0:10 9 ¼ 0:90 > = if G=Gcrit P 0:2 > ; ¼ 0:0
ð26Þ
At a value of 0.90 for Xlo, we follow the approach of Hoyer (1998), where calculations are started in the bubble region (typically at qualities near 0%) for high flow rates. 4. Assessment against tube data This exercise comprises comparisons of predicted CHF against selected tube data of Bennett et al. (1967) and Becker (1965). These data sets are considered to be reliable and are frequently used for CHF model validation. The data cover the following conditions: Tube diameter: 0.00922–0.01262 m Heated length: 1–5.6 m Pressure: 7–14 MPa Mass flux: 0.19–5.63 Mg m2 s1 Inlet subcooling: 33.2–1150 kJ kg1 Orientation: vertical Heat-flux profile: uniform To check the overall performance of the annular dryout model, the measured and calculated dryout heat fluxes are compared in Fig. 1 for Becker and Bennett et al.’s data. The average and RMS errors between measurements and calculations are also shown in the figure. The predictions show satisfactory variation with a small systematic deviation. Note that the model was iterated to predict CHF at the exit of the tube. This yielded a different CHF and hence a different critical quality from that of the experiment. As such, dryout heat-flux comparisons at the critical qualities of the experiments cannot be made.
Figs. 2–6 show the dependency of the prediction error on the inlet subcooling, pressure, mass flux, diameter, and heated length. Fig. 2 shows that prediction errors are almost independent of inlet subcooling, with the exception of some overprediction at very low subcoolings. The largest relative errors occurred at the lower subcoolings (approximately less than 200 kJ kg1). Fig. 3 shows that the prediction errors appear to be independent of pressure. There is, however, a slight overprediction bias at all pressures, except for the 8 MPa cases, where there is a slight underprediction bias. The largest relative errors are observed for the highest pressure of 14 MPa. Fig. 4 shows that the prediction errors are almost independent of mass flux. There is a slight over-prediction bias for mass fluxes below 3 Mg m2 s1. This is mostly caused by the higher scatter of the data at these conditions, and the shorter test sections from which these data points were obtained. A detailed explanation is provided at the end of this section. Fig. 5 shows that while there is a systematic overprediction bias of the dryout heat flux for the 9.22 mm, 9.4 mm and 10 mm diameter test sections, the dryout heat flux is on average underpredicted for the 12.62 mm test section. However, the limited number of data sets (only two are used here), as well as the small number of data points tested for a given test-section diameter, prohibit any general conclusions. Fig. 6 shows that there is no specific trend in the dryout heat flux errors with heated length. As discussed in Section 3.5, the main difficulty in applying the film dryout model is in determining the initial conditions at the start of annular flow. The parameters that are needed at the start of the annular flow are the initial deposition rate, the initial entrainment rate and the flow mass quality. These are neither measured nor easy to determine, and an ad hoc approach is usually required. For instance, in advance of model verification or validation, a series of calculations are usually performed by changing the mass quality and the initial entrained liquid mass fraction at the start of annular flow until predictions give the best agreement with the data (Govan et al., 1988; Sugawara, 1990). Figs. 7 and 8 show the impact of initial conditions on the predicted dryout for a long tube (5.563 m) and a short tube (1.829 m). For this comparison, three sets of initial conditions were used: 1. the default initial conditions as given by Eq. (26), 2. Xlo = 0.99 at X = 0.01, and 3. Xlo and X adjusted, for each case, to give best agreement with the data. Fig. 7 shows that for the long test section the dryout predictions are insensitive to the initial conditions at low mass-flow rates. As the mass flow increases, however, the predictions become more sensitive to initial conditions. Fig. 8 shows that, for the short test section, dryout predictions are very sensitive to initial conditions at all mass flow rates. In general, especially for short test sections, it is expected that the initial conditions must depend on mass flux, inlet subcooling and pressure. Therefore, it is desirable to develop a relationship that is able to reflect this dependency for a wide range of inlet conditions, to eliminate the need for random guessing. This should be further investigated in future work. Furthermore, in the short test sections, used in the Bennett et al. (1967) and Becker (1965) experiments, not only the phenomena of droplet deposition and entrainment do not have sufficient length to reach equilibrium conditions but also the CHF is not entirely a liquid film dryout type. In such short test sections, to reach CHF conditions, the test section local heat flux, at relatively higher subcooling, must be significantly higher than that in the longer test sections at equivalent conditions. This usually produces favourable conditions for departure from nucleate boiling (DNB) type CHF or
N. Hammouda et al. / Annals of Nuclear Energy 94 (2016) 313–324
Fig. 1. Comparison of predictions vs measured dryout heat flux (Becker and Bennett data).
Fig. 2. Prediction accuracy against inlet subcooling (Becker and Bennett data).
Fig. 3. Prediction accuracy against pressure (Becker and Bennett data).
319
320
N. Hammouda et al. / Annals of Nuclear Energy 94 (2016) 313–324
Fig. 4. Prediction accuracy against mass flux (Becker and Bennett data).
Fig. 5. Prediction accuracy against diameter (Becker and Bennett data).
transition type CHF (CHF that is not quite DNB type or film dryout type). DNB is the CHF point that marks the end of the nucleate boiling regime, where the most widely suggested mechanism for CHF is bubble crowding and coalescence near the wall to form a vapour film. This type of CHF is usually encountered at pressurised light water reactor (PWR) conditions, and requires fundamentally different modelling concepts and assumptions compared to the liquid film dryout of annular flow. The work described herein covers the liquid film dryout in annular flow, as this is the type of CHF encountered in CANDU heavy water pressurised reactors for the postulated design basis events, such as loss of flow and small loss of coolant accidents. Hoyer (1998) suggested using Eqs. (9) and (10) to account for heat-flux-induced entrainment and deposition inhibition through film evaporation. The deposition-inhibition correlation, however, is empirical, and is expected to be valid only for the conditions examined by Hoyer for flow in vertical tubes. The heat-fluxinduced entrainment correlation is valid only for the conditions specified in Section 2.2. Fig. 9 shows the impact of applying Eqs. (9) and (10) on the dryout heat flux. Fig. 9 shows that the
prediction accuracy deteriorates when using these equations, and hence they are not recommended for the ASSERT-PV code.
5. Assessment against bundle data The bundle experimental data sets selected for model prediction comparisons are the Stern Laboratories (SL) 28-element Fortman, 2010, 37-element (Fortman et al., 1997; Fortman, 2011) and 43-element (Dimmick et al., 1999; Leung et al., 2001) bundles. These particular datasets were chosen because (a) the bundle geometries were typical of CANDU designs, (b) the tests included axially non-uniform heat flux profiles, and (c) the data were obtained using the latest experimental techniques, and hence the data are relatively reliable. The test sections of the three bundles used for the comparison are shown in Fig. 10 (Rao et al., 2014b). This section compares the annular flow dryout model’s prediction of the onset of dryout power with the experimental measurements in the SL 28-element, 37-element and 43-element bundles. The data cover the following conditions (Rao et al., 2014b):
N. Hammouda et al. / Annals of Nuclear Energy 94 (2016) 313–324
321
Fig. 6. Prediction accuracy against heated length.
Fig. 7. Impact of initial conditions on the predictions using long test section (Becker and Bennett data).
Rod diameter: 15.3 mm (28e), 13.1 mm (37e), 11.5 and 13.5 mm (43e) Channel length: 12 bundles: 5.76 m (heated), 5.94 m (total) Exit pressure: 5–11 MPa Channel flow rate: 7–29 kg s1 Inlet temperature: 200–290 °C Applied power: 3.25–12.6 MW Orientation: horizontal Axial heat-flux profile: Cosine (28e), downstream skewed (37e, 43e) The detailed experiment data and ASSERT modelling for the bundle tests are described in reference (Rao et al., 2014b). The measured and calculated dryout powers are compared in Fig. 11 for the three types of SL bundles. It is clear that the model predicts most of the data within the ±20% error band. Given the scatter in the experimental data, the model’s results for dryout power are in reasonable agreement with the experiment. The CHF look-up table results are also in excellent agreement with the measurements. ASSERT predicts the CHF using the LW-T-
2005 CHF table (Groeneveld et al., 2005) and its associated correction factors for diameter, quality, gap size, flow-orientation and enhancement effects (Rao et al., 2014b). The annular flow dryout model, however, does not fully account for gravity effect, especially at low flow rates, or CHF enhancement due to appendages. The possible impact of these on dryout power predictions will be discussed in Section 6. Fig. 12–14 show the dependency of the prediction error on the inlet temperature, exit pressure and inlet mass flow rate. Fig. 12 shows a relatively large scattering in prediction errors at all inlet temperatures. Fig. 13 shows that the largest relative errors are observed for the higher pressure of 9–11 MPa. The prediction errors appear to be independent of exit pressure and inlet temperature. Fig. 14 shows that while the predictions are very satisfactory at high mass flow rates, the model clearly overpredicts the data at low mass flow rates. This is partly a consequence of not accounting for gravity effects on droplet concentration distribution in the subchannel at low flow rates, as discussed in section 3.4. This is an area that deserves further investigation in the future, to improve the dryout model’s performance at low-flow rates in horizontal channels.
322
N. Hammouda et al. / Annals of Nuclear Energy 94 (2016) 313–324
Fig. 8. Impact of initial conditions on the predictions using short test section (Becker and Bennett data).
Fig. 9. Effect of heat-flux controlled deposition and entrainment rates on the predictions (Becker and Bennett data).
Fig. 10. Stern Lab test section, axial view (top) and cross-section view (28-element, bottom-left, 37-element, bottom-middle, 43-element, bottom-right) (Dimmick et al., 1999).
N. Hammouda et al. / Annals of Nuclear Energy 94 (2016) 313–324
323
Fig. 11. Comparison predictions vs measured dryout power (Stern Lab data). Fig. 13. Prediction accuracy against exit pressure (Stern Lab data).
Fig. 12. Prediction accuracy against inlet temperature (Stern Lab data).
6. Conclusions and recommendations 6.1. Conclusions 1. An annular flow dryout model, based on that of Whalley et al. (1974) and using the entrainment and deposition correlations of Hewitt and Govan, 1990, was modified for rod bundle geometries and implemented in ASSERT-PV. 2. The ASSERT-PV predictions of dryout power using the mechanistic CHF model have been compared with data from tube and from CANDU-bundle geometries. The overall average and RMS errors of the predictions for the two sets of tube data of water are 0.15% and 5.46%, respectively. The overall average and RMS errors of the predictions for the SL 28-element, 37element and CANFLEX bundle data are 0.9% and 9.9%, respectively.
Fig. 14. Prediction accuracy against flow rate (Stern Lab data).
3. The model was found to be relatively sensitive to the initial conditions, especially at high mass fluxes and short heated length. 4. The model was found to be sensitive to flow orientation (e.g., horizontal flow), due to the increased importance of gravity on liquid film thickness, entrainment rate, deposition rate, and droplet concentration at lower mass-flow rates. 6.2. Recommendations From the analysis presented in this paper, it is recommended that both the annular flow dryout modes and the CHF look-up table be considered as primary methods for predicting CHF in tubes and bundles. More specifically, the CHF look-up table should remain the default method of predicting dryout-type CHF until the liquid film dryout model predictions are improved for low mass-flow rates in horizontal channels.
324
N. Hammouda et al. / Annals of Nuclear Energy 94 (2016) 313–324
The effects of appendages on the liquid film and droplet dynamics in annular flow have not been included in the present dryout model formulation. They are expected to have a significant impact on the CHF, and should be investigated thoroughly. It is well established experimentally that dryout occurs just upstream of the spacer in annular flow for BWR grid-type spacers, and the same occurrence has been confirmed for CANDU-type appendages, such as spacer pads, bearing pads and end plates. It should be noted that work is currently under way to qualify a modified annular flow dryout model that includes the effect of bundle appendages on the balance of the droplet entrainment and deposition in a new version of ASSERT-PV (v3.2.1). This study shows that while the annular flow dryout model predictions are satisfactory at high mass fluxes, the model clearly overpredicts the data at low mass fluxes. This is partly a consequence of not accounting for gravity effects on droplet concentration distribution in the subchannel at low flow rates, as discussed in Section 3.4. This is an area that deserves further investigation in the future, to improve the dryout power prediction at low flow rates in horizontal channels. As discussed in Section 3.5, one of the difficulties in applying the annular flow dryout model is in determining the initial conditions at the start of annular flow, and an ad hoc approach is usually required. In general, it is expected that the initial conditions must depend on mass flux, inlet subcooling and pressure. Therefore, it is desirable that a relationship be developed that can reflect this dependency for a wide range of inlet conditions, to minimise the uncertainty associated with the initial conditions. References Becker, K.M., 1965. An Analytical and Experimental Study of Burnout Conditions in Vertical Round Ducts Report AE-178. Aktiebolaget Atomenergi, Stockholm, Sweden. Bennett, A.W., Hewitt, G.F., Kearsey, H.A., Keeys, R.K.F., 1967. Heat Transfer to Steam-Water Flowing in Uniformly Heated Tubes in which the Critical Heat Flux has been Exceeded UKAEA Report, AERE- R5373. Biasi, L., Clerici, G.C., Garribba, S., Sala, R., Tozzi, A., 1967. Studies on burnout. Part 3: a new correlation for round ducts and uniform heating and its comparison with world data. Energia Nucl. 14, 530–536. Bowring, R.W., 1972. A Simple but Accurate Round Tube Uniform Heat Flux Dryout Correlation Over the Range 0.7–17 MN/m2 UKAEA Report AEEW R789. Butterworth, D., 1973. An Analysis of Film Flow for Horizontal Flow and Condensation in a Horizontal Tube Report AERE-R7575, pp. 1–18. Carver, M.B., Kiteley, J.C., Tahir, A., Banas, A.O., Rowe, D.S., 1990. Simulation of flow and phase distribution in vertical and horizontal bundles using the ASSERT subchannel code. Nucl. Eng. Des. 122, 413–424. Dimmick, G.R., Inch, W.W.R., Jun, J.S., Suk, H.C., Hadaller, G.I., Fortman, R.A., Hayes, R.C., 1999. Full scale water CHF testing of the CANFLEX bundle. In: The 6th International Conference on CANDU Fuel. Canadian Nuclear Society. Fortman, R.A., 2010. PDO and CHF Tests On CANDU 28-Element Fuel Stern Laboratories Report (internal report) SL-195. Fortman, R.A., 2011. Critical Heat Flux and Post-Dryout Tests Using the Reference 37-Element Fuel Simulation in Water Stern Laboratories Report (internal report) SL-212. Fortman, R.A., Hadaller, G.I., Hayes, R.C., Stern, F., 1997. Heat transfer studies with CANDU fuel simulators. In: The 5th International Conference on Nuclear Engineering (ICONE5), Nice, France.
Fukano, T., Ousaka, A., 1989. Prediction of the circumferential distribution of film thickness in horizontal and near-horizontal gas-liquid annular flows. Int. J. Multiphase Flow 15 (3), 403–419. Govan, A.H., Hewitt, G.F., Owen, D.G., Bott, T.R., 1988. An improved CHF modelling code. In: 2nd UK National Heat Transfer Conf., Glasgow. Griffith, P., Pearson, J.F., Lepkowski, R.J., 1977. Critical heat flux during a loss-ofcooling accident. Nucl. Saf. 18. Groeneveld, D.C., Shan, J.Q., Vasic, A.Z., Leung, L.K.H., Durmayaz, A., Yang, J., Cheng, S.C., Tanase, A., 2005. The 2005 CHF look-up table. In: 11th International Topical Meeting on Nuclear Reactor Thermal-Hydraulics (NURETH-11), Avignon, France, 2005 October 2–6, Paper No. 166. Hammouda, N., Carlucci, L.N., Waddington, G.M., 2001. ASSERT-PV IST V3R0: Software Theory Manual AECL Report, RC-2545, FFC-FCT-322. Hewitt, G.F., Govan, A.H., 1990. Phenomenological modelling of non-equilibrium flows with phase change. Int. J. Heat Mass Transfer 33 (2), 229–242. Hoyer, N., 1998. Calculation of dryout and post-dryout heat transfer for tube geometry. Int. J. Multiphase Flow 24 (2), 319–334. Ishii, M., 1975. Thermo-Fluid Dynamic Theory of Two-Phase Flow. Eyrolles, Paris. James, P.W., Wilkes, N.S., Conkie, W., Burns, A., 1987. Developments in the modelling of horizontal annular two-phase flow. Int. J. Multiphase Flow 13 (2), 173–198. Laurinat, J.E., Hanratty, T.J., Jepson, W.P., 1985. Film thickness distribution for gasliquid annular flow in horizontal pipes. Physicochem. Hydrodyn. 6, 179–195. Leung, L.K.H., Jun, J.S., Bullock, D.E., Inch, W.W.R., Suk, H.C., 2001. Dryout power in a CANFLEX bundle string with raised bearing pads. In: The 7th International Conference on CANDU Fuel. Canadian Nuclear Society. Lim, J.C., Weisman, J., 1988. A phenomenologically based prediction of rod-bundle dryout. Nucl. Eng. Des. 105, 363–371. Lin, T.F., Jones, O.C., Lahey, R.T., Block, R.C., Murase, M., 1985. Film thickness measurements and modeling in horizontal annular flows. Physicochem. Hydrodyn. 6, 197–206. Milashenko, V.I., Nigmatulin, B.I., Petukhov, V.V., Tribunkin, N.I., 1989. Burnout and distribution of liquid in evaporative channels of various lengths. Int. J. Multiphase Flow 15, 393–402. Paras, S.V., Karabelas, A.J., 1991. Droplet entrainment and deposition in horizontal annular flow. Int. J. Multiphase Flow 17 (4), 455–468. Rao, Y.F., Cheng, Z., Waddington, G.M., Nava-Dominguez, A., 2014a. ASSERT-PV 3.2: advanced subchannel thermalhydraulics code for CANDU fuel bundles. Nucl. Eng. Des. 275, 69–79. Rao, Y.F., Cheng, Z., Waddington, G.M., 2014b. Assessment of subchannel code ASSERT-PV for prediction of critical heat flux in CANDU bundles. Nucl. Eng. Des. 276, 216–227. Rowe, D.S., Carver, M.B., Kiteley, J.C., Tahir, A., Banas, A.O., 1988. ASSERT-IV Theory Manual AECL Report, CRNL-4221. Rowe, D.S., Carver, M.B., Banas, A.O., Kiteley, J.C., 1993. ASSERT-PV Theoretical Description and Numerical Solution ARD-TD-391, COG-92-371, RC-925. Saito, T., Ishizuka, T., Kimura, J., Kagawa, T., 1983. Application of multi-fluid model in dryout prediction. In: Merilo, M. (Ed.), Thermal Hydraulics of Nuclear Reactors, vol. I. American Nuclear Society, LaGrange Park, IL. Stewart, C.W. et al., 1977. COBRA-IV: The Model and the Method Battelle Pacific Northwest Laboratories Report, BNWL-2214. Sugawara, S., 1990. Droplet deposition and entrainment modeling based on the three-fluid model. Nucl. Eng. Des. 122, 67–84. Vaitekunas, D.A., Tain, R.M., Cheng, S.C., 1996. Assessment of Critical Heat Flux Models AECL Report COG-95-448 (also UO-MCG-TH-95002). Weisman, J., Pei, B.S., 1983. Prediction of critical heat flux in flow boiling at low qualities. Int. J. Heat Mass Transfer 26, 1463–1477. Whalley, P.B., 1977. The calculation of dryout in a rod bundle. Int. J. Mutiphase Flow 3, 501–515. Whalley, P.B., Hutchinson, P., Hewitt, G.F., 1974. The calculation of critical heat flux in forced convection boiling. In: Fifth Int. Heat Transfer Conference, Tokyo, Paper B6.11. Wheeler, C.L. et al., 1976. COBRA-IV-I: An Interim Version of COBRA for ThermalHydraulic Analysis of Rod Bundle Nuclear Fuel Elements and Cores Battelle Pacific Northwest Laboratories Report, BNWL-1962.