- Email: [email protected]
Simulations on LN2–VN2 flooding phenomenon in inclined tubes using a modified AIAD model
Cryogenics 97 (2019) 100–108
Contents lists available at ScienceDirect
Cryogenics journal homepage: www.elsevier.com/locate/cryogenics
Research paper
Simulations on LN2–VN2 flooding phenomenon in inclined tubes using a modified AIAD model
T
Liu Yua, Hong Chenb, Xu Gaob, Rui Zhoua, Huangjun Xiea, Limin Qiua, Xiaobin Zhanga,
⁎
a
Institute of Refrigeration and Cryogenics, Zhejiang University, The Key Laboratory of Refrigeration and Cryogenic Technology of Zhejiang Province, Hangzhou 310027, China b State Key Laboratory of Technologies in Space Cryogenic Propellants, Beijing 10028, China
ARTICLE INFO
ABSTRACT
Keywords: Flooding LN2 AIAD model CFD Inclined tube
Simulations basing on the Eulerian multiphase model and a modified algebraic interface area density (AIAD) model are carried out to study the flooding mechanisms with liquid nitrogen (LN2) and its vapor (VN2) in an inclined pipe. The original AIAD model predicts well for the water-air flooding velocity, which assumes three flow patterns in the whole hydraulic process, including the free surface flow, bubbly flow and droplet flow. Much different flooding phenomenon for LN2-VN2 is experimentally observed by the present authors, which indicates that some parameters have to be modified for LN2-VN2 occasion. The flooding onset is mainly identified by the sharp increase of pressure drop in this paper. As expected, the modified AIAD model predicts LN2-VN2 flooding velocity with a higher accordance with experiment than before. Both the experiments and calculations show that the flooding process of LN2-VN2 goes through three stages from stratified flow to wavy flow, and finally to mist flow with the appearance of spray entrainment, which is different from the final slug flow in water-air cases. Influences of geometric factors, such as tube length and diameter, are also investigated in this research.
1. Introduction Gas-liquid countercurrent flow is always hindered by the interfacial friction, and the greater is the relative velocity, the greater is the frictional force. When the geometry and working fluid pair are determined, there is a critical relative velocity, beyond which the flow gets into an unstable state. This critical velocity is the point at which flooding occurs [1]. Countercurrent two-phase flow is frequently encountered in a lot of industrial situations. It exists widely in nuclear reactors under abnormal operation, gravity heat pipes and distillation towers. Flooding encountered in industrial occasion often deteriorates the regular operation of devices [2]. For example, Flooding in gravity heat pipe tends to dry out the pipe, which means invalidation of equipment. The gas/ liquid system complete the heat and mass transfer through countercurrent contact in the structured packed distillation tower. If the gas velocity is too large, the liquid cannot continue to flow downward due to flooding. This will destroy the mass transfer and reduce the efficiency. So, the flooding velocity between phases is a vital parameter that we need to determine to maintain a safe operating environment. Large amount of experiments have been carried out in scaled models by Valléé [3,4], Kusunoki [5] and Issa [6,7], to investigate the flooding behaviors in the hot leg of a pressurized water reactor. Besides, Zapke
⁎
and Kröger [8,9], Suzuki [10], Jayanti and Hewitt [11] also studied mechanisms of flooding in vertical, inclined and horizontal tubes. Empirical correlations were concluded basing on these experimental data. However, formulas put forward by different scholars are quite different. The most widely used correlation for flooding is the wellknown Wallis equations, which is defined as follows: (1)
JG1/2 + m ·JL 1/2 = C where 1/2
JG = JG
G
gD (
(2)
G)
L
1/2
JL = JL
L
gD (
L
G)
(3)
JG and JL are the non-dimensional flow rates, more commonly known as Wallis parameter. JG and JL denote the gas superficial velocity and liquid superficial velocity, respectively. There are still many unresolved problems concerning the mechanism triggering flooding as well as correlations appropriate for practical applications. Computational Fluid Dynamics (CFD) has become a powerful tool
Corresponding author. E-mail address: [email protected] (X. Zhang).
https://doi.org/10.1016/j.cryogenics.2018.12.002 Received 28 May 2018; Received in revised form 4 December 2018; Accepted 5 December 2018 Available online 08 December 2018 0011-2275/ © 2018 Elsevier Ltd. All rights reserved.
Cryogenics 97 (2019) 100–108
L. Yu et al.
Nomenclature
t u p g v t
κ
G ,G , Yk , Y JL , JG JL ,JG D FD A NR e vr Qi rd
kinetic viscosity of continuous fluid, Pa.s drag coefficient sip velocity/relative velocity blending function blending coefficient volume fraction limit particle viscosity number surface tension coefficient, N/m
µc CD U f a
void fraction density, kg/m3 time, s velocity vector, m/s pressure, Pa gravity, m/s2 viscous shear force, Pa turbulent shear force, Pa turbulent kinetic energy turbulent diffusion rate production of κ, effective diffusion coefficient of κ, diffusion of κ, superficial velocity of liquid or gas, m/s Wallis parameter for liquid, gas diameter of tube, m interfacial drag force local interfacial area density,1/m particle Reynolds number relative velocity, m/s volume flow rate, m3/s radius of particle, m
limit
Nµ
Subscript
k, i B D FS L G c
liquid (L) or gas (G) bubble flow droplet flow free surface flow liquid phase gas phase lower index for continuous phase
Superscript laminar flow turbulent flow
v t
for studying complex two-phase flow. It is noted that the key to flooding simulation lies in the accurate modeling of interfacial momentum exchange, which is very complicated because it involves flow transformation between the continuous and separated phases. In flooding process of water-air, flow pattern experiences a dynamic transition from stratified flow to slug flow, and even to bubble or droplet flow. As an drag force modelling technique, the Algebraic Interface Area Density (AIAD) model has been validated in water-air cases to well predict complicated transitions between different flow patterns [12–15]. Thomas Höhne and Deendarlianto [14] performed three-dimensional CFD simulations on a 1/3 scale hot leg with rectangular cross section with ANSYS CFX12.1 basing on the Euler-Euler modeling approach and AIAD model. In terms of flooding velocity, they obtained good agreement with published experimental data. Murase and Kinoshita [16] carried out two-dimensional numerical simulations on a hot leg and pressuizer surge line, and successfully figured out the non-dimensional superficial velocity curve. However, since the original AIAD model was proposed, it has been mainly used for water-air cases. Attentions posed on LN2-VN2 are rather few [17–20]. Chen et al. [21] ever investigated the flooding behaviors of LN2-VN2 by experiment method. However, numerical calculations have not been performed yet. This paper devotes to the 2-D numerical study on flooding mechanism of LN2-VN2 in an inclined tube, basing on the Eulerian multiphase flow approach and a modified AIAD model. Some parameters in AIAD are modified basing on physical property differences between LN2-VN2 and water-air. Results are then used to analyze the effects of tube diameter and length. Usually the flooding is discussed basing on so-called flooding curves on JG JL graph [22,23]. But in the Refs. [24,25], the flooding velocity from experimental measurements are present in the form of JG JL [24,25]. In order to make direct comparison with experiment, the latter is adopted in this paper. Wallis parameter is considered when analyzing influences of tube dimeter. There are two ways to determine the flooding onset: one is finding the maximum liquid superficial velocity in a given gas superficial velocity [22], the other is the opposite [25,26]. JG - JL curves in this paper provides the information on the maximum JG which is possible for a given JL .
2. Numerical model 2.1. Governing equations The RANS approach is frequently used to model two-phase flow involving complicated momentum exchange. Compared with one-fluid model, the Eulerian two-fluid model predicts better for small-scale bubbles or droplets, and it tends to simulate a sharper interface. According to the phenomenon investigated in Chen’s experiments [21], LN2 liquid is inevitably crushed into small-scaled droplets by the highspeed stream when flooding occurs. Considering the complexity of flooding flow, the Eulerian multiphase approach is adopted. In the Eulerian multiphase framework, a single pressure is shared by all phases. The momentum and continuity equations are solved for each phase, and the interphase drag force is used to calculate the momentum exchange between the phases. Energy equations are not solved basing on the assumption that the two phases are always in thermal equilibrium. The continuity and momentum equations are written as follows, respectively: k k
t
+
k k uk
t
·(
+
k k uk )
·(
=0
k k u k uk )
(4)
=
k
pk +
k kg
+
k(
v
+
t)
+ FD (5)
The subscript k denotes liquid (L) or gas (G). v and t are laminar shear force and turbulent shear force. Except for the local volume fraction k and additional source or sink terms to account for the interphase interaction effects, they are in fact “very much alike” to their single-phase counterparts [27]. Through AIAD technique used for modelling drag force, closure to the interfacial exchange terms FD is generally attained and the coupling of two phases is obtained. Contrast with κ- turbulence closure, the standard κ- closure can better model low-Reynolds number flows near wall. The feasibility of κmodel for flooding has been validated by some researchers [14]. It is based on the transport equations for the turbulent kinetic energy κ and
101
Cryogenics 97 (2019) 100–108
L. Yu et al.
t
t
(
)+
(
)+
xi
xi
( ui) =
(
xj
ui) =
k
xj
xj
xj
+G
to κ:
2.3. Modified algebraic interface area density (AIAD) model
Y
+G
Fig. 1 plots the blending functions fD / fB (Eqs. (6)–(8)) against G (gas volume fraction) under different values of αB, limit (same as αD, limit) and aB (same as aD). It is found that the parameter aB determines the slope of curves, the greater is the value of aB, the larger is the slope, which implies a quicker transition between phase morphologies. While, whether a specific flow pattern exists in a broader range of G is mainly determined by the parameter αB, limit, the greater is the value of αB, limit, the smaller is the range of G corresponding to the free surface region. It is thus clear that the effects of these parameters are critical and significant. Fig. 2 shows typical snapshots of flooding phenomenon for water-air and LN2-VN2 at almost the same Reynolds number, which is obtained by present authors [21]. Contrast to water-air, LN2-VN2 hold a relatively smaller surface stress and bigger gas-liquid density ratio. These disparities result in a very different flooding morphology and triggering mechanism for LN2-VN2 [21]. In LN2-VN2 case, there are more small waves along the interface before flooding, then the interface suddenly disappears and the gas passage becomes filled with tiny droplets, forming the mist flow. This differs from the water-air case that one of the developed interfacial waves grows up with a distinct wave front, which is blocked by the gas and finally forms the slug flow. After flooding quenches, the whole interface becomes smooth like before. Refer to Ref. [21] for more details. From comparative observation, it is obvious that LN2-VN2 interface holds a poor stability, for that the interfacial wave is very easily blown into liquid droplets by the gas. Therefore, the G value corresponding to free surface region should be reduced. In addition, morphology transition takes place more rapidly for LN2-VN2 flow, which means the slope of curve in transition zone is larger. That is to say, αB, limit and aB (also αD, limit and aD) in AIAD model should be larger for LN2-VN2 case. Parameter studies are then carried out by using try-out method. This is the first time that AIAD model has ever been modified, and the coefficients of αD,limit = αB,limit = 0.45 and aD = aB = 150 are adopted for LN2-VN2.
(6)
Y
(7)
2.2. Original algebraic interface area density (AIAD) model The drag force FD in AIAD is defined in below:
FD =
1 A 2
LG CD
|U |2
(8)
where A is the local interfacial area density, CD is the drag coefficient, U is the slip velocity and LG is the density of continues phase for dispersed flow, or average density of two phases for free surface flow. The values of CD and A are related to local flow pattern, for that AIAD algorithm allows detections of morphological form and the corresponding switching of each correlation from one object pair to another. In the AIAD model proposed by Thomas et al. [12–14], flow morphologies are classified into three kinds according to the volume fraction: droplets, bubbles and free surface flow. The blending function f which enable switching between these morphologies are defined as follows:
[
L
D,limit )
fB = 1 + eaB (
[
G
B,limit )
fFS = 1
fB
fD = 1 + e aD (
fD
]1
(9)
]1
(10) (11)
where L and G are the volume fractions for liquid and gas, respectively; αD, limit and αB, limit are volume fraction limits for droplet and bubble, and aD and aB are blending coefficients. In previous studies for water-air [12–14], the values of αD,limit = αB,limit = 0.3 and aD = aB = 70 are adopted basing on the recognition that bubbly flow rarely exceeds a gaseous volume fraction of about 0.25–0.35 when the transition to resolved structure occurs [28]. Then, the local interfacial area density A and drag coefficient CD are calculated as the sum of Aj (CD, j) with fj as the weighting factor (j = FS, B, D):
CD = fFS CD,FS + fB CD,B + fD CD,D
(12)
A = fFS AFS + fB AB + fD AD
(13)
2.4. Drag coefficient for different flow morphology In previous studies, the drag coefficients for bubble (CD,B ) and droplet flow (CD,D ) are both designated to be 0.44, independent of the local morphology and actual working fluids. It’s not consistent with the real situation and will deteriorate the accuracy of simulation. Ishii and Hibiki [29] stated in their book that the particle viscosity number and Reynolds number are decisive in drag correlation. They are defined as follows:
where D, B, FS denote droplet, bubble and free surface flow, respectively.
1.00
fB
0.50
B
-------
0.25 0.00 0.00
0.25
0.50
0.75
D,limit
1.00
0.75
0.75
fD
— — —
D
0.75
1.00
0.50
B/D
B/D
0.50
B/D,limit
0.25
0.25
0.25
0.50
0.75
1.00
G G
under different values of αB,limit (same as αD,limit) and aB (same as aD). 102
=0.45
0.00 0.00
1.00
=0.3
B/D B/D,limit
0.00
=0.3
B/D,limit
B,limit
G
Fig. 1. Curves of blending function against
fFS
the turbulent frequency , considered as the ratio of
Cryogenics 97 (2019) 100–108
L. Yu et al.
(a) water-air, = 30°, ReL=2878
(b) VN2-LN2, = 30°, ReL=2825
Fig. 2. Experimental observations of flooding phenomenon [21].
Nµ =
µc (
NR e =
1
c
2rd
)2
g
c
Table 1 Used physical properties of LN2 and VN2 (1 bar).
(14)
|vr |
µc
(15)
Most of the time, the drag coefficient can be presented by a function of Reynolds number, which is called viscous regime. Beyond this range is Newton’s regime, in which the drag coefficient is independent of Reynolds number but increases linearly with the size of particle since the distortion and irregular motions happened. Physical properties of LN2 and VN2 are listed in Table 1. According to Eq. (14), it is easy to know that LN2-VN2 flooding flow belongs to the viscous regime because of its Nµ value on the order of 10−3, so the drag coefficient could be calculated by following correlation [29]:
Density m3)
LN2 VN2
806.08 4.8655
(kg/
Dynamic viscosity µ (kg/ m·s)
Surface tension coefficient (N/m)
1.6065× 10−4 5.4741× 10−6
0.009
10
5% +2
8
JG, sim /m s-1
24 0.75 CD = (1 + 0.1NRe ) NR e
Fluids
(16)
6
-25
4
%
30 45 60
2
2.5. Simulation conditions Ousaka et al. [25,30–31] performed a lot of experiments on waterair, with different tube lengths, diameters, angels and fluid properties considered, and obtained massive data. Chen [24] took some of them as counterparts and carried out numerical calculations basing on the original AIAD model. Fig. 3 illustrates the comparison of calculated flooding velocity JG with the experimental data at different JL . It can be seen that the two are in fairly good agreement with an error of about 25%, which shows the reliability of original AIAD model for water-air flooding. It has been tried by the present authors that using the original AIAD model to simulate LN2-VN2 flooding. Computational domain incorporates the investigated area of a long-narrow tube in the middle and two ends having a gradual expansion structure to eliminate disturbance, as depicted in Fig. 4 Inclination is realized by setting gravity. Pressures are monitored near both ends. This 2D-model was also adopted by other researchers to simulate flooding flow [32,33]. However, the solver frequently crashes down during calculations and the discrepancy of flooding velocity between simulations and experiments is quite large (seen in Fig. 5). Therefore, modifications are performed on the original AIAD model as described above (Sections 2.3 and 2.4). Then, numerical calculations on LN2-VN2 are re-launched. Three sets of grids are used for grid independence validation, which consist of 38,640, 113,060 and 142,736 cells, respectively. We verified the grid independence basing on a simple numerical experiment: LN2 and VN2 enter the tube from gas-in and liquid-in boundary, respectively, flowing in the opposite direction at JL = 0.05 m/s and JG = 1.0 m/s. After the steady liquid film is formed in the tube, we plotted the value of area averaged gas pressure against x
0
0
2
4
6
JG, exp /m s-1
8
10
Fig. 3. Comparison of flooding velocities between simulations and experiments [24].
coordinate. The results are as given in Fig. 5. It’s is found that the second and third grid show very close predictions, while the result of first grid is not that good. So, the second grid which has 113,060 cells are adopted for the following simulations. The computational domain is initialized to be full of gas with the velocity of 1 m/s. LN2 comes into tube from the liquid-in boundary until a stable liquid film is formed. Then, gas velocity is increased by 0.05 m/s in a time step size of 1 s. Liquid/gas inlet are both set as velocity-in-boundary, and liquid/gas outlet are outflow-boundary. Liquid velocity ranges from 0.02 m/s to 0.07 m/s with an interval of 0.01 m/s. The film behaviors of interfacial wave is also investigated [32]. 3. Results and discussions 3.1. Model validation Two important superficial velocities, JL and JG, are defined as:
Ji =
103
4Q i (i = L, G) D2
(17)
Cryogenics 97 (2019) 100–108
L. Yu et al.
Fig. 4. Computational domain for flooding simulations.
error of 20% by modified AIAD model (red dot). Both simulated and experimental JG decrease monotonously with an increasing JL. This exactly validates the rationality of modified AIAD model for LN2-VN2 cases. Fig. 7 illustrates the numerically predicted fluctuations of flooding for LN2-VN2 case. It can be seen that the interfacial waves are concentrated in areas close to the liquid exit when flooding initiates. This is called lower flooding by Deendarlianto [31]. He experimented with water-air and proposed two mechanisms of flooding flow, i.e. the lower flooding and upper flooding. The term of “lower and upper” indicates the location where flooding initiates is at the lower (x/L ranged from 0.0 to 0.5, L is the length of test section) or the upper part (x/L ranged from 0.5 to 1.0) of the tube. He also pointed out that the lower flooding often occurred at lower liquid velocity and higher tube inclination, which is exactly the setup in this research.
Relative Pressure /Pa
10
38640 cells 113060 cells 142736 cells
0
-10
-20
-30
-0.2
0.0
0.2
0.4
0.6
x /m
0.8
3.2. Observation of flooding phenomenon and related features
Fig. 5. The area averaged gas pressure against x coordinate with different grid resolutions.
4
12
a
10
LN2-VN2, modified AIAD
JG /m s-1
-1
water-air, experiment [22] water-ari, simulation LN2-VN2, fitting line
4
LN2-VN2, modified AIAD LN2-VN2, fitting line
2
1
2 0
LN2-VN2, original AIAD
3
LN2-VN2, original AIAD
6
LN2-VN2, experiment
b
LN2-VN2, experiment
8
JG /m s
Fig. 8 depicts the comparison of flooding phenomenon between experiments and simulations. The contours are calculated by modified AIAD model, with LN2-VN2 as the working fluid at an inclination of 30°.
20%
20% 0.0
0.1
0.2
JL /m s
-1
0.3
0.4
0.00
0.02
0.04
0.06
JL /m s-1
0.08
0.10
Fig. 6. Comparisons of simulated flooding velocities with experimental data of LN2-VN2 and water-air.
It’s found that when gas velocity is relatively small, the two-phase flow holds in a stable stratified state as presented in Fig. 8(a). As gas velocity increases, the liquid film begins to fluctuate and the interfacial wave starts to grow up. The gradually growing wave blocks some of the channel, resulting in a decrease of effective cross-sectional area for gas flow. So, the gas flow begins speeding at the wave crest, which will promote the interfacial wave to flow backward (Fig. 8(b)). Finally, when a certain velocity is reached, the local accelerated stream will suddenly tear up the interface wave, causing the appearance of spray entrainment, sometimes even mist flow (Fig. 8(c)). The evolution of interface wave calculated in this paper is basically the same as that observed in experiment [21], with the only difference that the growing wave will finally contact the upper wall, forming slug flows for waterair situation [21]. A careful observation on flooding phenomenon of LN2-VN2 reveals
where Qi is the volume flow rate of liquid or gas, D is the tube diameter. The present authors plotted flooding velocity curve on JG - JL figure. Fig. 6(a) depicts experimental and simulated flooding curves, and the lower three are coming from LN2-VN2. The enlarged details about LN2-VN2 curve are demonstrated in Fig. 6(b). As can be seen from Fig. 6(a), flooding velocity of water-air is much larger than that of LN2VN2. This is primarily caused by the differences of physical properties [21]. When the original AIAD model is used to model flooding process of water, satisfactory agreement with experimental data can be obtained [24]. It indicates that the original AIAD model is well qualified for flooding calculations for water. However, when it is used for LN2VN2, the predictions are found to be excessively high. As shown in Fig. 6(b), the error between original AIAD and experiments is so high as nearly 50% (seen in Fig. 6(a)). As for the modified AIAD model, the measured flooding gas velocity JG (blue 1dot) is predicted within an
1
(footnote continued) version of this article.
For interpretation of color in Fig. 6 and 10, the reader is referred to the web 104
Cryogenics 97 (2019) 100–108
L. Yu et al.
Fig. 7. Interfacial waves at flooding.
Fig. 8. Comparisons of flooding process between experiment and simulation [21].
that its occurrence is always accompanied by several distinct features, one of which is the spray entrainment. The entrainment of liquid droplets is usually initiated at the lower part near liquid outlet, but it can also take place in upper area as flooding progresses. After some droplets are entrained for a certain distance, their moving speeds will decrease due to the attenuation of kinetic energy and those droplets may deposit onto the liquid film. It increase the perturbation of liquid film at that location, thus accelerating the transition to flooding. In simulated dynamic process, it is clearly observed that after droplets land at a certain point, an interface wave is generated there immediately and then gradually grow up during its downward movement. As demonstrated in Fig. 9, it is an example of entrainment phenomenon obtained by numerical method ( JL = 0.05m/s ). At first, an interface wave is caused by tiny disturbances (a). Due to the shear forces exerted by gas flow, some liquid at a crest will be deformed into a long ribbonlike lump (b). Then, this lump moves forward and finally lands at a certain point on liquid film (c). External intruding momentum destroys the flow stability at that location and thus the second interface wave is generated (d) there. During this process, the first interface wave appears to “stay on” the surface of liquid film. This is because the inertia of wave’s downward motion could substantially cancel out the carrying force from gas flow, i.e. this wave is on the onset of reversing its’ direction at this moment. This entrainment process coincides well with experimental observations (Fig. 2(b)). Dukler and Smith [34] interpreted this phenomenon as a change of flow rate of liquid. They observed that the droplets falling into liquid film caused a circulating flow inside the point, so the actual liquid down-flow rate at this point was higher than that injected into the tube. Therefore, a relatively lower gas velocity may trigger the occurrence of flooding. Another important characteristic of flooding phenomenon is the large-scale disturbance wave. The entrained spray is exactly originated
from here. After the large-scale disturbance waves are shredded by strong shear stress of ascending gas flow, they are turned into numerous tiny droplets, which will be carried forward by the upward gas, then forming the spray entrainment. The appearance of large-scale wave will lead to an increase in interfacial roughness, which is also an increase in shear stress. In turn, the increased shear stress will hinder the downward movement of liquid film and thus stimulate the occurrence of flooding. When gas velocity approaches the flooding onset, a steep increase in pressure drop is also an important feature of flooding phenomenon. In order to obtain the pressure-drop data, detection surfaces are respectively disposed near the liquid inlet and outlet. Fig. 10 shows the pressure-drop curves detected in this research. The inlet gas velocities are also plotted, indicated as blue lines. It can be seen that before reaching the flooding onset, the pressure drop between two ends shows a gradual increasing trend due to the increasing gas velocity. And each time the gas velocity is lifted by a small step, the pressure-drop curves show a small jump correspondingly. Whereas, when the flooding onset occurs, the pressure-drop abruptly increases in a great amplitude. For this reason, this has been taken by many researchers [4,26], also the present authors, as a tool to accurately predict the occurrence of flooding onset. What is slightly different in this research is that: after the flooding point is reached, the gas velocity will still be lifted as a function of time, so a steep increase in pressure drop may occur several times. Therefore, this paper takes the first significant rise in pressure drop as the judgment basis for flooding. The occurrence of flooding takes away a large volume of liquid. This makes liquid film become very thin after flooding and the interface will be stable again, so the pressure drop also correspondingly falls back. Bankoff [1] ever stated that the pressure drop of multiphase countercurrent flow consists of three components: the momentum reduction term, frictional resistance term 105
Cryogenics 97 (2019) 100–108
L. Yu et al.
400
1.5 1.0
200
20
1.5 1.0
400
30
Time /s
0
50 0.0
40
1600
0
10
20
30
Time /s
JL = 0.07m/s
Gas veloctiy Pressure drop
1200
0.5
small jump
3.0 2.5 2.0 -1
10
800
1.5
800
1.0 400
0.5
small jump 0
0
10
JG /m s
0
2.0
0.5
small jump
0
1200
JL = 0.05m/s Gas velocity 2.5 Pressure drop -1
2.0
JG /m s
600
Pressure drop /Pa
2.5 -1
800
3.0
1600
3.0
JG /m s
JL = 0.04m/s Gas velocity Pressure drop
Pressure drop /Pa
Pressure drop /Pa
Fig. 9. Entrainment phenomenon on the onset of flooding.
20
30
Time /s
40
50
0.0
Fig. 10. Characteristics of pressure-drop in an inclined tube (L = 1.0 m, D =20 mm).
106
40
50
0.0
Cryogenics 97 (2019) 100–108
L. Yu et al.
3.0
JG m s-1
2.5
configuration. In this paper, the effects of tube diameter and length are taken into consideration. The effects of tube diameter have been widely concerned because it is related to possible scaling laws in a single-channel flow. Five groups, whose diameters are 10, 15, 20, 25, 30 mm, respectively, are simulated in this paper, with superficial velocity of LN2 to be 0.05 m/s. The critical VN2 velocities on the onset of flooding are recorded in Fig. 11(a). The superficial gas velocities are found to be largely dependent on the diameter, which tells that a higher JG is required to trigger flooding in larger diameter tubes. Some earlier studies [25] also predicted a similar trend of JG for water-air flow when diameter rises at a fixed JL . For the pipes whose diameter is smaller than 10 mm, Mouza [35] have conducted a series of valuable and pioneering experiments and found that, the flooding mechanism is different for small diameter pipes. So the problem of whether the modified AIAD model can be extended to small diameter cases needs more research work in the future. As can be seen from the mathematical expression (Eqs. (1)–(3)), JG* and JL* are dimensionless parameters containing the diameter impact. Some analytical models and correlations [36,37] predict that the gas flooding velocity JG monotonically rises with the diameter for water-air pair, however, the trend is opposite when it comes to the Wallis parameter JG*. For LN2-VN2 cases, a similar trend is found for JG as shown in Fig. 11(a), while the discipline of Wallis parameter JG* is slightly different. According to the Fig. 11(b), JG* keeps a slow increasing trend with the diameter under a given liquid superficial velocity. It means that the gas flooding superficial velocity is not proportional to the square root of D as indicated in JG*. The Wallis number consider the effect of the tube diameter explicitly in the denominator, but the tube diameter effect seems to be more complicated than the Wallis parameter indicates [1]. The effects of tube length on flooding for water-air have also been studied by many researchers [25,26]. However, the discussions on the same issue for LN2-VN2 cases are still lacking. So here, the effect of length is discussed by investigating the trend of JG when altering tube length from 750 to 1050 mm. The tube diameter is set at 20 mm and the liquid superficial velocity (JL ) is kept as 0.05 m/s. Results are listed in Fig. 12. The critical gas velocities appear to be inversely proportional to the tube length, and the slope is about −2.3 s−1. Ousaka [25] also predicted a monotonic decreasing trend of JG against the tube length when investigating flooding phenomenon of water-air, while the slope was slightly different as −1.6 s−1.
(a)
L = 1.0 m JL=0.05 m/s
/
2.0 1.5 1.0 0.5 0.0 10
JG* /Wallis Parameter
10-3 0.5
0.4
15
20 25 D /mm
30
(b)
L=1.0 m JL=0.05 m/s
0.3 0.2 0.1 0.0
10
15
20 D /mm
25
30
Fig. 11. Effects of tube diameter on flooding velocity-(a): vs. D; (b): (Wallis Parameter) vs. D.
3.0
D= 20 mm JL=0.05m/s
2.5
JG /m s
-1
2.0 1.5
4. Conclusions
1.0
The original AIAD model, which was well qualified for flooding calculations for water, was modified on some coefficients to predict the flooding velocity of LN2-VN2. These modifications were verified by published experimental data from LN2-VN2 cases. Then, with this tool, numerical simulations on flooding mechanism of LN2-VN2 were carried out, and the following results were drawn:
0.5 0.0
750 800 850 900 950 1000 1050 L /mm
(1) The flooding curve were plotted on JG - JL figure, and good agreement with the experimental data was obtained. Besides, it showed that the interfacial fluctuations were almost concentrated in tube’s “lower part”, which was the same as Deendarlianto’s observation. These modifications were thus appropriately verified. (2) It appeared that the flooding process of LN2-VN2 went through three stages from stratified flow to wavy flow, and finally to mist flow with the appearance of spray entrainment. And further careful observations revealed that its occurrence was always accompanied by several features, which were the spray entrainment, the largescale disturbance wave and the steep increase in pressure drop. (3) Effects of the tube diameter and length were investigated in this research. And it’s found that the flooding velocity would be higher
Fig. 12. Effects of tube length on flooding velocity.
and gravity term. The first term could be omitted when flow regime in tube is relatively stable, accompanied by a negligible variation of liquid film thickness. However, after flooding occurs, the frictional resistance term becomes dominant due to obstruction of large-scale disturbance wave at the interface. 3.3. Influences of geometric factors Unlike water-air fluid pair [9,25], there is still a shortage of exploration for LN2-VN2 cases about the topic of influences of geometry 107
Cryogenics 97 (2019) 100–108
L. Yu et al.
with larger tube diameter and smaller tube length. (4) The modified AIAD model proposed in the paper is only valid for LN2-VN2 system because the coefficients is determined according to the properties of LN2-VN2 pair. And the model is validated by the experimental data of LN2-VN2 system.
2018;124:269–78. [16] Murase M, Kinoshita I, Kusunoki T, Lucas D, Tomiyama A. Countercurrent flow limitation in slightly inclined tubes with elbows. J Nucl Eng Radiation Sci 2015;1(4):041009. [17] Fiedler S, Auracher H. Experimental and theoretical investigation of reflux condensation in an inclined small diameter tube. Int J Heat Mass Transf 2004;47(19–20):4031–43. [18] Barnea D, Yoseph NB, Taitel Y. Flooding in inclined tubes — effect of entrance section. Can J Chem Eng 2010;64(2):177–84. [19] Kinoshita I, Murase M, Nariai T, Tomiyama A. Countercurrent gas-liquid flow in a rectangular channel simulating a PWR hot leg (3): evaluation of effects of fluid properties using VOF method. Japan J Multiph Flow 2010;24(4):445–52. [20] Deendarlianto, Höhne T, Apanasevich P, Lucas D, Vallée C, Beyer M. Application of a new drag coefficient model at CFD-simulations on free surface flows relevant for the nuclear reactor safety analysis. Ann Nucl Energy 2012;39(1):70–82. [21] Chen J, Xu L, Xiong W, Qiu L, Zhang X. Experimental results of flooding experiments in an inclined tube with liquid nitrogen and its vapor. Cryogenics 2014;62:1–6. [22] Wang ZW, Tian WX, Yu JT, Zhang DL, Su GH, Qiu SZ, Chen RH, Dong B. Experimental investigation of CCFL in pressurizer surge line with air-water. International Topical Meeting on Nuclear Reactor Thermal Hydraulics Nureth. 2015. [23] Wongwises S. Effect of inclination angles and upper end conditions on the countercurrent flow limitation in straight circular pipes. Int Commun Heat Mass Transfer 1998;25(1):117–25. [24] Chen J, Tang Y, Zhang W, Wang Y, Qiu L, Zhang X. Computational fluid dynamic simulations on liquid film behaviors at flooding in an inclined tube ☆. Chin J Chem Eng 2015;23(9):1460–8. [25] Ousaka A, Deendarlianto, Kariyasaki A, Fukano T. Prediction of flooding gas velocity in gas–liquid counter-current two-phase flow in inclined tubes. Nucl Eng Des 2006;236(12):1282–92. [26] Guo T, Ji HJ. Experimental study on flooding and flow reversal in small diameter tubes with various inclinations and horizontal lengths. Int J Refrig 2014;38(1):290–8. [27] Yeoh GH, Tu J. Computational Techniques for Multiphase Flows. Elsevier; 2009. [28] Höhne T, Mehlhoop JP. Validation of closure models for interfacial drag and turbulence in numerical simulations of horizontal stratified gas–liquid flows. Int J Multiph Flow 2014;62(2):1–16. [29] Ishii M, Hibiki T. Thermo-fluid Dynamics of Two-phase Flow. US: Springer; 2006. [30] Deendarlianto, Ousaka A, Indarto, Kariyasaki A, Lucas D, Vierow K, Valléé C, Hogan K. The effects of surface tension on flooding in counter-current two-phase flow in an inclined tube. Exp Therm Fluid Sci 2010;34(7):813–26. [31] Deendarlianto, Ousaka A, Kariyasaki A, Fukano T. Investigation of liquid film behavior at the onset of flooding during adiabatic counter-current air–water twophase flow in an inclined tube. Nucl Eng Des 2005;235(21):2281–94. [32] Trifonov YY. Flooding in two-phase counter-current flows: numerical investigation of the gas–liquid wavy interface using the Navier-Stokes equations. Int J Multiph Flow 2010;36(7):549–57. [33] Valluri P, Spelt PDM, Lawrence CJ, Hewitt GF. Numerical simulation of the onset of slug initiation in laminar horizontal channel flow. Int J Multiph Flow 2008;34(2):206–25. [34] Dukler AE, Smith L. Two phase interactions in counter-current flow: studies of the flooding mechanism. Annual report, November 1975–October 1977 [PWR]. 1977. [35] Mouza AA, Paras SV, Karabelas AJ. Incipient flooding in inclined tubes of small diameter. Int J Multiph Flow 2003;29(9):1395–412. [36] Chung SS, Tien LP, Liu CL. Flooding in two-phase countercurrent flows. II. Experimental investigation. Physicochem Hydrodyn 1980;1:209–20. [37] Richter HJ. Flooding in tubes and annuli. Int J Multiph Flow 1981;7:647–58.
Acknowledgments This work is supported by the National Natural Science Foundation of China (No.51576169; 51636007) and the Research Fund of State Key Laboratory of Technologies in Space Cryogenic Propellants and National key R&D plan of China (2017YFB0603702). References [1] Bankoff SG, Sang CL. A Critical Review of the Flooding Literature. Berlin Heidelberg: Springer; 1986. [2] Issa AA, Macian R. A review of CCFL phenomenon. Ann Nucl Energy 2011;38(9):1795–819. [3] Vallée C, Nariai T, Futatsugi T, Tomiyama A, Lucas D, Murase M. Experimental characterisation of the interfacial structure during counter-current flow limitation in a model of the hot leg of a PWR. Sci Technol Nuclear Instal 2012;2012:339–42. [4] Vallée C, Seidel T, Lucas D, Beyer M, Prasser HM, Pietruske H, et al. Counter-current flow limitation in a model of the hot leg of a PWR—comparison between air/water and steam/water experiments. Nucl Eng Des 2012;245(3):113–24. [5] Kusunoki T, Nozue T, Hayashi K, Hosokawa S, Tomiyama A, Murase M. Condensation experiments for counter-current flow limitation in an inverted U-tube. J Nucl Sci Technol 2016;53(4):486–95. [6] Issa AA, Macian R. Experimental investigation of countercurrent flow limitation (CCFL) in a large-diameter hot-leg geometry: a detailed description of CCFL mechanisms, flow patterns and high-quality HSC imaging of the interfacial structure in a 1/3.9 scale of PWR geometry. Nucl Eng Des 2014;280:550–63. [7] Issa AA, Macian R. Experimental investigation and CFD validation of countercurrent flow limitation (CCFL) in a large-diameter PWR hot-leg geometry. J Nucl Sci Technol 2016;53(5):647–55. [8] Zapke A, Kröger DG. The influence of fluid properties and inlet geometry on flooding in vertical and inclined tubes. Int J Multiph Flow 1996;22(3):461–72. [9] Zapke A, Kröger DG. Countercurrent gas–liquid flow in inclined and vertical ducts — I: flow patterns, pressure drop characteristics and flooding. Int J Multiph Flow 2000;26(9):1439–55. [10] Suzuki S, Ueda T. Behaviour of liquid films and flooding in counter-current twophase flow—part 1. Flow in circular tubes. Int J Multiph Flow 1977;3(6):517–32. [11] Jayanti S, Hewitt GF. Prediction of the slug-to-churn flow transition in vertical twophase flow. Int J Multiph Flow 1992;18(6):847–60. [12] Apanasevich P, Lucas D, Beyer M, Szalinski L. CFD based approach for modeling direct contact condensation heat transfer in two-phase turbulent stratified flows. Int J Therm Sci 2015;95:123–35. [13] Höhne T, Geissler T, Bieberle A, Hampel U. Numerical modeling of a horizontal annular flow experiment using a droplet entrainment model. Ann Nucl Energy 2015;77(77):351–60. [14] Höhne T, Deendarlianto A. Numerical simulations of counter current flow experiments using a morphology detection algorithm. J Comput Multiph Flows 2012;4(3). [15] Chen J, Zeng R, Zhang X, Qiu L, Xie J. Numerical modeling of flow film boiling in cryogenic chilldown process using the AIAD framework. Int J Heat Mass Transf
108
Contents lists available at ScienceDirect
Cryogenics journal homepage: www.elsevier.com/locate/cryogenics
Research paper
Simulations on LN2–VN2 flooding phenomenon in inclined tubes using a modified AIAD model
T
Liu Yua, Hong Chenb, Xu Gaob, Rui Zhoua, Huangjun Xiea, Limin Qiua, Xiaobin Zhanga,
⁎
a
Institute of Refrigeration and Cryogenics, Zhejiang University, The Key Laboratory of Refrigeration and Cryogenic Technology of Zhejiang Province, Hangzhou 310027, China b State Key Laboratory of Technologies in Space Cryogenic Propellants, Beijing 10028, China
ARTICLE INFO
ABSTRACT
Keywords: Flooding LN2 AIAD model CFD Inclined tube
Simulations basing on the Eulerian multiphase model and a modified algebraic interface area density (AIAD) model are carried out to study the flooding mechanisms with liquid nitrogen (LN2) and its vapor (VN2) in an inclined pipe. The original AIAD model predicts well for the water-air flooding velocity, which assumes three flow patterns in the whole hydraulic process, including the free surface flow, bubbly flow and droplet flow. Much different flooding phenomenon for LN2-VN2 is experimentally observed by the present authors, which indicates that some parameters have to be modified for LN2-VN2 occasion. The flooding onset is mainly identified by the sharp increase of pressure drop in this paper. As expected, the modified AIAD model predicts LN2-VN2 flooding velocity with a higher accordance with experiment than before. Both the experiments and calculations show that the flooding process of LN2-VN2 goes through three stages from stratified flow to wavy flow, and finally to mist flow with the appearance of spray entrainment, which is different from the final slug flow in water-air cases. Influences of geometric factors, such as tube length and diameter, are also investigated in this research.
1. Introduction Gas-liquid countercurrent flow is always hindered by the interfacial friction, and the greater is the relative velocity, the greater is the frictional force. When the geometry and working fluid pair are determined, there is a critical relative velocity, beyond which the flow gets into an unstable state. This critical velocity is the point at which flooding occurs [1]. Countercurrent two-phase flow is frequently encountered in a lot of industrial situations. It exists widely in nuclear reactors under abnormal operation, gravity heat pipes and distillation towers. Flooding encountered in industrial occasion often deteriorates the regular operation of devices [2]. For example, Flooding in gravity heat pipe tends to dry out the pipe, which means invalidation of equipment. The gas/ liquid system complete the heat and mass transfer through countercurrent contact in the structured packed distillation tower. If the gas velocity is too large, the liquid cannot continue to flow downward due to flooding. This will destroy the mass transfer and reduce the efficiency. So, the flooding velocity between phases is a vital parameter that we need to determine to maintain a safe operating environment. Large amount of experiments have been carried out in scaled models by Valléé [3,4], Kusunoki [5] and Issa [6,7], to investigate the flooding behaviors in the hot leg of a pressurized water reactor. Besides, Zapke
⁎
and Kröger [8,9], Suzuki [10], Jayanti and Hewitt [11] also studied mechanisms of flooding in vertical, inclined and horizontal tubes. Empirical correlations were concluded basing on these experimental data. However, formulas put forward by different scholars are quite different. The most widely used correlation for flooding is the wellknown Wallis equations, which is defined as follows: (1)
JG1/2 + m ·JL 1/2 = C where 1/2
JG = JG
G
gD (
(2)
G)
L
1/2
JL = JL
L
gD (
L
G)
(3)
JG and JL are the non-dimensional flow rates, more commonly known as Wallis parameter. JG and JL denote the gas superficial velocity and liquid superficial velocity, respectively. There are still many unresolved problems concerning the mechanism triggering flooding as well as correlations appropriate for practical applications. Computational Fluid Dynamics (CFD) has become a powerful tool
Corresponding author. E-mail address: [email protected] (X. Zhang).
https://doi.org/10.1016/j.cryogenics.2018.12.002 Received 28 May 2018; Received in revised form 4 December 2018; Accepted 5 December 2018 Available online 08 December 2018 0011-2275/ © 2018 Elsevier Ltd. All rights reserved.
Cryogenics 97 (2019) 100–108
L. Yu et al.
Nomenclature
t u p g v t
κ
G ,G , Yk , Y JL , JG JL ,JG D FD A NR e vr Qi rd
kinetic viscosity of continuous fluid, Pa.s drag coefficient sip velocity/relative velocity blending function blending coefficient volume fraction limit particle viscosity number surface tension coefficient, N/m
µc CD U f a
void fraction density, kg/m3 time, s velocity vector, m/s pressure, Pa gravity, m/s2 viscous shear force, Pa turbulent shear force, Pa turbulent kinetic energy turbulent diffusion rate production of κ, effective diffusion coefficient of κ, diffusion of κ, superficial velocity of liquid or gas, m/s Wallis parameter for liquid, gas diameter of tube, m interfacial drag force local interfacial area density,1/m particle Reynolds number relative velocity, m/s volume flow rate, m3/s radius of particle, m
limit
Nµ
Subscript
k, i B D FS L G c
liquid (L) or gas (G) bubble flow droplet flow free surface flow liquid phase gas phase lower index for continuous phase
Superscript laminar flow turbulent flow
v t
for studying complex two-phase flow. It is noted that the key to flooding simulation lies in the accurate modeling of interfacial momentum exchange, which is very complicated because it involves flow transformation between the continuous and separated phases. In flooding process of water-air, flow pattern experiences a dynamic transition from stratified flow to slug flow, and even to bubble or droplet flow. As an drag force modelling technique, the Algebraic Interface Area Density (AIAD) model has been validated in water-air cases to well predict complicated transitions between different flow patterns [12–15]. Thomas Höhne and Deendarlianto [14] performed three-dimensional CFD simulations on a 1/3 scale hot leg with rectangular cross section with ANSYS CFX12.1 basing on the Euler-Euler modeling approach and AIAD model. In terms of flooding velocity, they obtained good agreement with published experimental data. Murase and Kinoshita [16] carried out two-dimensional numerical simulations on a hot leg and pressuizer surge line, and successfully figured out the non-dimensional superficial velocity curve. However, since the original AIAD model was proposed, it has been mainly used for water-air cases. Attentions posed on LN2-VN2 are rather few [17–20]. Chen et al. [21] ever investigated the flooding behaviors of LN2-VN2 by experiment method. However, numerical calculations have not been performed yet. This paper devotes to the 2-D numerical study on flooding mechanism of LN2-VN2 in an inclined tube, basing on the Eulerian multiphase flow approach and a modified AIAD model. Some parameters in AIAD are modified basing on physical property differences between LN2-VN2 and water-air. Results are then used to analyze the effects of tube diameter and length. Usually the flooding is discussed basing on so-called flooding curves on JG JL graph [22,23]. But in the Refs. [24,25], the flooding velocity from experimental measurements are present in the form of JG JL [24,25]. In order to make direct comparison with experiment, the latter is adopted in this paper. Wallis parameter is considered when analyzing influences of tube dimeter. There are two ways to determine the flooding onset: one is finding the maximum liquid superficial velocity in a given gas superficial velocity [22], the other is the opposite [25,26]. JG - JL curves in this paper provides the information on the maximum JG which is possible for a given JL .
2. Numerical model 2.1. Governing equations The RANS approach is frequently used to model two-phase flow involving complicated momentum exchange. Compared with one-fluid model, the Eulerian two-fluid model predicts better for small-scale bubbles or droplets, and it tends to simulate a sharper interface. According to the phenomenon investigated in Chen’s experiments [21], LN2 liquid is inevitably crushed into small-scaled droplets by the highspeed stream when flooding occurs. Considering the complexity of flooding flow, the Eulerian multiphase approach is adopted. In the Eulerian multiphase framework, a single pressure is shared by all phases. The momentum and continuity equations are solved for each phase, and the interphase drag force is used to calculate the momentum exchange between the phases. Energy equations are not solved basing on the assumption that the two phases are always in thermal equilibrium. The continuity and momentum equations are written as follows, respectively: k k
t
+
k k uk
t
·(
+
k k uk )
·(
=0
k k u k uk )
(4)
=
k
pk +
k kg
+
k(
v
+
t)
+ FD (5)
The subscript k denotes liquid (L) or gas (G). v and t are laminar shear force and turbulent shear force. Except for the local volume fraction k and additional source or sink terms to account for the interphase interaction effects, they are in fact “very much alike” to their single-phase counterparts [27]. Through AIAD technique used for modelling drag force, closure to the interfacial exchange terms FD is generally attained and the coupling of two phases is obtained. Contrast with κ- turbulence closure, the standard κ- closure can better model low-Reynolds number flows near wall. The feasibility of κmodel for flooding has been validated by some researchers [14]. It is based on the transport equations for the turbulent kinetic energy κ and
101
Cryogenics 97 (2019) 100–108
L. Yu et al.
t
t
(
)+
(
)+
xi
xi
( ui) =
(
xj
ui) =
k
xj
xj
xj
+G
to κ:
2.3. Modified algebraic interface area density (AIAD) model
Y
+G
Fig. 1 plots the blending functions fD / fB (Eqs. (6)–(8)) against G (gas volume fraction) under different values of αB, limit (same as αD, limit) and aB (same as aD). It is found that the parameter aB determines the slope of curves, the greater is the value of aB, the larger is the slope, which implies a quicker transition between phase morphologies. While, whether a specific flow pattern exists in a broader range of G is mainly determined by the parameter αB, limit, the greater is the value of αB, limit, the smaller is the range of G corresponding to the free surface region. It is thus clear that the effects of these parameters are critical and significant. Fig. 2 shows typical snapshots of flooding phenomenon for water-air and LN2-VN2 at almost the same Reynolds number, which is obtained by present authors [21]. Contrast to water-air, LN2-VN2 hold a relatively smaller surface stress and bigger gas-liquid density ratio. These disparities result in a very different flooding morphology and triggering mechanism for LN2-VN2 [21]. In LN2-VN2 case, there are more small waves along the interface before flooding, then the interface suddenly disappears and the gas passage becomes filled with tiny droplets, forming the mist flow. This differs from the water-air case that one of the developed interfacial waves grows up with a distinct wave front, which is blocked by the gas and finally forms the slug flow. After flooding quenches, the whole interface becomes smooth like before. Refer to Ref. [21] for more details. From comparative observation, it is obvious that LN2-VN2 interface holds a poor stability, for that the interfacial wave is very easily blown into liquid droplets by the gas. Therefore, the G value corresponding to free surface region should be reduced. In addition, morphology transition takes place more rapidly for LN2-VN2 flow, which means the slope of curve in transition zone is larger. That is to say, αB, limit and aB (also αD, limit and aD) in AIAD model should be larger for LN2-VN2 case. Parameter studies are then carried out by using try-out method. This is the first time that AIAD model has ever been modified, and the coefficients of αD,limit = αB,limit = 0.45 and aD = aB = 150 are adopted for LN2-VN2.
(6)
Y
(7)
2.2. Original algebraic interface area density (AIAD) model The drag force FD in AIAD is defined in below:
FD =
1 A 2
LG CD
|U |2
(8)
where A is the local interfacial area density, CD is the drag coefficient, U is the slip velocity and LG is the density of continues phase for dispersed flow, or average density of two phases for free surface flow. The values of CD and A are related to local flow pattern, for that AIAD algorithm allows detections of morphological form and the corresponding switching of each correlation from one object pair to another. In the AIAD model proposed by Thomas et al. [12–14], flow morphologies are classified into three kinds according to the volume fraction: droplets, bubbles and free surface flow. The blending function f which enable switching between these morphologies are defined as follows:
[
L
D,limit )
fB = 1 + eaB (
[
G
B,limit )
fFS = 1
fB
fD = 1 + e aD (
fD
]1
(9)
]1
(10) (11)
where L and G are the volume fractions for liquid and gas, respectively; αD, limit and αB, limit are volume fraction limits for droplet and bubble, and aD and aB are blending coefficients. In previous studies for water-air [12–14], the values of αD,limit = αB,limit = 0.3 and aD = aB = 70 are adopted basing on the recognition that bubbly flow rarely exceeds a gaseous volume fraction of about 0.25–0.35 when the transition to resolved structure occurs [28]. Then, the local interfacial area density A and drag coefficient CD are calculated as the sum of Aj (CD, j) with fj as the weighting factor (j = FS, B, D):
CD = fFS CD,FS + fB CD,B + fD CD,D
(12)
A = fFS AFS + fB AB + fD AD
(13)
2.4. Drag coefficient for different flow morphology In previous studies, the drag coefficients for bubble (CD,B ) and droplet flow (CD,D ) are both designated to be 0.44, independent of the local morphology and actual working fluids. It’s not consistent with the real situation and will deteriorate the accuracy of simulation. Ishii and Hibiki [29] stated in their book that the particle viscosity number and Reynolds number are decisive in drag correlation. They are defined as follows:
where D, B, FS denote droplet, bubble and free surface flow, respectively.
1.00
fB
0.50
B
-------
0.25 0.00 0.00
0.25
0.50
0.75
D,limit
1.00
0.75
0.75
fD
— — —
D
0.75
1.00
0.50
B/D
B/D
0.50
B/D,limit
0.25
0.25
0.25
0.50
0.75
1.00
G G
under different values of αB,limit (same as αD,limit) and aB (same as aD). 102
=0.45
0.00 0.00
1.00
=0.3
B/D B/D,limit
0.00
=0.3
B/D,limit
B,limit
G
Fig. 1. Curves of blending function against
fFS
the turbulent frequency , considered as the ratio of
Cryogenics 97 (2019) 100–108
L. Yu et al.
(a) water-air, = 30°, ReL=2878
(b) VN2-LN2, = 30°, ReL=2825
Fig. 2. Experimental observations of flooding phenomenon [21].
Nµ =
µc (
NR e =
1
c
2rd
)2
g
c
Table 1 Used physical properties of LN2 and VN2 (1 bar).
(14)
|vr |
µc
(15)
Most of the time, the drag coefficient can be presented by a function of Reynolds number, which is called viscous regime. Beyond this range is Newton’s regime, in which the drag coefficient is independent of Reynolds number but increases linearly with the size of particle since the distortion and irregular motions happened. Physical properties of LN2 and VN2 are listed in Table 1. According to Eq. (14), it is easy to know that LN2-VN2 flooding flow belongs to the viscous regime because of its Nµ value on the order of 10−3, so the drag coefficient could be calculated by following correlation [29]:
Density m3)
LN2 VN2
806.08 4.8655
(kg/
Dynamic viscosity µ (kg/ m·s)
Surface tension coefficient (N/m)
1.6065× 10−4 5.4741× 10−6
0.009
10
5% +2
8
JG, sim /m s-1
24 0.75 CD = (1 + 0.1NRe ) NR e
Fluids
(16)
6
-25
4
%
30 45 60
2
2.5. Simulation conditions Ousaka et al. [25,30–31] performed a lot of experiments on waterair, with different tube lengths, diameters, angels and fluid properties considered, and obtained massive data. Chen [24] took some of them as counterparts and carried out numerical calculations basing on the original AIAD model. Fig. 3 illustrates the comparison of calculated flooding velocity JG with the experimental data at different JL . It can be seen that the two are in fairly good agreement with an error of about 25%, which shows the reliability of original AIAD model for water-air flooding. It has been tried by the present authors that using the original AIAD model to simulate LN2-VN2 flooding. Computational domain incorporates the investigated area of a long-narrow tube in the middle and two ends having a gradual expansion structure to eliminate disturbance, as depicted in Fig. 4 Inclination is realized by setting gravity. Pressures are monitored near both ends. This 2D-model was also adopted by other researchers to simulate flooding flow [32,33]. However, the solver frequently crashes down during calculations and the discrepancy of flooding velocity between simulations and experiments is quite large (seen in Fig. 5). Therefore, modifications are performed on the original AIAD model as described above (Sections 2.3 and 2.4). Then, numerical calculations on LN2-VN2 are re-launched. Three sets of grids are used for grid independence validation, which consist of 38,640, 113,060 and 142,736 cells, respectively. We verified the grid independence basing on a simple numerical experiment: LN2 and VN2 enter the tube from gas-in and liquid-in boundary, respectively, flowing in the opposite direction at JL = 0.05 m/s and JG = 1.0 m/s. After the steady liquid film is formed in the tube, we plotted the value of area averaged gas pressure against x
0
0
2
4
6
JG, exp /m s-1
8
10
Fig. 3. Comparison of flooding velocities between simulations and experiments [24].
coordinate. The results are as given in Fig. 5. It’s is found that the second and third grid show very close predictions, while the result of first grid is not that good. So, the second grid which has 113,060 cells are adopted for the following simulations. The computational domain is initialized to be full of gas with the velocity of 1 m/s. LN2 comes into tube from the liquid-in boundary until a stable liquid film is formed. Then, gas velocity is increased by 0.05 m/s in a time step size of 1 s. Liquid/gas inlet are both set as velocity-in-boundary, and liquid/gas outlet are outflow-boundary. Liquid velocity ranges from 0.02 m/s to 0.07 m/s with an interval of 0.01 m/s. The film behaviors of interfacial wave is also investigated [32]. 3. Results and discussions 3.1. Model validation Two important superficial velocities, JL and JG, are defined as:
Ji =
103
4Q i (i = L, G) D2
(17)
Cryogenics 97 (2019) 100–108
L. Yu et al.
Fig. 4. Computational domain for flooding simulations.
error of 20% by modified AIAD model (red dot). Both simulated and experimental JG decrease monotonously with an increasing JL. This exactly validates the rationality of modified AIAD model for LN2-VN2 cases. Fig. 7 illustrates the numerically predicted fluctuations of flooding for LN2-VN2 case. It can be seen that the interfacial waves are concentrated in areas close to the liquid exit when flooding initiates. This is called lower flooding by Deendarlianto [31]. He experimented with water-air and proposed two mechanisms of flooding flow, i.e. the lower flooding and upper flooding. The term of “lower and upper” indicates the location where flooding initiates is at the lower (x/L ranged from 0.0 to 0.5, L is the length of test section) or the upper part (x/L ranged from 0.5 to 1.0) of the tube. He also pointed out that the lower flooding often occurred at lower liquid velocity and higher tube inclination, which is exactly the setup in this research.
Relative Pressure /Pa
10
38640 cells 113060 cells 142736 cells
0
-10
-20
-30
-0.2
0.0
0.2
0.4
0.6
x /m
0.8
3.2. Observation of flooding phenomenon and related features
Fig. 5. The area averaged gas pressure against x coordinate with different grid resolutions.
4
12
a
10
LN2-VN2, modified AIAD
JG /m s-1
-1
water-air, experiment [22] water-ari, simulation LN2-VN2, fitting line
4
LN2-VN2, modified AIAD LN2-VN2, fitting line
2
1
2 0
LN2-VN2, original AIAD
3
LN2-VN2, original AIAD
6
LN2-VN2, experiment
b
LN2-VN2, experiment
8
JG /m s
Fig. 8 depicts the comparison of flooding phenomenon between experiments and simulations. The contours are calculated by modified AIAD model, with LN2-VN2 as the working fluid at an inclination of 30°.
20%
20% 0.0
0.1
0.2
JL /m s
-1
0.3
0.4
0.00
0.02
0.04
0.06
JL /m s-1
0.08
0.10
Fig. 6. Comparisons of simulated flooding velocities with experimental data of LN2-VN2 and water-air.
It’s found that when gas velocity is relatively small, the two-phase flow holds in a stable stratified state as presented in Fig. 8(a). As gas velocity increases, the liquid film begins to fluctuate and the interfacial wave starts to grow up. The gradually growing wave blocks some of the channel, resulting in a decrease of effective cross-sectional area for gas flow. So, the gas flow begins speeding at the wave crest, which will promote the interfacial wave to flow backward (Fig. 8(b)). Finally, when a certain velocity is reached, the local accelerated stream will suddenly tear up the interface wave, causing the appearance of spray entrainment, sometimes even mist flow (Fig. 8(c)). The evolution of interface wave calculated in this paper is basically the same as that observed in experiment [21], with the only difference that the growing wave will finally contact the upper wall, forming slug flows for waterair situation [21]. A careful observation on flooding phenomenon of LN2-VN2 reveals
where Qi is the volume flow rate of liquid or gas, D is the tube diameter. The present authors plotted flooding velocity curve on JG - JL figure. Fig. 6(a) depicts experimental and simulated flooding curves, and the lower three are coming from LN2-VN2. The enlarged details about LN2-VN2 curve are demonstrated in Fig. 6(b). As can be seen from Fig. 6(a), flooding velocity of water-air is much larger than that of LN2VN2. This is primarily caused by the differences of physical properties [21]. When the original AIAD model is used to model flooding process of water, satisfactory agreement with experimental data can be obtained [24]. It indicates that the original AIAD model is well qualified for flooding calculations for water. However, when it is used for LN2VN2, the predictions are found to be excessively high. As shown in Fig. 6(b), the error between original AIAD and experiments is so high as nearly 50% (seen in Fig. 6(a)). As for the modified AIAD model, the measured flooding gas velocity JG (blue 1dot) is predicted within an
1
(footnote continued) version of this article.
For interpretation of color in Fig. 6 and 10, the reader is referred to the web 104
Cryogenics 97 (2019) 100–108
L. Yu et al.
Fig. 7. Interfacial waves at flooding.
Fig. 8. Comparisons of flooding process between experiment and simulation [21].
that its occurrence is always accompanied by several distinct features, one of which is the spray entrainment. The entrainment of liquid droplets is usually initiated at the lower part near liquid outlet, but it can also take place in upper area as flooding progresses. After some droplets are entrained for a certain distance, their moving speeds will decrease due to the attenuation of kinetic energy and those droplets may deposit onto the liquid film. It increase the perturbation of liquid film at that location, thus accelerating the transition to flooding. In simulated dynamic process, it is clearly observed that after droplets land at a certain point, an interface wave is generated there immediately and then gradually grow up during its downward movement. As demonstrated in Fig. 9, it is an example of entrainment phenomenon obtained by numerical method ( JL = 0.05m/s ). At first, an interface wave is caused by tiny disturbances (a). Due to the shear forces exerted by gas flow, some liquid at a crest will be deformed into a long ribbonlike lump (b). Then, this lump moves forward and finally lands at a certain point on liquid film (c). External intruding momentum destroys the flow stability at that location and thus the second interface wave is generated (d) there. During this process, the first interface wave appears to “stay on” the surface of liquid film. This is because the inertia of wave’s downward motion could substantially cancel out the carrying force from gas flow, i.e. this wave is on the onset of reversing its’ direction at this moment. This entrainment process coincides well with experimental observations (Fig. 2(b)). Dukler and Smith [34] interpreted this phenomenon as a change of flow rate of liquid. They observed that the droplets falling into liquid film caused a circulating flow inside the point, so the actual liquid down-flow rate at this point was higher than that injected into the tube. Therefore, a relatively lower gas velocity may trigger the occurrence of flooding. Another important characteristic of flooding phenomenon is the large-scale disturbance wave. The entrained spray is exactly originated
from here. After the large-scale disturbance waves are shredded by strong shear stress of ascending gas flow, they are turned into numerous tiny droplets, which will be carried forward by the upward gas, then forming the spray entrainment. The appearance of large-scale wave will lead to an increase in interfacial roughness, which is also an increase in shear stress. In turn, the increased shear stress will hinder the downward movement of liquid film and thus stimulate the occurrence of flooding. When gas velocity approaches the flooding onset, a steep increase in pressure drop is also an important feature of flooding phenomenon. In order to obtain the pressure-drop data, detection surfaces are respectively disposed near the liquid inlet and outlet. Fig. 10 shows the pressure-drop curves detected in this research. The inlet gas velocities are also plotted, indicated as blue lines. It can be seen that before reaching the flooding onset, the pressure drop between two ends shows a gradual increasing trend due to the increasing gas velocity. And each time the gas velocity is lifted by a small step, the pressure-drop curves show a small jump correspondingly. Whereas, when the flooding onset occurs, the pressure-drop abruptly increases in a great amplitude. For this reason, this has been taken by many researchers [4,26], also the present authors, as a tool to accurately predict the occurrence of flooding onset. What is slightly different in this research is that: after the flooding point is reached, the gas velocity will still be lifted as a function of time, so a steep increase in pressure drop may occur several times. Therefore, this paper takes the first significant rise in pressure drop as the judgment basis for flooding. The occurrence of flooding takes away a large volume of liquid. This makes liquid film become very thin after flooding and the interface will be stable again, so the pressure drop also correspondingly falls back. Bankoff [1] ever stated that the pressure drop of multiphase countercurrent flow consists of three components: the momentum reduction term, frictional resistance term 105
Cryogenics 97 (2019) 100–108
L. Yu et al.
400
1.5 1.0
200
20
1.5 1.0
400
30
Time /s
0
50 0.0
40
1600
0
10
20
30
Time /s
JL = 0.07m/s
Gas veloctiy Pressure drop
1200
0.5
small jump
3.0 2.5 2.0 -1
10
800
1.5
800
1.0 400
0.5
small jump 0
0
10
JG /m s
0
2.0
0.5
small jump
0
1200
JL = 0.05m/s Gas velocity 2.5 Pressure drop -1
2.0
JG /m s
600
Pressure drop /Pa
2.5 -1
800
3.0
1600
3.0
JG /m s
JL = 0.04m/s Gas velocity Pressure drop
Pressure drop /Pa
Pressure drop /Pa
Fig. 9. Entrainment phenomenon on the onset of flooding.
20
30
Time /s
40
50
0.0
Fig. 10. Characteristics of pressure-drop in an inclined tube (L = 1.0 m, D =20 mm).
106
40
50
0.0
Cryogenics 97 (2019) 100–108
L. Yu et al.
3.0
JG m s-1
2.5
configuration. In this paper, the effects of tube diameter and length are taken into consideration. The effects of tube diameter have been widely concerned because it is related to possible scaling laws in a single-channel flow. Five groups, whose diameters are 10, 15, 20, 25, 30 mm, respectively, are simulated in this paper, with superficial velocity of LN2 to be 0.05 m/s. The critical VN2 velocities on the onset of flooding are recorded in Fig. 11(a). The superficial gas velocities are found to be largely dependent on the diameter, which tells that a higher JG is required to trigger flooding in larger diameter tubes. Some earlier studies [25] also predicted a similar trend of JG for water-air flow when diameter rises at a fixed JL . For the pipes whose diameter is smaller than 10 mm, Mouza [35] have conducted a series of valuable and pioneering experiments and found that, the flooding mechanism is different for small diameter pipes. So the problem of whether the modified AIAD model can be extended to small diameter cases needs more research work in the future. As can be seen from the mathematical expression (Eqs. (1)–(3)), JG* and JL* are dimensionless parameters containing the diameter impact. Some analytical models and correlations [36,37] predict that the gas flooding velocity JG monotonically rises with the diameter for water-air pair, however, the trend is opposite when it comes to the Wallis parameter JG*. For LN2-VN2 cases, a similar trend is found for JG as shown in Fig. 11(a), while the discipline of Wallis parameter JG* is slightly different. According to the Fig. 11(b), JG* keeps a slow increasing trend with the diameter under a given liquid superficial velocity. It means that the gas flooding superficial velocity is not proportional to the square root of D as indicated in JG*. The Wallis number consider the effect of the tube diameter explicitly in the denominator, but the tube diameter effect seems to be more complicated than the Wallis parameter indicates [1]. The effects of tube length on flooding for water-air have also been studied by many researchers [25,26]. However, the discussions on the same issue for LN2-VN2 cases are still lacking. So here, the effect of length is discussed by investigating the trend of JG when altering tube length from 750 to 1050 mm. The tube diameter is set at 20 mm and the liquid superficial velocity (JL ) is kept as 0.05 m/s. Results are listed in Fig. 12. The critical gas velocities appear to be inversely proportional to the tube length, and the slope is about −2.3 s−1. Ousaka [25] also predicted a monotonic decreasing trend of JG against the tube length when investigating flooding phenomenon of water-air, while the slope was slightly different as −1.6 s−1.
(a)
L = 1.0 m JL=0.05 m/s
/
2.0 1.5 1.0 0.5 0.0 10
JG* /Wallis Parameter
10-3 0.5
0.4
15
20 25 D /mm
30
(b)
L=1.0 m JL=0.05 m/s
0.3 0.2 0.1 0.0
10
15
20 D /mm
25
30
Fig. 11. Effects of tube diameter on flooding velocity-(a): vs. D; (b): (Wallis Parameter) vs. D.
3.0
D= 20 mm JL=0.05m/s
2.5
JG /m s
-1
2.0 1.5
4. Conclusions
1.0
The original AIAD model, which was well qualified for flooding calculations for water, was modified on some coefficients to predict the flooding velocity of LN2-VN2. These modifications were verified by published experimental data from LN2-VN2 cases. Then, with this tool, numerical simulations on flooding mechanism of LN2-VN2 were carried out, and the following results were drawn:
0.5 0.0
750 800 850 900 950 1000 1050 L /mm
(1) The flooding curve were plotted on JG - JL figure, and good agreement with the experimental data was obtained. Besides, it showed that the interfacial fluctuations were almost concentrated in tube’s “lower part”, which was the same as Deendarlianto’s observation. These modifications were thus appropriately verified. (2) It appeared that the flooding process of LN2-VN2 went through three stages from stratified flow to wavy flow, and finally to mist flow with the appearance of spray entrainment. And further careful observations revealed that its occurrence was always accompanied by several features, which were the spray entrainment, the largescale disturbance wave and the steep increase in pressure drop. (3) Effects of the tube diameter and length were investigated in this research. And it’s found that the flooding velocity would be higher
Fig. 12. Effects of tube length on flooding velocity.
and gravity term. The first term could be omitted when flow regime in tube is relatively stable, accompanied by a negligible variation of liquid film thickness. However, after flooding occurs, the frictional resistance term becomes dominant due to obstruction of large-scale disturbance wave at the interface. 3.3. Influences of geometric factors Unlike water-air fluid pair [9,25], there is still a shortage of exploration for LN2-VN2 cases about the topic of influences of geometry 107
Cryogenics 97 (2019) 100–108
L. Yu et al.
with larger tube diameter and smaller tube length. (4) The modified AIAD model proposed in the paper is only valid for LN2-VN2 system because the coefficients is determined according to the properties of LN2-VN2 pair. And the model is validated by the experimental data of LN2-VN2 system.
2018;124:269–78. [16] Murase M, Kinoshita I, Kusunoki T, Lucas D, Tomiyama A. Countercurrent flow limitation in slightly inclined tubes with elbows. J Nucl Eng Radiation Sci 2015;1(4):041009. [17] Fiedler S, Auracher H. Experimental and theoretical investigation of reflux condensation in an inclined small diameter tube. Int J Heat Mass Transf 2004;47(19–20):4031–43. [18] Barnea D, Yoseph NB, Taitel Y. Flooding in inclined tubes — effect of entrance section. Can J Chem Eng 2010;64(2):177–84. [19] Kinoshita I, Murase M, Nariai T, Tomiyama A. Countercurrent gas-liquid flow in a rectangular channel simulating a PWR hot leg (3): evaluation of effects of fluid properties using VOF method. Japan J Multiph Flow 2010;24(4):445–52. [20] Deendarlianto, Höhne T, Apanasevich P, Lucas D, Vallée C, Beyer M. Application of a new drag coefficient model at CFD-simulations on free surface flows relevant for the nuclear reactor safety analysis. Ann Nucl Energy 2012;39(1):70–82. [21] Chen J, Xu L, Xiong W, Qiu L, Zhang X. Experimental results of flooding experiments in an inclined tube with liquid nitrogen and its vapor. Cryogenics 2014;62:1–6. [22] Wang ZW, Tian WX, Yu JT, Zhang DL, Su GH, Qiu SZ, Chen RH, Dong B. Experimental investigation of CCFL in pressurizer surge line with air-water. International Topical Meeting on Nuclear Reactor Thermal Hydraulics Nureth. 2015. [23] Wongwises S. Effect of inclination angles and upper end conditions on the countercurrent flow limitation in straight circular pipes. Int Commun Heat Mass Transfer 1998;25(1):117–25. [24] Chen J, Tang Y, Zhang W, Wang Y, Qiu L, Zhang X. Computational fluid dynamic simulations on liquid film behaviors at flooding in an inclined tube ☆. Chin J Chem Eng 2015;23(9):1460–8. [25] Ousaka A, Deendarlianto, Kariyasaki A, Fukano T. Prediction of flooding gas velocity in gas–liquid counter-current two-phase flow in inclined tubes. Nucl Eng Des 2006;236(12):1282–92. [26] Guo T, Ji HJ. Experimental study on flooding and flow reversal in small diameter tubes with various inclinations and horizontal lengths. Int J Refrig 2014;38(1):290–8. [27] Yeoh GH, Tu J. Computational Techniques for Multiphase Flows. Elsevier; 2009. [28] Höhne T, Mehlhoop JP. Validation of closure models for interfacial drag and turbulence in numerical simulations of horizontal stratified gas–liquid flows. Int J Multiph Flow 2014;62(2):1–16. [29] Ishii M, Hibiki T. Thermo-fluid Dynamics of Two-phase Flow. US: Springer; 2006. [30] Deendarlianto, Ousaka A, Indarto, Kariyasaki A, Lucas D, Vierow K, Valléé C, Hogan K. The effects of surface tension on flooding in counter-current two-phase flow in an inclined tube. Exp Therm Fluid Sci 2010;34(7):813–26. [31] Deendarlianto, Ousaka A, Kariyasaki A, Fukano T. Investigation of liquid film behavior at the onset of flooding during adiabatic counter-current air–water twophase flow in an inclined tube. Nucl Eng Des 2005;235(21):2281–94. [32] Trifonov YY. Flooding in two-phase counter-current flows: numerical investigation of the gas–liquid wavy interface using the Navier-Stokes equations. Int J Multiph Flow 2010;36(7):549–57. [33] Valluri P, Spelt PDM, Lawrence CJ, Hewitt GF. Numerical simulation of the onset of slug initiation in laminar horizontal channel flow. Int J Multiph Flow 2008;34(2):206–25. [34] Dukler AE, Smith L. Two phase interactions in counter-current flow: studies of the flooding mechanism. Annual report, November 1975–October 1977 [PWR]. 1977. [35] Mouza AA, Paras SV, Karabelas AJ. Incipient flooding in inclined tubes of small diameter. Int J Multiph Flow 2003;29(9):1395–412. [36] Chung SS, Tien LP, Liu CL. Flooding in two-phase countercurrent flows. II. Experimental investigation. Physicochem Hydrodyn 1980;1:209–20. [37] Richter HJ. Flooding in tubes and annuli. Int J Multiph Flow 1981;7:647–58.
Acknowledgments This work is supported by the National Natural Science Foundation of China (No.51576169; 51636007) and the Research Fund of State Key Laboratory of Technologies in Space Cryogenic Propellants and National key R&D plan of China (2017YFB0603702). References [1] Bankoff SG, Sang CL. A Critical Review of the Flooding Literature. Berlin Heidelberg: Springer; 1986. [2] Issa AA, Macian R. A review of CCFL phenomenon. Ann Nucl Energy 2011;38(9):1795–819. [3] Vallée C, Nariai T, Futatsugi T, Tomiyama A, Lucas D, Murase M. Experimental characterisation of the interfacial structure during counter-current flow limitation in a model of the hot leg of a PWR. Sci Technol Nuclear Instal 2012;2012:339–42. [4] Vallée C, Seidel T, Lucas D, Beyer M, Prasser HM, Pietruske H, et al. Counter-current flow limitation in a model of the hot leg of a PWR—comparison between air/water and steam/water experiments. Nucl Eng Des 2012;245(3):113–24. [5] Kusunoki T, Nozue T, Hayashi K, Hosokawa S, Tomiyama A, Murase M. Condensation experiments for counter-current flow limitation in an inverted U-tube. J Nucl Sci Technol 2016;53(4):486–95. [6] Issa AA, Macian R. Experimental investigation of countercurrent flow limitation (CCFL) in a large-diameter hot-leg geometry: a detailed description of CCFL mechanisms, flow patterns and high-quality HSC imaging of the interfacial structure in a 1/3.9 scale of PWR geometry. Nucl Eng Des 2014;280:550–63. [7] Issa AA, Macian R. Experimental investigation and CFD validation of countercurrent flow limitation (CCFL) in a large-diameter PWR hot-leg geometry. J Nucl Sci Technol 2016;53(5):647–55. [8] Zapke A, Kröger DG. The influence of fluid properties and inlet geometry on flooding in vertical and inclined tubes. Int J Multiph Flow 1996;22(3):461–72. [9] Zapke A, Kröger DG. Countercurrent gas–liquid flow in inclined and vertical ducts — I: flow patterns, pressure drop characteristics and flooding. Int J Multiph Flow 2000;26(9):1439–55. [10] Suzuki S, Ueda T. Behaviour of liquid films and flooding in counter-current twophase flow—part 1. Flow in circular tubes. Int J Multiph Flow 1977;3(6):517–32. [11] Jayanti S, Hewitt GF. Prediction of the slug-to-churn flow transition in vertical twophase flow. Int J Multiph Flow 1992;18(6):847–60. [12] Apanasevich P, Lucas D, Beyer M, Szalinski L. CFD based approach for modeling direct contact condensation heat transfer in two-phase turbulent stratified flows. Int J Therm Sci 2015;95:123–35. [13] Höhne T, Geissler T, Bieberle A, Hampel U. Numerical modeling of a horizontal annular flow experiment using a droplet entrainment model. Ann Nucl Energy 2015;77(77):351–60. [14] Höhne T, Deendarlianto A. Numerical simulations of counter current flow experiments using a morphology detection algorithm. J Comput Multiph Flows 2012;4(3). [15] Chen J, Zeng R, Zhang X, Qiu L, Xie J. Numerical modeling of flow film boiling in cryogenic chilldown process using the AIAD framework. Int J Heat Mass Transf
108