Computational Fluid Dynamics Simulations of the AirâWater Two

Subscriber access provided by TUFTS UNIV

Thermodynamics, Transport, and Fluid Mechanics

CFD simulations of the air-water twophase vertically upward bubbly flow in pipes Dhiraj Ashokrao Lote, Vadakanchery Vinod, and Ashwin Wasudeo Patwardhan Ind. Eng. Chem. Res., Just Accepted Manuscript • DOI: 10.1021/acs.iecr.8b01579 • Publication Date (Web): 13 Jul 2018 Downloaded from http://pubs.acs.org on July 13, 2018

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

CFD simulations of the air-water two-phase vertically upward bubbly flow in pipes Dhiraj A. Lote1, Vadakanchery Vinod2, Ashwin W. Patwardhan1*

1

Department of Chemical Engineering, Institute of Chemical Technology, Matunga, Mumbai-400019, India

2

Indira Gandhi Center for Atomic Research, Kalpakkam, Tamil Nadu-603102, India [email protected], [email protected]

*Corresponding author: Tel: 91-22-33612018; Fax: 91-22-33611020 Email address: [email protected]

1 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Abstract In this study, CFD simulation of air-water two-phase flow in vertical pipe has been carried out. A wide range of superficial liquid and gas velocity 0.376 – 3.53 m/s and 0.0036 – 1.275 m/s and dispersed phase holdup range from 0.3 – 38.74 % is used in the simulation. In this study, for drag, lift, wall lubrication and turbulent dispersion force the model of Grace, Tomiyama, Hosokawa and Burns are used respectively. Eulerian-Eulerian multiphase approach with k-ω SST turbulence model was used. Bubble diameter is an essential parameter to perform the Eulerian simulation. In the present study, we developed new correlation for bubble diameter and with this predicted bubble diameter simulations were performed. The developed CFD model is well capable of predicting gas void fraction, interfacial area concentration, pressure drop, gas/liquid velocity and turbulent kinetic energy. Keywords: Bubble size, interfacial forces, Two phase flow, CFD simulation, Core peak, Wall peak.


Page 2 of 58

Page 3 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1. Introduction Two phase flow is common phenomena which encountered in many industrial processes like nuclear reactors, steam generators, chemical and petrochemical reactors and pipelines. In gas – liquid two phase flow system, gas flows as a dispersed phase in a continuous liquid phase in the form of bubbles. This kind of two phase flow occurs in a vertical pipe flow. Gas and liquid flowing together in a pipe, it takes the possible regimes of bubbly, slug, churn, annular and froth. The distribution of gas void fraction in bubbly flow regime in a pipe takes the possible regimes of wall peak and core peak. These regimes depend on different variables like dispersed phase diameter, gas, and liquid superficial velocity, dispersed phase volume fraction etc. Characteristics of co-current vertical bubbly upflow in pipes such as pressure drop and heat transfer mechanism depend on the distribution of gas void fraction in a pipe. Many researchers have given insight into two phase flow to understand the fundamentals of the radial distribution of gas hold up1,2,3. In the gas-liquid two-phase bubbly flow system, bubbles with smaller size tend to move towards the wall and shows the wall peak phenomenon. This is because the lift coefficient gets the positive value hence, the direction of lift force acting on the bubble is towards the pipe wall. Further, bubbles with larger size tend to move towards center of pipe and shows the core peak phenomenon. In this case, the lift coefficient gets the negative value hence, the direction of lift force acting on the bubble is towards the pipe center. The bubble gets larger in size as the coalescence of smaller bubbles happen. Hence, smaller bubble move towards the wall and bigger bubble move towards the center of the pipe3. From the literature, it is clear that the small bubbles remain near wall region (high void fraction), resulting in wall peak phenomenon and bigger bubbles move toward the pipe centre (high void fraction) showing core peak phenomenon3. Thus, the dynamic pressure profile in a vertical pipe will be non-uniform. If 3 ACS Paragon Plus Environment


Page 4 of 58

we know the void fraction distribution across the pipe, corresponding non-uniform pressure distribution can be predicted. Two-phase density is a function of void fraction. Total pressure drop consist of summation of hydrostatic, frictional and momentum pressure drop. Within these three pressure drop component, the frictional pressure drop is strongly dependent on two phase density. Hence, for accurate prediction of pressure drop, accurate prediction of two phase density is mandatory. For the prediction of two phase density, accurate prediction of gas void fraction is required. Further, from the book “Frontiers and progress in multiphase flow” by Ghajar and Bhagwat4, Two-phase frictional pressure drop is very much sensitive to the gas-liquid distribution across the pipe cross section as well as along the pipe length.

Void fraction

distribution in bubbly flow regime depends on drag and non drag forces. The dispersed gas bubble in a continuous liquid phase experiences the skin drag and form drag. Along with the drag, non drag force also play an important role in the prediction of radial distribution of gas void fraction. Only the drag force cannot give the exact distribution of gas void fraction. As we know the radial distribution of gas void fraction takes the possible regimes of wall peak to core peak. This confirms the lateral movement of the bubble in the vertical pipe. This lateral movement of gas bubble occurs because of the non drag forces which acts perpendicular to the flow direction. The spherical bubble generally experiences the positive lift coefficient. If this is the case then lift force acts towards lower velocity region in a pipe. The situation gets complex when bubble deforms; because of deformation, lift coefficient may be positive or negative. This positive or negative lift coefficient decides the radial movement of the bubble in a pipe. The visualization and comprehensive information through experiment are very difficult. Hence, CFD can be one of the tools to address this difficulty. CFD has a capability to provide detailed, information which might be difficult to get from the experiment. 4 ACS Paragon Plus Environment

Page 5 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Eulerian simulation is performed in this work. The gas void fraction distribution in the vertical pipe is function of bubble size. Smaller bubbles show wall peak phenomenon and bigger bubbles show core peak phenomenon. The shape of a bubble may be spherical or non-spherical depending on Eötvös number (Eo). From the literature, it confirms that the radial movement of bubble occurs as the bubble changes its size3. The bubble size can be measured experimentally or it can be fitted in the simulation to match with the experimental data. Hence, the effort has been made to predict the gas distribution, interfacial area concentration, pressure drop, gas and liquid velocities and turbulent kinetic energy without fitting the bubble diameter. In this study, the correlation is developed to find out the bubble size which is dependent on gas hold up and gas and liquid superficial velocity. The novelty of this work is predictions can be obtained without fitting the bubble diameter. Population balance modelling is one of the ways for the prediction of radial distribution of gas void fraction and gas and liquid velocity. Performing population balance modelling, it is computationally expensive. Without compromising much accuracy, the similar predictions for the radial distribution of gas void fraction, interfacial area concentration, pressure drop, gas and liquid velocity and turbulent kinetic energy can be drawn by a constant bubble size. In this study, some of the recent developments in simulation methods are reviewed. Peña-Monferrer et al.,5 Shang6, and Wang7 used the Eulerian approach with k-ε turbulence model. Islam et al.8 also implemented Eulerian approach with low Re correction turbulence model. Kriebitzsch and Rzehak9, Pellacani et al.,10 Rzehak and Krepper11,12,13, Rzehak et al.14, Rzehak and Kriebitzsch15, and Yamoah et al.16 performed simulations with Eulerian approach using shear stress transport (SST) k-ω turbulence model. Marfaing et al.17 used Neptune code for air-water simulation with k-ε turbulence model.

Montoya et al.18 carried out Eulerian



Page 6 of 58

simulations with SST model. Duan et al.19 in their simulation also implemented k-ω SST model. Wang and Sun20 performed the Eulerian simulation of two phase flow.

In the recent

development of CFD simulation O Marfaing et al.21,22 performed the simulation with small pipe diameter (i.e. 14 mm) by using NEPTUNE code. Recent developments in simulations for airwater two phase flow and present work are tabulated in Table 1. The interfacial force models such as drag, lift, wall lubrication and turbulent dispersion used by various authors are also tabulated in Table 1. Lote et al.53 performed the Eulerian simulation for gas-liquid two-phase flow system in vertical pipe. Author implemented Grace drag model, Tomiyama model for lift and wall lubrication and Burns model for turbulent dispersion. Frank et al.23 performed steady state simulations and concluded that the virtual mass force does not significantly alter the results. We have also tested the influence of virtual mass force on several test cases and found that, the virtual mass force does not influence the results. Hence, virtual mass force is neglected in this study. From the literature, for the prediction of the radial distribution of gas void fraction, researchers have reported the bubble diameters. Either they have experimentally measured or fitted the bubble diameter in simulation. In this study, we developed the new correlation for the prediction of bubble diameter. This predicted bubble diameter is then given as an input to the CFD simulation and simulations were performed. 2. Two fluid model (Euler-Euler approach) In this work, the Eulerian approach is used where each phase treats as interpenetrating continua with the inclusion of gas phase volume fraction concept. In Eulerian approach, it solves momentum and continuity equation separately for both the phases. This is the widely used formulation for gas-liquid two phase flow system24,25. A proper set of conservation equations are required to describe the basic flow, phase interaction, heat, and mass transfer processes. In the 6 ACS Paragon Plus Environment

Page 7 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


momentum equation, the coupling between both the phases is achieved through pressure and interfacial term.

The general conservation equation of mass and momentum are given in

equations (1) and (2) respectively. (q = l for liquid and q = g for gas) Continuity equation

∂ (α q ρ q ) ∂t

+ ∇ ⋅ (α q ρ q v q ) = 0

(1)

Momentum equation

∂ (α q ρ q vq ) ∂t where, ρq,

T + ∇ ⋅ (α q ρ q vq vq ) = −α q ∇P + ∇ ⋅ α G µe ,q (∇vq + ∇vq ) + Fq + α q ρ q g  

(2)

vq and µe,q is the density, velocity and effective viscosity of phase q respectively, g is

the acceleration due to gravity, P is the pressure and, t is the time and Fq is the total interface force. 2.1. Turbulence model Turbulence is a very complex phenomena, prediction of turbulence in a continuous phase is more difficult. The situation gets worse when there is an addition of extra turbulence in a continuous phase (water liquid) due to the presence of bubbles. The shear stress transport (SST) model is used in this work which is developed by Menter26. This model is a blend of two equations the k-ε model and two equation the k-ω model. This model attributes application of the k-ε model in the bulk and application of k-ω model near the wall. In multiphase flow, the turbulence in the continuous phase gets affected by the existence of dispersed phase.

Hence, the fluctuations in the turbulent velocity causes by inherent

turbulence of continuous phase as well as the presence of dispersed phase. assumption, Sato et al.27 derived the turbulent dynamic viscosity of liquid phase.


From this


µ lt = µ lts + µ ltd

Page 8 of 58

(3)

where, the shear induced turbulence µ lts can be expressed as.

µlts = Cµ ρl

kl2

(4)

εl

where, kl is the turbulent kinetic energy and ε l is the dissipation rate of the turbulent kinetic energy. The bubble induced turbulence is given by equation (5).

µ ltd = C µ , g ρ lα g Db v g − vl

(5)

For the gas phase, turbulent viscosity can be expressed as.

µ gt =

ρ g µ lt ρ l Prtg

(6)

where, Prtg is the turbulent Prandtl number of the gas phase. Its value is set to unity.

Cµ and Cµ , g are the model constants which take the value of 0.09 and 0.6 respectively. 2.2. Interface force model For accurate predictions of the radial distribution of gas void fraction, interfacial force play an important role. These forces are drag, lift, wall lubrication, turbulent dispersion and virtual mass force. The virtual mass force is not considered in this study because of the steady state flow condition. The expression for the total interphase force Fq is given in equation (7).


Page 9 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Fq = FD + FL + FWL + FTD + Fv

(7)

In the present work, Grace model for drag, Tomiyama model for lift, Hosokawa model for wall lubrication and Burns model for turbulent dispersion is used.

2.2.1. Drag force The drag force is the resistive force for the motion of bubble, offered by surrounding continuous liquid. The drag force expression for a single bubble is given by equation (8).

π  ρ FD = CD  db2  l (vg − vl ) vg − vl 4  2

(8)

where, CD is the drag coefficient, db is the bubble diameter, and (vg − vl ) is the slip velocity. The drag coefficient CD is different for single bubble and the swarm of the bubble. The drag coefficient CD is Reynolds number dependent in case of the spherical bubble and also Eötvös number (Eo) dependent in case of non-spherical bubbles. In this study, the drag coefficient model of Grace28 was used which takes into account spherical as well as non-spherical shapes. Grace developed the correlation for drag coefficient for different bubble shapes by conducting a dimensionless analysis for individual bubbles rising in a fluid. This model demonstrates the dependency of terminal velocity of a bubble on bubble Reynolds number Reb, Eötvös number (Eo) and Morton number (Mo). The expression for drag coefficient is given by equation (9). At low bubble Reynolds number, bubble remains in the spherical shape and CD can be calculated by CD, sphere given by equation (9). At higher bubble Reynolds number, the bubble may not remain in

spherical shape and it gets deform. This deformation in shape results in the cap bubbly shape or elliptical shape and corresponding drag coefficient is calculated by equation (9).



C D = max(min(C D ,ellipse , C D ,cap ), C D ,sphere )

C D ,sphere

 24  Re  b = 0.687  24 1 + 0.15 Re b  Re b

(

C D ,ellipse =

C D ,cap =

Page 10 of 58

(9)

Re b < 0.01 (10)

)

Re b ≥ 0.01

4 gd b ( ρl − ρ g )

(11)

3U t2 ρl

8 3

(12)

2.2.2. Lift Force A particle flowing in shear flow, experiences transverse force which is called the lift force29. The lift force plays a very important role in the radial distribution of gas void fraction across the pipe. The main cause of lift force is the vortex effect and the velocity gradient in the main flow. The observation of Kariyasaki30 was a spherical particle experience positive lift force coefficient which acts in the direction of decreasing liquid velocity and a deformed bubble influences negative lift force coefficient which acts in the opposite direction of decreasing velocity. Drew and Lahey1 gave a general expression for lift force equation (13). FL = −C L ρ lα g (v g − vl ) × (∇ × vl )

(13)

where, CL is the lift force coefficient. For the ellipsoid and cap bubbly flow regime, Tomiyama et al.31 proposed a lift coefficient model. The lift coefficient dependence on bubble Reynolds number (Reb) and Eötvös number (Eo) (equation (14)). At the low Eötvös number (Eo) bubble maintains the spherical


Page 11 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


shape and moves towards the wall, while at high Eötvös number (Eo) it gets deformed and moves towards the centre of the pipe. The reason is vortex behind deformed bubble becomes slanted and thereby, bubble migrates to the pipe center for high Eötvös number (Eo). Tomiyama studied the sensitivity of contamination in the system (air-water in glycerol solution). Tomiyama performed the experiments with the air-water system and developed the lift coefficient correlation which predicts the lift coefficient for spherical, cap bubbly and elliptical shape bubbles. From the equation (14), it is clear that the lift coefficient can take positive as well as negative value depending on Eötvös number (Eo). At lower Eötvös number (Eo) value, the lift coefficient is positive but at higher Eötvös number (Eo) it changes the sign to negative. More insight on lift force acting on a bubble is demonstrated by Kulkarni Amol.32

[

min 0.288 tanh(0.121Reb ), f ( Eo' )  CL =  f ( Eo' )  − 0.27

]

Eo' ≤ 4 4 < Eo' ≤ 10 10 < Eo

f ( Eo ' ) = 0 .00105 Eo '3 − 0 .0159 Eo ' 2 − 0 .0204 Eo ' + 0 .474

Eo ' =

(14)

'

∆ρgd h2

(15)

(16)

σ

d h = d b (1 + 0.16 Eo 0.757 )1/ 3

(17)

2.2.3. Wall lubrication force There is asymmetric drainage of the liquid surrounding bubble near the wall. At wall side of the bubble, drainage of liquid is less and at opposite side of bubble drainage of liquid is



Page 12 of 58

more. This asymmetric drainage around the bubble generates hydrodynamic pressure difference which drives the bubble away from the wall33. This force is called as a wall lubrication force. 2

FWL = CWL ρlα g (vg − vl )|| nW

(18)

where, CWL is the wall lubrication force coefficient, (v g − vl )|| is the phase relative velocity and

nW is the unit normal pointing away from the wall. Hosokawa et al.34 proposed the correlation for wall lubrication force coefficient (equation (19, 20)). This wall lubrication coefficient depends on bubble Reynolds number Reb as well as on Eötvös number (Eo). Author measured the lateral movement of the bubble near the wall of the vertical plate. The proposed correlation is valid for bubbles as well as for particles. This correlation of wall lubrication coefficient is applicable to both high and low viscosity system. Hosokawa assumes that, the bubble should not colloid with the wall.  7  C w = max  1.9 ,0.0217 Eo   Re b 

CWL

(19)

y   1− w   1 C wc d b  = C w max 0,  Cwd  y y w ( w ) m−1   C wc d b  

(20)

2.2.4. Turbulent dispersion force models The turbulent fluctuations in the liquid velocity and its effect on the disperse phase are accounted by the turbulent dispersion force. Turbulent dispersion force is driven by the void fraction gradient. Turbulent dispersion force flattens the void fraction distribution. Burns et al.35 modeled the turbulent dispersion by performing a time average of the interphase drag term. 12 ACS Paragon Plus Environment

Page 13 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


After double averaged momentum equation it gives additional term which account for the turbulent dispersion force. This is called the Favre Averaged Drag (FAD) model for turbulent dispersion. This Favre averaged drag model for turbulent dispersion can be used with any Reynolds averaged turbulence model and for any arbitrary number of phases with arbitrary morphologies. The expression is given by equation (21).

FTD = −CTD K gl

µlt  ∇α g ∇α l  − ρ l Sc  α g α l 

(21)

here, CTD = 1 and Sc = 0.9

3. Simulation conditions In this work, numerical simulations were performed on adiabatic bubbly flow regime. The simulation conditions given in Table 2 are taken from the available literature. These cases are (1) "wall peak" case where a steep increase in the void fraction is observed near the wall. (2) Less steep or "flattened wall peak" case where a flatten wall peak is observed near the wall. (3) "core peak" case where gas volume fraction is more at the centre of pipe and zero at the wall. Wall peak cases reported by Hibiki et al.,36 Liu,37 Lucas et al.,38 Sun et al.,39 flattened wall peak case measured by Liu,37 Lucas et al.38 and core peak case measured by Fu40 are considered in this study. The details of the experimental conditions for wall peak, flattened wall peak and core peak cases are given in Table 2. An air-water system was used in the simulation with constant material properties. respectively.

The density and viscosity of air was 1.185 kg/m3, 0.0183×10-3 Pa.s

Similarly, the density and viscosity of water was 997 kg/m3, 0.89×10-3 Pa.s

respectively. The constant value of 0.072 N/m of surface tension was used in the simulation. A



fixed bubble diameter was given in the simulation for each case at the pipe inlet i.e. from 1.9 mm to 9.3 mm. These values of bubble diameter were already provided in the literature. The pipe diameters vary from 48.3 mm to 57.2 mm. The length of tubes varies from 2.922 m to 8 m. The superficial velocity of the gas and liquid phases ranges from 0.0111 to 1.275 m/s and 0.19 to 2.607 m/s respectively. Measurement locations for all these cases considered here are also given in Table 2. Velocity inlet, pressure outlet and no slip boundary condition are imposed at inlet, outlet and pipe wall respectively. Bubble diameter (db) is a property of a dispersed phase and used as a constant value in the simulations. Simulations were conducted on the platform of Ansys Fluent (V17.0). Simulations were carried out at the steady state and fully developed flow conditions. Phase coupled SIMPLE algorithm is used for the pressure velocity coupling. The QUICK method was used in the spatial discretization of volume fraction and the second order upwind scheme was used for the rest of the equations. We have identified the cases from the literature and categorized based on wall peak and core peak. From the literature, it has been observed that the wall peak case where the gas void fraction is more towards the wall and for core peak case it is more towards the center of the pipe. The cases which are listed in Table 2, it can be observed that for wall peak cases the bubble diameters are small corresponding gas volume fraction and gas and liquid superficial velocity are also small. On the other hand for core peak case the bubble diameter is large and corresponding gas volume fraction and gas and liquid superficial velocity are also large. Hence, for the development of correlation, the bubble diameter is correlated to gas volume fraction and gas and liquid superficial velocity.


Page 14 of 58

Page 15 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


4. Results and Discussion 4.1. Grid sensitivity study Three dimensional uniform grids are used for the simulation of the adiabatic bubbly flow. Cross sectional view of mesh, as well as flow domain along with boundary conditions are shown in Figure 1 (a, b). To get mesh independence, the numbers of cells in cross section as well as in axial position were varied and corresponding simulations were performed. The details of mesh count and mesh size distribution are given in Table 3. The effect of grid size was studied using 20, 25 and 30 grids in the radial direction and 400, 500 and 600 grids in the axial direction. This grid independence study was carried out for the case of Lucas et al.38 Radial distribution of gas volume fraction and gas velocity was considered for grid sensitivity. Figure 2 (a-b). shows the comparison of different grids on gas volume fraction and gas velocity profile at L/D = 60. It can be seen from Figure 2a. that, predictions obtained with the grid resolution of 1200×400 (480000) (grid size radially and axially 1.28 mm and 8.45 mm, respectively) over predict the gas void fraction at the centre and at the wall. After increasing the grid numbers in radial as well as in axial direction i.e. 1600×500 (800000) (grid size radially and axially 0.85 mm and 6.76 mm respectively) it also overpredicts gas void fraction at the centre and at the wall of a pipe. Another mesh was generated of 2000×600 (1200000) (grid size radially and axially 0.64 mm and 5.63 mm respectively) it also over predicts the gas void fraction. From the investigation, prediction of 1600 × 500 and 2000×600 grids show similar results. Hence, grid of 25 × 500 has been chosen because it satisfies the grid independence study. For the prediction of gas velocity, all the mesh considered here show slightly overprediction (Figure 2b). The grid of 1600 × 500 results in similar prediction of grid 2000×600. Considering the cost of computation time and accuracy, the simulations were performed with grid of 1600 × 500. The simulated result does not 15 ACS Paragon Plus Environment


vary in core region but a small variation is observed near the pipe wall. From Figure 2. it can be seen that grid resolution does not affect much of the radial distribution of void fraction as well as gas velocity.

4.2. Influence of the turbulence model Turbulence is very complex phenomena in two phase flow. The additional turbulence causes by presence of dispersed phase in two phase flow makes it complex. Till date there is no standard in the selection of turbulent model for adiabatic bubbly flow. Hence, the comparison of simulated results for radial distribution of gas void fraction by k-ε, k-ω and RSM model is presented in this study. For the comparisons of these different turbulent models, the wall peak case by Hibiki et al.36 is considered. From the Figure 3. a little difference in the prediction of gas void fraction is seen by the two equation k-ε and k-ω model. The RSM model underpredicts the gas void fraction at the center as well as at the wall. In measurement of radial distribution of gas void fraction, a peak is observed near the wall. On comparing different turbulent models, only the k-ω SST model predicts the peak at wall while other models (k-ε std, k-ε RNG, k-ε Realizable, k-ω std, k-ω BSL and RSM) overpredict the gas fraction at the wall. From the Figure 3. k-ε model overpredict the radial distribution of gas void fraction at the wall which is in well agreement with the simulation results of Frank et al.41 and Cheung et al.42,43. The realistic prediction of gas void fraction is shown by k-ω SST model at close to the wall as well as at the center. Hence, for further investigation the k-ω SST model is considered.

4.3. Interfacial force models In this present work, the cases are selected keeping in mind that it should cover wall peak to core peak phenomena. The behaviour of interfacial force coefficients such as drag, lift, and


Page 16 of 58

Page 17 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


ρ (v − vl )db    . For wall lubrication are different for different bubble Reynolds number  Reb = l g µl   high Reb region, bubbles may get distorted to cap bubbly and ellipsoidal shape. Simple drag correlations are not suitable for such deformed bubbles. Grace drag model considers the bubble distortion and corresponding drag coefficient is calculated. The behaviour of drag coefficient with Reb is given in Figure 4 (a). With an increase in Reb, drag coefficient decreases up to Reb = 400 and then remains constant. With increasing bubble Reynolds number (Reb), initially lift force coefficient increases up to Reb = 30 then it decreases up to Reb = 70 and above Reb = 70 it remains constant (Figure. 4b). In Tomiyama model, lift coefficient depends on bubble Reynolds number (Reb) and Eötvös number (Eo) only. Wall peak to core peak cases with the wide range of bubble Reynolds number (Reb) have been taken to investigate the dependence of lift coefficient on bubble Reynolds number (Reb). In core peak case with high bubble Reynolds number (Reb) and Eötvös number (Eo), the lift coefficient calculated by Tomiyama model changes its sign to negative as seen from Figure. 4(b), hence the bubble moves towards pipe centre. The study of lift force coefficient on bubble diameter has been carried out. It is observed that the reversal of sign of lift coefficient, i.e. from positive (+) to negative (-) occurs at around 5.8 mm bubble diameter (see Figure. 4 (c)). For comparisons of different wall lubrication force coefficient models, three geometries are selected based on pipe diameters (48.3 mm (Fu), 50.8 mm (Hibiki) and 51.2 mm (Lucas)). From the correlation of wall lubrication coefficient CWL,,it is a function of distance from the wall; it only takes effect within a region close enough to the wall. CWL is expressed in term of bubble Reynolds number (Reb), Eötvös number (Eo) and bubble diameter (db). The comparison of wall 17 ACS Paragon Plus Environment


lubrication coefficient is done for Fu, Hibiki and Lucas cases (see Figure. 4 (d-f)). A significant increase of CWL value can be observed near wall region (yw/D ≈ 0). The value of CWL almost decreases to zero near the pipe centre (yw/D ≈ 0.1). In flattened wall peak case and core peak case where the Reb value is high, the effect of CWL values on phase distribution of gas is minimal. If we compare the wall peak and core peak case, in wall peak case the most of the bubbles are concentrated near the wall while at core peak case it is concentrated at the center of the pipe. The effect of wall lubrication force is less for less void fraction near the wall (for core peak).

4.4. Predictions for wall peak cases For the prediction of the radial distribution of gas void fraction, gas and liquid velocity, the experimental data of Hibiki et al.36, Liu,37 Lucas et al.38 and Sun et al.39 are considered. These cases cover gas volume fraction at the inlet from 2.64 % i.e. low gas volume fraction to 10.6 % i.e. high gas volume fraction. The bubble diameter also varies from 1.9 mm to 4.5 mm. In all the cases the bubble diameters taken for simulation was already reported in the literature. The comparison of radial distribution of gas void fraction experimental as well simulations is shown in Figure 5 (a, c, e, g). For the case of Hibiki, a good agreement with the measurement of the radial distribution of gas void fraction is observed at pipe centre and slight under prediction is observed towards the wall (Figure 5(a)). Prediction for the case of Sun is depicted in Figure 5 (c). Prediction of the radial distribution of gas void fraction matches well with experimental data. Prediction for L21B case is shown in Figure 5 (e). Good agreement for the radial distribution of gas void fraction with experimental data is observed in this case. For the case of MT063, predictions are slightly deviating (Figure 5 (g)). This is the case where bigger bubble diameter (4.5 mm) is reported. In the case MT063, small deviation is observed at the wall and at the core. In this case, CFD model under predict the gas void fraction at the center 18 ACS Paragon Plus Environment

Page 18 of 58

Page 19 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


of the pipe while it over predicts gas void fraction at the wall. The root mean square error (RMSE) for Hibiki, Sun, L21B, and MT063 are 0.0164, 0.0037, 0.0143 and 0.0099, respectively. A comparison between liquid and gas axial velocity profiles obtained from simulation with experimental data is shown in Figure 5 (b, d, f, h). For the validation of liquid velocity, cases of Sun and Liu and for the gas velocity cases of Hibiki and MT063 are considered. The experimental result shows that maximum velocity of the liquid is at the centre of the pipe. As we go far from the centre i.e. towards the pipe wall, a steep decrease in velocity is observed (see Figure 5 (b, d, f, h)). The predictions of gas and liquid velocity by using CFD model at the wall as well as at the center is found in good agreement with the experiments. Good agreement with the measurement is observed for the case of Hibiki at the wall as well as at the centre. In case of Sun CFD model slightly underpredicts the liquid velocity at the center. Prediction for L21B and MT063 cases result in good agreement with the experimental data. The RMSE for the cases of Hibiki, Sun, L21B, and MT063 are 0.0199 m/s, 0.0592 m/s, 0.0322 m/s and 0.0650 m/s, respectively.

4.5. Predictions for Flatten wall peak cases For the flattened wall peak phenomenon, the cases of Lucas, L11A, and MT041 are considered. The comparison of gas velocity, liquid velocity and the radial distribution of gas void fraction are shown in Figure 6. Figure 6 (a, c, e) show the comparison of experimental and simulated gas volume fraction distribution radially. For the case of Lucas where the bubble diameter is 4.5 mm, the prediction is in good agreement with experimental data. In this case, CFD model marginally overpredicts gas void fraction at the centre but good agreement with the measurement is observed at the wall. In the case of L11A prediction of a gas void fraction are in good agreement with experimental data in the core region but slightly under prediction is 19 ACS Paragon Plus Environment


observed at the pipe wall. For this case, the bubble diameter is 2.94 mm. For the case of MT041, where the bubble diameter is 4.5 mm distribution of gas void fraction results in good agreement with experimental data at the centre of the pipe but it slightly over predicts at the wall. This may be because of non drag forces acting at the vicinity of the wall. For accurate prediction of the radial distribution of gas volume fraction special attention should be given towards modelling of interfacial forces and more care should be taken in the vicinity of the wall. The RMSE values for Lucas, L11A, and MT041 are 0.0029, 0.0153, and 0.0071, respectively. A comparison between liquid and gas axial velocity profiles obtained from simulation with experimental data is shown in Figure 6 (b, d, f,). For the validation of liquid velocity, the case of Liu and for the gas velocity cases of Lucas and MT041 are considered. Predictions of gas velocity by CFD model elucidate the same pattern what we observed in experimental data. For the case of Lucas, prediction of gas velocity profile is slightly over predicted at the centre as well as at the pipe wall. In the case of L11A, predictions result is in satisfactory agreement with experimental data.

For the case of MT041, predictions are in good agreement with the

experimental data. The RMSE of Lucas, L11A, and MT041 cases are 0.0763 m/s, 0.0833 m/s and 0.0378 m/s, respectively.

4.6. Predictions for core peak case The core peak phenomenon is observed for Fu Run 1 case. In this case, bubble diameter is bigger (9.3 mm) as compared to that observed in wall peak and flattened wall peak case. The superficial velocity is also high in this case. The prediction for gas void fraction and gas/liquid velocity is shown in Figure 7. It can be seen from Figure 4 (c) that lift coefficient changes its sign from positive to negative at 5.8 mm bubble diameter. This lift force on bubble acts in the


Page 20 of 58

Page 21 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


direction of increasing liquid velocity (towards pipe centre) showing core peak phenomenon. From the Figure 7(a), it can be said that for larger bubble the lift coefficient changes its sign to negative, the corresponding void fraction also shows more at the centre and it shows very less at the wall. The prediction for the radial distribution of gas void fraction is in good agreement with the experimental data. The RMSE in Fu Run 1 case is 0.0486. With the available models of interfacial forces, the developed CFD model slightly under predicts gas volume fraction towards the pipe centre and slightly over predict near the wall. Figure 7 (b) shows the comparison of prediction of gas velocity with experimental data. The predictions are in good agreement with the experimental data. Slight overprediction is observed at the center as well as at the pipe wall. The RMSE for the gas velocity is 1.2738 m/s. For the liquid velocity, experimental data was not available in the literature hence only prediction is drawn in Figure 7c. Low slip velocity between the two phases in the wall peak case has little effect on the distortion of the small bubbles. Thus the bubble keeps the spherical shape. As we go from wall peak case to core peak case the bubble diameter keeps on increasing. Correspondingly phase distribution pattern across the pipe changes. For all these cases simulated, results are in good agreement with experimental data at the centre of the pipe but the small deviation is observed near the wall. This is because at near wall; the non drag forces such as lift, wall lubrication, and turbulent dispersion play an important role and its coefficients are very much sensitive to gas void fraction distribution near wall region. Looking at the predictions by developed CFD model, the prediction of the radial distribution of gas void fraction is a strong function of bubble size. As the bubble size changes in the two phase flow in pipe corresponding radial distribution also changes. After summarizing the



results for the prediction of gas void fraction and liquid and gas velocities, the developed CFD model is well capable of predicting the results.

5. Development of correlation Radial distribution of gas void fraction across a pipe strongly depends on bubble size. For smaller bubbles, distribution of gas in the pipe shows wall peak phenomena while with bigger bubble, distribution of gas void fraction across the pipe shows core peak phenomena. A correlation pertaining to the bubble size has been developed in the present study. The developed correlation accounts the effect of gas void fraction (hold up) and gas and liquid superficial velocity. Wall peak and core peak phenomena are dependent on bubble diameter in a pipe. At smaller volume fraction bubbles keep spherical shape but as the volume fraction in a pipe increases, the coalescence frequency also increases which results in the formation of a bigger bubble. The larger bubbles show the core peak phenomena. If there is an increase in gas and liquid velocity, the turbulence in the system increases which results in breakage of bubbles. The breakage of bigger bubble results in the formation of small bubbles. In a vertical pipe with gas liquid flow, small as well as large bubbles are present. Hence the correlation is developed to find out the bubble size which depends on the gas and liquid superficial velocity (JG, JL) and gas void fraction (αg) (see equation 22). In the gas-liquid two-phase system in the vertical pipe the flow regime map is plotted JL V/s JG. Hence, increased in JL and/or JG the flow regime changes from bubbly to slug, slug to churn and churn to annular. In all these flow regimes, the void fraction of dispersed phase (gas) is different. Hence, for the development of correlation the dependence of db on JL, JG and α is considered. For the development of correlation, the experimental data was taken from the available literature. These literature includes the study of Hosokawa and Tomiyama44, Fu40, Lucas et al.38, 22 ACS Paragon Plus Environment

Page 22 of 58

Page 23 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Wang et al.45 and Liu et al.37. The pipe diameter range was considered from 25 mm to 57.2 mm, the superficial gas and liquid velocities are 0.405 to 5.1 m/s and 0.006 to 1.275 m/s respectively and gas void fraction from 0.007 to 0.383. Multiple linear regression analysis was performed to obtain the dependence of JL, JG and αg on db. Multiple linear regression analysis gives the value of multiple R, R square and standard error of 0.8826, 0.7791 and 02189, respectively. The unit of bubble diameter db is in meter and JL and JG are in m/s.

d b = 0.0132 × J L0.9398 × J G−0.6981 × α g0.9719

(22)

5.1. Validation of correlation and predictive nature of CFD model The dispersed phase diameter was calculated by using newly developed correlation and is given input to the CFD model. Iterative method was employed for the prediction of bubble diameter (see Figure 8). At the beginning of the simulation, we have to define the αg at the inlet along with the JL and JG. After the converged solution, we get new αg which is the overall void fraction in the whole domain. Based on this αg value we get the new JL and JG. There is not much significant difference in the new values of JL and JG. Normally, within 3 to 4 runs of simulation we get the constant bubble diameter. The computational time required for one simulation is around 11 hrs. Whereas, for population balance the computational time is around 16 hrs. For the validation of proposed correlations, the cases of Monrós-Andreu et al.,46 Liu and Bankoff47,48 and Hosokawa and Tomiyama49 were considered. The experimental flow conditions for these cases are given in Table 4. Table 2 summarizes the experimental flow conditions which were already available in the literature. In Table 2, the bubble diameter reported is either measured in experiment or it is fitted in the simulations by the authors. The Table 4 demonstrates the experimental flow conditions that are selected for the validation of newly proposed correlation of 23 ACS Paragon Plus Environment


bubble diameter. In this Table all the bubble diameters reported are predicted from the newly proposed correlation by iterative method. As the proposed correlation for bubble diameter is valid for pipe diameter 25 to 57.2 mm ID. The above mention three cases has been selected where in Monrós-Andreu case the pipe diameter is 52 mm, in Liu and Bankoff case pipe diameter is 38 mm and in Hosokawa and Tomiyama case pipe diameter is 25 mm. These three cases intentionally has been selected which covers the pipe diameter from 25 mm to 52 mm. The Monrós Andreu cases with flow condition of VL05JG005TA and VL05JG010TA are considered for the validation of proposed correlation for wall peak case. In these two cases, the range of superficial liquid and gas velocity are 0.5 m/s and 0.05 to 0.1 m/s respectively. The comparisons of prediction by CFD model with the experimental data are shown in Figure 9. for the case of (VL05JG005TA). Figure 9 (a) shows the comparison of the experimental radial distribution of gas void fraction with prediction. A little under prediction is observed for the prediction of gas volume fraction at the center as well as at the pipe wall. Figure 9 (b) shows the comparison of experimental interfacial area concentration (IAC) with prediction by CFD model. The interfacial area concentration profile follows the similar trend of gas void fraction. As little under prediction is observed in the prediction of gas volume fraction same result is replicated in the prediction of interfacial area concentration. CFD model slightly underpredict the interfacial area concentration at the wall as well as at the center. The comparison of gas axial velocity with the prediction is shown in Figure 9 (c). The prediction of gas velocity is in excellent agreement with the experimental data. Figure 9 (d). shows the comparison of experimental liquid axial velocity with prediction. A good agreement is observed between CFD model and experimental data for liquid velocity. The RMSE for radial distribution of gas void fraction, axial gas and liquid velocity are 0.0375, 0.1420 m/s and 0.0221 m/s, respectively. 24 ACS Paragon Plus Environment

Page 24 of 58

Page 25 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 9 (e, f, g, h) show the comparisons of experimental data of radial distribution of gas void fraction, interfacial area concentration, axial gas and liquid velocity for the case of (VL05JG010TA). For the prediction of gas volume fraction, CFD model slightly underpredicts the gas volume fraction at the center as well as at the wall (Figure 9 e). The similar results are replicated in the prediction of interfacial area concentration.

The developed CFD model

underpredicts the interfacial area concentration at the center as well as at the wall (Figure 9 f). A good agreement with the measurement is observed in the prediction of gas velocity. CFD model slightly under predict the gas velocity at the center and slightly over predict at the wall (Figure 9 g). An excellent agreement is observed in the prediction of liquid velocity (Figure 9 h). The RMSE for the radial distribution of gas void fraction, axial gas and liquid velocity are 0.0193, 0.0921 m/s and 0.0551 m/s, respectively. The validity of the above developed correlation is also checked for the core peak case measured by Monrós Andreu (VL10JG030TA). The flow conditions for this case are given in Table 4. In this case, the liquid and gas superficial velocity are 1 and 0.3 m/s respectively. For this case, only the radial distributions of gas void fraction, interfacial area concentration and gas velocity are available in the literature. Hence, for the liquid velocity only predictions are drawn (see Figure (10d)). With the newly developed correlation (equation (22)) bubble diameter was calculated, then simulations were performed by iterative method and predictions are drawn. The developed CFD model results in little under prediction of gas void fraction and interfacial area concentration (Figure 10 a, b). Good agreement with the measurement is observed in the prediction of gas velocity (Figure 10 c). The experimental data for the liquid velocity was not reported in this case hence, only the prediction of liquid velocity is plotted (Figure 10 d). The



RMSE for the radial distribution of gas void fraction and gas velocity are 0.0969 and 0.1934 m/s, respectively. The validity of the newly developed correlation and predictability of the developed CFD model has been tested for the lower pipe diameter cases which include the case of Liu and Bankoff and Hosokawa and Tomiyama. The experimental flow conditions for these cases are given in Table 4. The predictions for the radial distribution of gas void fraction, gas and liquid velocity for the case of Liu and Bankoff is shown in Figure 11. In this case the JL and JG are 0.376 m/s and 0.347 m/s respectively. From Figure 11 a. CFD model slightly underpredict the gas void fraction at the center as well as at the wall. For the prediction of liquid velocity, the developed model slightly overpredicts the liquid velocity at the wall and it slightly under predict at the center of pipe (Figure 11b). The predictions are in good agreement with the experimental data of gas velocity (Figure11c). The RMSE for the radial distribution of gas void fraction, gas velocity and liquid velocity are 0.1039, 0.0996 m/s and 0.1295 m/s, respectively. In the next case of the Liu and Bankoff, the JL is 1.087 m/s and JG is 0.027 m/s. In this case the superficial liquid velocity is somewhat high and superficial gas velocity is somewhat less as compare to first case. Predictions for the radial distribution of gas volume fraction, gas and liquid velocity is shown in Figure 11 (d-f). The prediction of gas volume fraction and liquid velocity matches well with the experimental data (Figure 11 d-e). The experimental data of gas velocity was not reported in the literature, hence only the predictions are drawn in Figure 11 f. The RMSE for the radial distribution of gas void fraction and liquid velocity are 0.0082 and 0.0477 m/s, respectively. The prediction for gas volume fraction, liquid velocity and turbulent kinetic energy is also drawn for lowest pipe diameter case i.e. 25 mm ID (Hosokawa and Tomiyama case). Two 26 ACS Paragon Plus Environment

Page 26 of 58

Page 27 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


cases are selected from the experimental data measured by Hosokawa and Tomiyama. The experimental flow conditions are also given in Table 4. In the first case, JL is 1 m/s and JG is 0.02 m/s. In the prediction of radial distribution of gas void fraction, the developed CFD model results in good agreement with the experimental data at the center of pipe and a little over prediction is observed near wall (Figure 12a). The developed CFD model shows good results for liquid velocity with slight overprediction (Figure 12 b).

An excellent agreement with the

experimental data of turbulent kinetic energy is observed towards the center of the pipe while little under prediction is observed at the wall (Figure 12 c). The RMSE for radial distribution of gas void fraction, liquid velocity and the turbulent kinetic energy are 0.0124, 0.1370 m/s and 0.0011 (m/s)2, respectively. In the next case of Hosokawa and Tomiyama, the JL and JG are 0.5 m/s and 0.018 m/s respectively. The predictions for radial distribution of gas void fraction, liquid velocity and turbulent kinetic energy is plotted in Figure 12 d-f. The developed CFD model under predict the gas void fraction at the center and it over predict at the wall (Figure 12 a). The under prediction is also observed for liquid velocity at the center as well as at the wall (Figure 12 e). The developed CFD model under predict the turbulent kinetic energy (Figure 12 f). The RMSE for the radial distribution of gas void fraction, liquid velocity and the turbulent kinetic energy are 0.0183, 0.4004 m/s and 0.0025 (m/s)2, respectively.

5.2. Pressure drop predictions The applicability of the developed correlation and the predictive nature of CFD model have been checked for the pressure drop. As the developed correlation is valid for the pipe diameter ranges from 25 mm to 57.2 mm ID. The cases for the prediction of pressure drop is considered which fall in the above mention range of pipe diameter. Three cases have been selected from the available literature. These three cases are the measurement of pressure drop by 27 ACS Paragon Plus Environment


Page 28 of 58

Shiba and Yamazaki50 (pipe ID - 24.5 mm), measurement by Hernandez Perez51 (Pipe ID – 38 mm) and measurement by Descamps et al.52 (pipe ID – 50 mm).

The experimental flow

conditions are mentioned in the Table S1 of supporting information. In the case of the pressure drop measured by Descamps, the prediction and the measured pressure drop values are shown in Figure 13 a. In this case the pipe diameter is 50 mm and the total length is 7 m. The differential pressure was measured by Descamps at the length of 2 m i.e. ∆L = 2 m. For the prediction of pressure drop same process is followed as mention in Figure 8. From the Figure 13 a. it can be seen that a good agreement of pressure drop with experimental data is observed at lower gas superficial velocity while at higher gas superficial velocity a little deviation is observed. Prediction of pressure drop for 38 mm pipe ID is shown in Figure 13 b. Hernandez Parez measured the pressure drop in a pipe having 6.5 m in length. The pressure drop was measured at the length of 0.75 m i.e. ∆L = 0.75 m. In this case pressure drop was measured at constant liquid superficial velocity of 0.73 m/s and at varying gas superficial velocity. The comparison of predicted and measured pressure drop is shown in Figure 13 b. The developed CFD model slightly under predict the pressure drop. Figure 13 c. shows the prediction of pressure drop for the smallest pipe ID (24.5 mm). The pressure drop was measured by Masayoshi Shiba in a pipe having length 1.35 m. The differential pressure was measured at the length of 0.5 m i.e. ∆L = 0.5 m. The pressure drop was measured at constant liquid superficial velocity of 1.94 and 3.53 m/s and at varying gas superficial velocity. A good agreement is observed between the prediction and measurement of pressure drop data.


Page 29 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


In order to validate the newly developed correlation for bubble diameter and CFD model, the simulations were performed for the prediction of pressure drop. The predictions for pressure drop are carried out for smaller pipe diameter (ID – 24.5 mm), moderate pipe diameter (ID – 38 mm) and higher pipe diameter (ID – 50 mm). The comparison of prediction of pressure drop and measurement of pressure drop is shown in Figure 13 d. for all the pipe diameter. The prediction of pressure drop results in ±15 % error with the measured pressure drop data.

6. Conclusions The applicability of interphase force models depends on bubbly flow regime (wall peak and core peak). The interfacial forces; drag, lift, wall lubrication and turbulent dispersion are considered in this study. The simulation validations with these interfacial forces are in good agreement with experimental data. Grace drag model, Tomiyama lifts model, Hosokawa wall lubrication model and Burns model for turbulent dispersion was used in this study.

The

simulation predictions of the radial distribution of gas void fraction, axial liquid/gas velocity, interfacial area concentration, turbulent kinetic energy and pressure drop are compared with experimental data. In the transverse movement of a bubble, non drag forces such as lift, wall lubrication, and turbulent dispersion play the significant role. With this consent "wall peak" and "core peak" phenomenon is observed experimentally from literature as well as predictions in the simulation. For the simulation validation, the cases with gas volume fraction range from low gas void fraction i.e. 0.3% to high gas void fraction i.e. 38.74% were considered. The simulation results are in good agreement with these low and high gas volume fraction experimental data. A correlation has been developed to find out bubble diameter which is gas and liquid superficial velocity (JG, JL) and gas hold up (αg) dependent. The predictive ability of this correlation was tested for a wall as well as for a core peak case. For both the cases, the developed 29 ACS Paragon Plus Environment


correlation gives satisfactory results. The predictions of the radial distribution of gas void fraction, axial velocity of liquid and gas phase, interfacial area concentration and turbulent kinetic energy were checked by developed CFD model.

It shows a good agreement with

experimental data in all the cases. The applicability of the newly developed correlation for bubble diameter and predictive nature of developed CFD model is also checked for the pressure drop predictions. The prediction of pressure drop results in ± 15 % error with measured pressure drop.

Acknowledgment: The authors would like to thank Indira Gandhi Center for Atomic Research for funding the research project (Sanction number: IGCAR/FRTG/SE & HD/ETHS/SG/06/2015). One of the authors Dhiraj A. Lote would like to thank the University Grant Commission (UGC) for providing financial support and DST-FIST for providing high performance computing system.

Nomenclature v

Velocity vector (m/s)

P

Pressure (Pa)

g

Acceleration due to gravity (m/s2)

F

Interphase force (N)

t

Time (s)

k

Turbulence kinetic energy (m2/s2)

C

Coefficient

db

Bubble diameter (m)

Reb

Reynolds number

Eo

Eötvös number

Mo

Morton number

Ut

Terminal velocity (m/s)

Sc

Schmidt number


Page 30 of 58

Page 31 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


J

Superficial velocity (m/s)

L

Length (m)

yw

Distance to the nearest wall (m)

Z

Measurement location (m)

RSM

Reynolds Stress Model

Greek Symbols

α

Volume fraction

ρ

Density (kg/m3)

ε

Turbulence dissipation rate (m2/s3)

ω µ

Specific turbulence dissipation rate (1/s)

σ

Surface tension coefficient (N/m)

q

Phase symbol

g

Gas phase

l

Liquid phase

gl

Transfer of quantities between gas and liquid

dr

Drift

b

Bubble

D

Drag

L

Lift

Viscosity (m2/s)

Subscripts

WL

Wall lubrication

TD

Turbulent dispersion

w

Wall

h

Horizontal

t

Terminal



Supporting information: Table S1. Flow conditions for the validation of correlations for the prediction of pressure drop.

References (1)

Drew, D. A.; Lahey, R. T. The Virtual Mass and Lift Force on a Sphere in Rotating and Straining Inviscid Flow. Int. J. Multiph. Flow 1987, 13 (1), 113–121.

(2)

Lucas, D.; Krepper, E.; Prasser, H.-M. Use of Models for Lift, Wall and Turbulent Dispersion Forces Acting on Bubbles for Poly-Disperse Flows. Chem. Eng. Sci. 2007, 62 (15), 4146–4157.

(3)

Tomiyama, A. Struggle with Computational Bubble Dynamics. Multiphase Science and Technology. 1998, pp 369–405.

(4)

Ghajar, A. J.; Bhagwat, S. M. Flow Patterns, Void Fraction and Pressure Drop in GasLiquid Two Phase Flow at Different Pipe Orientations. In Frontiers and Progress in Multiphase Flow I; Springer, 2014; pp 157–212.

(5)

Peña-Monferrer, C.; Passalacqua, A.; Chiva, S.; Muñoz-Cobo, J. L. CFD Modelling and Validation of Upward Bubbly Flow in an Adiabatic Vertical Pipe Using the Quadrature Method of Moments. Nucl. Eng. Des. 2016, 301, 320–332.

(6)

Shang, Z. A Novel Drag Force Coefficient Model for Gas--Water Two-Phase Flows under Different Flow Patterns. Nucl. Eng. Des. 2015, 288, 208–219.

(7)

Wang, X. Simulations of Two-Phase Flows Using Interfacial Area Transport Equation. PhD 2010, 225.

(8)

Islam, A. S. M. A.; Adoo, N. A.; Bergstrom, D. J. Prediction of Mono-Disperse Gas– liquid Turbulent Flow in a Vertical Pipe. Int. J. Multiph. Flow 2016, 85, 236–244.


Page 32 of 58

Page 33 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(9)

Kriebitzsch, S.; Rzehak, R. Baseline Model for Bubbly Flows: Simulation of Monodisperse Flow in Pipes of Different Diameters. Fluids 2016, 1 (3), 29.

(10)

Pellacani, F.; Chiva, S.; Peña, C.; Macián, R. Implementation of a One-Group Interfacial Area Transport Equation in a CFD Code for the Simulation of Upward Adiabatic Bubbly Flow.

(11)

Rzehak, R.; Krepper, E. Bubble-Induced Turbulence: Comparison of CFD Models. Nucl. Eng. Des. 2013, 258, 57–65.

(12)

Rzehak, R.; Krepper, E. CFD Modeling of Bubble-Induced Turbulence. Int. J. Multiph. Flow 2013, 55, 138–155.

(13)

Rzehak, R.; Krepper, E. Closure Models for Turbulent Bubbly Flows: A CFD Study. Nucl. Eng. Des. 2013, 265, 701–711.

(14)

Rzehak, R.; Krepper, E.; Liao, Y.; Ziegenhein, T.; Kriebitzsch, S.; Lucas, D. Baseline Model for the Simulation of Bubbly Flows. Chem. Eng. Technol. 2015, 38 (11), 1972– 1978.

(15)

Rzehak, R.; Kriebitzsch, S. Multiphase CFD-Simulation of Bubbly Pipe Flow: A Code Comparison. Int. J. Multiph. Flow 2015, 68, 135–152.

(16)

Yamoah, S.; Martínez-Cuenca, R.; Monrós, G.; Chiva, S.; Macián-Juan, R. Numerical Investigation of Models for Drag, Lift, Wall Lubrication and Turbulent Dispersion Forces for the Simulation of Gas-Liquid Two-Phase Flow. Chem. Eng. Res. Des. 2015, 98, 17– 35.

(17)

Marfaing, O.; Guingo, M.; Laviéville, J.; Bois, G.; Méchitoua, N.; Mérigoux, N.; Mimouni, S. An Analytical Relation for the Void Fraction Distribution in a Fully Developed Bubbly Flow in a Vertical Pipe. Chem. Eng. Sci. 2016, 152, 579–585.



(18)

Montoya, G.; Baglietto, E.; Lucas, D.; Krepper, E.; Hoehne, T. Comparative Analysis of High Void Fraction Regimes Using an Averaging Euler-Euler Multi-Fluid Approach and a Generalized Two-Phase Flow (GENTOP) Concept. In Proceedings of the 2014 22st International Conference on Nuclear Engineering (ICONE22). Prague, Czech Republic; 2014.

(19)

Duan, X. Y.; Cheung, S. C. P.; Yeoh, G. H.; Tu, J. Y.; Krepper, E.; Lucas, D. Gas--Liquid Flows in Medium and Large Vertical Pipes. Chem. Eng. Sci. 2011, 66 (5), 872–883.

(20)

Wang, X.; Sun, X. Effects of Non-Uniform Inlet Boundary Conditions and Lift Force on Prediction of Phase Distribution in Upward Bubbly Flows with Fluent-IATE. Nucl. Eng. Des. 2011, 241 (7), 2500–2507.

(21)

Marfaing, O.; Guingo, M.; Laviéville, J.-M.; Mimouni, S. Analytical Void Fraction Profile Near the Walls in Low Reynolds Number Bubbly Flows in Pipes: Experimental Comparison and Estimate of the Dispersion Coefficient. Oil Gas Sci. Technol. d�IFP Energies Nouv. 2017, 72 (1), 4.

(22)

Marfaing, O.; Laviéville, J. Bubble Force Balance Formula for Low Reynolds Number Bubbly Flows in Pipes and Channels: Comparison of Wall Force Models. Int. J. Chem. React. Eng. 2017.

(23)

Frank, T.; Zwart, P. J.; Krepper, E.; Prasser, H.-M.; Lucas, D. Validation of CFD Models for Mono-and Polydisperse Air--Water Two-Phase Flows in Pipes. Nucl. Eng. Des. 2008, 238 (3), 647–659.

(24)

Drew, D. A.; Lahey, R. T. Application of General Constitutive Principles to the Derivation of Multidimensional Two-Phase Flow Equations. Int. J. Multiph. Flow 1979, 5 (4), 243–264.


Page 34 of 58

Page 35 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(25)

Ishii, M.; Hibiki, T. Thermo-Fluid Dynamics of Two-Phase Flow; Springer Science & Business Media, 2010.

(26)

Menter, F. R. Two-Equation Eddy-Viscosity Turbulence Models for Engineering Applications. AIAA J. 1994, 32 (8), 1598–1605.

(27)

Sato, Y.; Sadatomi, M.; Sekoguchi, K. Momentum and Heat Transfer in Two-Phase Bubble Flow-I. Theory. Int. J. Multiph. Flow 1981, 7 (2), 167–177.

(28)

Grace, J. R. Shapes and Velocities of Single Drops and Bubbles Moving Freely through Immiscible Liquids. Chem. Eng. Res. Des. 1976, 54, 167–173.

(29)

Auton, T. R. The Lift Force on a Spherical Body in a Rotational Flow. J. Fluid Mech.

1987, 183, 199–218. (30)

Kariyasaki, A. Behavior of a Single Gas Bubble in a Liquid Flow with a Linear Velocity Profile. In Proceedings of the 1987 ASME-JSME Thermal Engineering Joint Conference; 1987; Vol. 5, pp 261–267.

(31)

Tomiyama, A.; Tamai, H.; Zun, I.; Hosokawa, S. Transverse Migration of Single Bubbles in Simple Shear Flows. Chem. Eng. Sci. 2002, 57 (11), 1849–1858.

(32)

Kulkarni, A. A. Lift Force on Bubbles in a Bubble Column Reactor: Experimental Analysis. Chem. Eng. Sci. 2008, 63 (6), 1710–1723.

(33)

Antal, S. P.; Lahey, R. T.; Flaherty, J. E. Analysis of Phase Distribution in Fully Developed Laminar Bubbly Two-Phase Flow. Int. J. Multiph. Flow 1991, 17 (5), 635– 652.

(34)

Hosokawa, S.; Tomiyama, A.; Misaki, S.; Hamada, T. Lateral Migration of Single Bubbles due to the Presence of Wall. ASME 2002 Jt. US-European Fluids Eng. Div. Conf.

2002, 4 (4), 855–860.



(35)

Burns, A. D.; Frank, T.; Hamill, I.; Shi, J.-M.; others. The Favre Averaged Drag Model for Turbulent Dispersion in Eulerian Multi-Phase Flows. In 5th international conference on multiphase flow, ICMF; 2004; Vol. 4, pp 1–17.

(36)

Hibiki, T.; Ishii, M.; Xiao, Z. Axial Interfacial Area Transport of Vertical Bubbly Flows. Int. J. Heat Mass Transf. 2001, 44 (10), 1869–1888.

(37)

Liu, T. J. The Role of Bubble Size on Liquid Phase Turbulent Structure in Two-Phase Bubbly Flow. In Proc. Third International Congress on Multiphase Flow ICMF; 1998; Vol. 98, pp 8–12.

(38)

Lucas, D.; Krepper, E.; Prasser, H.-M. Development of Co-Current Air--Water Flow in a Vertical Pipe. Int. J. Multiph. Flow 2005, 31 (12), 1304–1328.

(39)

Sun, X.; Paranjape, S.; Kim, S.; Ozar, B.; Ishii, M. Liquid Velocity in Upward and Downward Air--Water Flows. Ann. Nucl. Energy 2004, 31 (4), 357–373.

(40)

Fu, X. Interfacial Area Measurement and Transport Modeling in Air-Water Two-Phase Flow. 2001.

(41)

Frank, T.; Shi, J.-M.; Burns, A. D. Validation of Eulerian Multiphase Flow Models for Nuclear Safety Application. 3rd Int. Symp. Two-Phase Flow Model. Exp. Pisa, 22-24 Sept. 2004 2004, No. September, 22–24.

(42)

Cheung, S. C. P.; Yeoh, G. H.; Tu, J. Y. On the Numerical Study of Isothermal Vertical Bubbly Flow Using Two Population Balance Approaches. Chem. Eng. Sci. 2007, 62 (17), 4659–4674.

(43)

Cheung, S. C. P.; Yeoh, G. H.; Tu, J. Y. On the Modelling of Population Balance in Isothermal Vertical Bubbly Flows�average Bubble Number Density Approach. Chem. Eng. Process. Process Intensif. 2007, 46 (8), 742–756.


Page 36 of 58

Page 37 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(44)

Hosokawa, S.; Tomiyama, A. Multi-Fluid Simulation of Turbulent Bubbly Pipe Flows. Chem. Eng. Sci. 2009, 64 (24), 5308–5318.

(45)

Wang, S. K.; Lee, S. J.; Jones Jr, O. C.; Lahey Jr, R. T. 3-D Turbulence Structure and Phase Distribution Measurements in Bubbly Two-Phase Flows. Int. J. Multiph. flow 1987, 13 (3), 327–343.

(46)

Monrós-Andreu, G.; Chiva, S.; Mart’\inez-Cuenca, R.; Torró, S.; Juliá, J. E.; Hernández, L.; Mondragón, R. Water Temperature Effect on Upward Air-Water Flow in a Vertical Pipe: Local Measurements Database Using Four-Sensor Conductivity Probes and LDA. In EPJ Web of Conferences; 2013; Vol. 45, p 1105.

(47)

Liu, T. J.; Bankoff, S. G. Structure of Air-Water Bubbly Flow in a Vertical pipe�II. Void Fraction, Bubble Velocity and Bubble Size Distribution. Int. J. Heat Mass Transf. 1993, 36 (4), 1061–1072.

(48)

Liu, T. J.; Bankoff, S. G. Structure of Air-Water Bubbly Flow in a Vertical pipe�I. Liquid Mean Velocity and Turbulence Measurements. Int. J. Heat Mass Transf. 1993, 36 (4), 1049–1060.

(49)

Hosokawa, S.; Tomiyama, A. Multi-Fluid Simulation of Turbulent Bubbly Pipe Flows. Chem. Eng. Sci. 2009, 64 (24), 5308–5318.

(50)

SHIBA, M.; YAMAZAKI, Y. A Comparative Study on the Pressure Drop of Air-Water Flow. Bull. JSME 1967, 10 (38), 290–298.

(51)

Hernandez Perez, V. Gas-Liquid Two-Phase Flow in Inclined Pipes, University of Nottingham, 2008.

(52)

Descamps, M. N.; Oliemans, R. V. A.; Ooms, G.; Mudde, R. F. Air-Water Flow in a Vertical Pipe: Experimental Study of Air Bubbles in the Vicinity of the Wall. Exp. Fluids



2008, 45 (2), 357–370. (53)

Lote, D. A.; V. Vinod; Patwardhan. A. W. Comparison of Models for Drag and NonDrag Forces for Gas-Liquid two-phase bubbly flow. Multiphase Science and Technology

2018 (DOI: 10.1615/MultScienTechn.2018025983)

List of Figures: Figure 1. (a) Cross sectional view of mesh and (b) flow domain with boundary condition of a vertical pipe. Figure 2. Effect of grid size on radial distribution of gas void fraction (a) and gas velocity (b). Figure 3. Prediction of radial distribution of gas void fraction by different turbulence models. Figure 4. (a) Drag coefficient (CD) Vs bubble Reynolds number (Reb), (b) variation of lift coefficient (CL) with bubble Reynolds number (Reb), (c) Lift coefficient Vs bubble diameter. Distribution of wall lubrication force coefficient (CWL) Vs yw/D for (d) wall peak case (e) flattened wall peak case and (f) core peak case. Figure 5. Simulated results of radial distribution of gas void fraction case for (a) Hibiki, (c) Sun, (e) L21B and (g) MT063 case and velocity for (b) Hibiki, (d) Sun, (f) L21B and (h) MT063 case. Figure 6. Simulated results of radial distribution of gas void fraction for (a) Lucas, (c) L11A, and (e) MT041 case and velocity for (b) Lucas, (d) L11A, and (f) MT041 case. Figure 7. Simulated results of radial distribution of gas void fraction (a), gas velocity (b) and liquid velocity (c) for "core peak" case of Fu Run 1 Figure 8. Iterative procedure for the calculation of bubble diameter. Figure 9. Prediction of radial distribution of gas void fraction (a), IAC (b), axial gas (c) and liquid (d) velocity for experimental data of Andreu case (VL05JG005TA) and gas void fraction (e), IAC (f), axial gas (g) and axial liquid (h) velocity for Andrew case (VL05JG010TA). Figure 10. Prediction of radial distribution of gas void fraction (a), IAC (b), and axial gas velocity (c) and liquid velocity (d) by new proposed correlation for experimental data of Andreu case (VL10JG030TA). Figure 11. Prediction of radial distribution of gas void fraction (a), axial liquid(b) and gas velocity (c) for JL-0.376, JG-0.347 and radial distribution of gas void fraction (d), axial liquid (e) and gas (f) velocity for JL-1.087, JG-0.027 for experimental data of Liu and Bankoff. Figure 12. Prediction for radial distribution of gas void fraction (a), liquid velocity (b) and turbulent kinetic energy (c) for JL-1 and JG-0.02 and radial distribution of gas void fraction (d), 38 ACS Paragon Plus Environment

Page 38 of 58

Page 39 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


liquid velocity (e) and turbulent kinetic energy (f) for JL-0.5 and JG-0.018 for the experimental data of Hosokawa and Tomiyama.

Figure 13. Gas-liquid pressure drop as a function of gas superficial velocity for (a) pipe ID – 50 mm, (b) pipe ID 38 mm, (c) pipe ID 24.5 mm, and (d) comparison of total pressure drop between measured and predicted values in vertical gas-liquid two phase flow.

List of Tables:

Table 1. Summery of previous literature on simulation flow conditions for bubbly, slug, churn and annular flow regimes. Table 2. Flow conditions for the wall peak, flattened wall peak and core peak case. Letters "MT" denotes the MT Loop test of (Prasser et al., 2003; Lucas et al., 2005) while letter "L" denotes the test of Liu (1998). Table 3. Different grids used in the simulation of two phase pipe flow. Table 4. Flow conditions for the validation of correlations for wall peak and core peak cases.



Figure 1. (a) Cross sectional view of mesh and (b) flow domain with boundary condition of a vertical pipe.


Page 40 of 58

Page 41 of 58

0.018

480000mesh 800000mesh 1200000mesh

Gas volume fraction (-)


0.04

0.015 0.65

0.02

Radial position r/R (-)

(a)

0.66

Expt

480000mesh 800000mesh 1200000mesh 0 0

0.2

0.4 0.6 Radial position r/R (-)

0.8

1

1.6 Gas velocity (m/s)

(b) 1.2

Gas velocity (m/s)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


0.633

0.631 0.88

480000mesh 800000mesh 1200000mesh r/R

0.882

0.8 Expt

0.4

480000mesh 800000mesh 1200000mesh

0 0

0.2

0.4 0.6 Radial position r/R

0.8

1

Figure 2. Effect of grid size on radial distribution of gas void fraction (a) and gas velocity (b). 41 ACS Paragon Plus Environment


0.2

Expt k-ε-std k-ε-RNG Gas void fraction (-)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 42 of 58

k-ε-Realizable k-ω-std k-ω-sst k-ω-BSL RSM

0.1

0 0

0.2

0.4

0.6

0.8

1

Radial Position (r/R)

Figure 3. Prediction of radial distribution of gas void fraction by different turbulence models.


Page 43 of 58

2000

100

Grace

(a)

Hibiki et al.

(d)

CD (-)

10

CWL (-)

0.1

0 1

10

100 Reb

1000

0.5

Tomiyama

0.25

yw / D = 0.019

0 10000 2000

0.02

0.04

0.06 yw/D (-)

0.08

0.1

Lucas et al.

(b)

1000

0

-0.25

(e)

First node's position

CWL (-)

CL (-)

First node's position

1000 1

yw / D = 0.021

0

-0.5 1

10

100 Reb

0.5

1000

0

10000

0.02

2000

Tomiyama

0.04 0.06 yw/D (-)

0.08

Fu et al.

0.1

(f)

(c)

0.25 1000

0

-0.25

CWL (-)

CL (-)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


-0.5

0

0

0.003

0.006

0.009

0.012

First node's position yw / D = 0.016

0

Bubble diameter (m)

0.02

0.04 0.06 yw/D (-)

0.08

Figure 4. (a) Drag coefficient (CD) Vs bubble Reynolds number (Reb), (b) variation of lift coefficient (CL) with bubble Reynolds number (Reb), (c) Lift coefficient Vs bubble diameter. Distribution of wall lubrication force coefficient (CWL) Vs yw/D for (d) wall peak case (e) flattened wall peak case and (f) core peak case.


0.1


2 (a)

Expt

Gas vrelocity (m/s)


0.2

Predicted 0.1

0 0.2

0.4 0.6 0.8 Radial Position (r/R)

(b)

Predicted 1

Expt

1

0.2

2 (c)

Predicted

0.1

0 Liquid velocity (m/s)

0.15


Expt

0 0

0.05

0.4 0.6 Radial position (r/R)

0.8

Expt

1

(d)

Predicted 1

0

0 0

0.2

0.4 0.6 Radial Position (r/R)

0.8

0

1

0.2


0.8

1

Expt

0.4

Liquid velocity (m/s)


2 (e)

Predicted 0.2

0

(f)

Expt Predicted 1

0 0

0.2

0.15

0.4 0.6 0.8 Radial position (r/R)

Expt

1

0

(g) Gas velocity (m/s)


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 44 of 58

Predicted 0.1

0.05

0

0.2


0.8

Expt 2

1

(h)

Predicted

1

0 0

0.2


0.8

1

0

0.2


0.8

Figure 5. Simulated results of radial distribution of gas void fraction case for (a) Hibiki, (c) Sun, (e) L21B and (g) MT063 case and velocity for (b) Hibiki, (d) Sun, (f) L21B and (h) MT063 case. 44 ACS Paragon Plus Environment

1

Page 45 of 58

1.5 (a)

Gas velocity (m/s)


0.09 Expt 0.06 Predicted 0.03

1

Predicted

0.5

0 0

0.2


0.8

0

1

0.4

0.2


0.8

1

1.5 (c)

Expt



(b)

Expt

0

Predicted 0.2

(d)

Expt 1

Predicted

0.5

0

0 0

0.2


0.8

1

0

0.09

0.2


0.8

1

3 (e)

Expt

Expt

Gas velocity (m/s)


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


0.06 Predicted 0.03

0

(f)

Predicted

2

1

0 0

0.2


0.8

1

0

0.2


0.8

Figure 6. Simulated results of radial distribution of gas void fraction for (a) Lucas, (c) L11A, and (e) MT041 case and velocity for (b) Lucas, (d) L11A, and (f) MT041 case.


1



1 Expt

(a)

Predicted 0.5

0 0

0.2


0.8

1

10 (b)

Gas velocity (m/s)

Expt Predicted 5

0 0

0.2


0.8

1

10


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 46 of 58

Predicted

(c)

5

0 0

0.2


0.8

1

Figure 7. Simulated results of radial distribution of gas void fraction (a), gas velocity (b) and liquid velocity (c) for "core peak" case of Fu Run 1. 46 ACS Paragon Plus Environment

Page 47 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Get the JL, JG and αg from experiments

Calculate bubble diameter (db) by putting JL, JG and αg in the equation (21)

Perform the simulation with calculated bubble diameter (db) and JL, JG and αg

After the converged solution, get the αg, JL, JG. Calculate the (db)

Is the bubble diameter same as previous (db)

No

Yes Plot the αg distribution, IAC, Pressure drop, ug, ul and k.

Figure 8. Iterative procedure for the calculation of bubble diameter.



0.6 (a)

Expt Predicted

0.2



0.3

0.1

(e) Expt Predicted

0.4

0.2

0

0 0

0.01 0.02 Radial position r, (m)

0


400 (b) Expt

Predicted

IAC (1/m)

IAC (1/m)

200

(f)

Expt

600

Predicted

300

0

0 0

0

0.01 0.02 Radial position r, (m) (c)

Expt Predicted

1

0.5

Gas vvelocity (m/s)

Gas vvelocity (m/s)


1.5

1.5

Expt

(g)

Predicted

1

0.5

0

0 0

0


1


1 (d)

Expt

Predicted 0.5



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 48 of 58

0

Expt

(h)

Predicted 0.5

0 0


0


Figure 9. Prediction of radial distribution of gas void fraction (a), IAC (b), axial gas (c) and liquid (d) velocity for experimental data of Andreu case (VL05JG005TA) and gas void fraction (e), IAC (f), axial gas (g) and axial liquid (h) velocity for Andrew case (VL05JG010TA).


1.6 1.2

Expt

0.8

Predicted

(a)

0.4 0 0


1500

(b)

Expt

IAC (1/m)

1000

Predicted 500

0 0


Gas velocity (m/s)

3

(c) 2

Expt 1

Predicted 0 0 3


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60



Page 49 of 58


(d)

2

Predicted

1

0 0


Figure 10. Prediction of radial distribution of gas void fraction (a), IAC (b), and axial gas velocity (c) and liquid velocity (d) by new proposed correlation for experimental data of Andreu case (VL10JG030TA). 49 ACS Paragon Plus Environment


0.1 (a)

Expt

0.9

Predicted 0.6

Expt Predicted

0.05

0.3 0

0 0

0.2


0.8

1

0

0.2


0.8

1

2 Expt

1.2


1.6 Liquid velocity (m/s)

(d)



1.2

(b)

Predicted

0.8 0.4 0

(e)

1

Expt Predicted

0 0

0.2


0.8

1

0

0.2

2

1.2


0.8

1 (f)

Gas velocity (m/s)

Gas velocity (m/s)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 50 of 58

0.9 0.6

Expt

0.3

(c)

Predicted

0

1 Predicted

0 0

0.2


0.8

1

0

0.2


0.8

Figure 11. Prediction of radial distribution of gas void fraction (a), axial liquid(b) and gas velocity (c) for JL-0.376, JG-0.347 and radial distribution of gas void fraction (d), axial liquid (e) and gas (f) velocity for JL-1.087, JG-0.027 for experimental data of Liu and Bankoff.


1

Page 51 of 58

0.2



0.2 (a) 0.15

Expt

Predicted

0.1 0.05

Predicted 0.1 0.05 0

0

0.2

0.4 0.6 0.8 Radial position (r/R)

0

1

0.2


0.8

1

2

2 (b) 1.5 1

Expt 0.5

Predicted

(e)

Expt



(d)

Expt

0.15

0

1.5

Predicted

1 0.5

0

0 0

0.2


0.8

1

0

0.8

1 (f)

Predicted

k (m2/s2)

Expt Predicted

0.02


Expt

(c) 0.03

0.2

0.02

0.04

k (m2/s2)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


0.01

0.01 0

0 0

0.2


0.8

1

0

0.2

0.4 0.6 0.8 Radial positions (r/R)

Figure 12. Prediction for radial distribution of gas void fraction (a), liquid velocity (b) and turbulent kinetic energy (c) for JL-1 and JG-0.02 and radial distribution of gas void fraction (d), liquid velocity (e) and turbulent kinetic energy (f) for JL-0.5 and JG-0.018 for the experimental data of Hosokawa and Tomiyama.


1


(a)

(b)

16000

(c)

Predicted Pressure drop (Pa/m)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 52 of 58

(d)

+ 15%

12000 - 15% 8000

4000

0 0

4000 8000 12000 Measured Pressure drop (Pa/m)

Figure 13. Gas-liquid pressure drop as a function of gas superficial velocity for (a) pipe ID – 50 mm, (b) pipe ID 38 mm, (c) pipe ID 24.5 mm, and (d) comparison of total pressure drop between measured and predicted values in vertical gas-liquid two phase flow.


16000

Page 53 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47


Table 1. Summery of previous literature on simulation flow conditions for bubbly, slug, churn and annular flow regimes. Author/year

Geometry (mm)

Flow Patterns

Superficial Velocities (m/s)

Volume fraction (%)

Reported bubble diameter (mm)

ASM.Atiqul Islam (2016)

D = 51.2

Bubbly

JL = 1.0170 - 1.6110 JG = 0.0574 - 0.2190

α = 4.5 - 20

dB = 4 - 6

C. PenaMonferrera (2016)

D = 52 L = 5500 D = 50.8 L = 3000 D = 25 L = 2000 D = 38

Bubbly

JL = 0.491 -2 JG = 0.036 0.2

α = 1.46 - 9.2

dB = 2.26 3.81

Open FOAM

dB = 2.5

NEPTU NE code

dB = 3 - 5

Open FOAM

O. Marfaing (2016)

Bubbly

Sebastian Kriebitzsch (2016)

D = 25 200

Bubbly

Zhi Shang (2015)

D = 5 0.8 L = 3810 D = 51.2 L = 4000 D = 50 L = 6000

Bubbly, Slug, Churn, Annular, Mist

Stephen Yamoah (2015)

D = 52 L = 5500

Bubbly

JL = 0.45 1 JG = 0.015 0.2 JL = 0.04 3.458 JG = 0.615 - 12.2

α = 1.7 - 15.7

JL = 0.5 - 1 JG = 0.05 0.3

α=522

CFD Code

α = 1090

dB = 3.02 6.5

CFX

Approach

Turbulence model

Drag force coefficient

Lift force coefficient

Wall Lubrication force coefficient

Eulerian

Low Re k-ε

Moonahan and fox

Tomiyama

Tomiyama

Eulerian

k-ε

SchillerNaumann, Roghair, MagnaudetLegendre

Wang

Antal

0.1

k-ε

0.1

0.03

Antal

0.003

Eulerian

k-⍵ SST

Ishi Zuber

Tomiyama

Hosokawa

Burns

Eulerian

k-ε

Shang

Tomiyama

Antal

Burns

Eulerian

k-⍵

Tomiyama, Grace, Simonnin, Ishi-Zuber

Tomiyama, Saffman, Magnaudet

Antal Frank

Lopez. Burns


Turbulent Dispersion force dispersion coefficient

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47

Roland Rzehak (2015)

D = 51.2 D = 57.2

Montoya Gustavo (2014)

D = 51.2 D = 195.3

Roland Rzehak (2013 (a))

D = 51.2 L = 3030 D = 57.2

Roland Rzehak (2013 (b))

JL = 0.4050 - 1.0167 JG = 0.0111 - 0.22

α = 115.7

Churn turbulent, Bubbly to churn turbulent Bubbly

JL = 1.017 JG = 0.219 0.342

α = 20 80

JL = 0.4050 - 1.0167 JG = 0.0096 - 0.2

α= 1.89 15.7

D = 57.2

Bubbly

JL = 0.5 - 1 JG = 0.12 0.22

Roland Rzehak (2013 (c))

D = 25, 57.2, 200

Bubbly and transition to slug

JL = 0.45 1 JG = 0 0.14

X.Y. Duan (2011)

Bubbly

JL = 0 - 4 JG = 0- 14

F. Pellacani (2011)

D = 51.2 L=3500 D= 195 L= 9000 D = 52 L = 3340

Bubbly, Bubbly to slug

JL = 0.5 - 2

Xia Wang (2011)

D = 48.3 L=3

Bubbly

JL = 0.064 - 2.34 JG = 0.043 0.506

dB = 2.94 4.5

Page 54 of 58

Open FOAM, Ansys CFX

Eulerian

k-⍵ SST

Ishi and Zuber

Tomiyama

Tomiyama/ Hosokawa

Burns

CFX

Eulerian

SST

Ishi and Zuber

Tomiyama

Hosokawa

Burns

dB = 2.94 4.5

Ansys CFX

Eulerian

k-⍵ SST

Ishi and Zuber

Tomiyama

Tomiyama/ Hosokawa

Burns

α = 9.6 - 15.7

dB = 2.94 4.22

Ansys CFX

Eulerian

k-⍵ SST

Ishi and Zuber

Tomiyama

Antal

Burns

α = 1.7 - 15.7

dB = 2.94 6.6

Ansys CFX, CATHA RE

Eulerian

k-⍵ SST

Ishi and Zuber

Tomiyama

Antal

Burns, Carrica, Lopez

dB = 0 - 60

CFX

k-⍵ SST

Zuber

Tomiyama

Antal

Burns

α=515

IATE approach

CFX

Eulerian

k-⍵ SST

Grace

Tomiyama

Antal

FAD

α= 11.7 12. 8

dB = 2.4 4.5

Fluent

Eulerian

Tomiyama

0, 0.1 and Tomiyama

MUSIC model


Page 55 of 58 1 2 3 Xia Wang 4 (2010) 5 6 7 8 Present work 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47


D = 48.3

Bubbly

JL = 0.064 - 5.1 JG = 0.039 1.275

α= 2.2125.7

dB = 1.8 9.3

Fluent

Eulerian

k-ε

Tomiyama

Tomiyama

Antal

D= 24.5 57.2 L = 29228000

Bubbly

JL = 0.376 – 3.53 JG = 0.0036 - 1.275

α = 0.338.74

Predicted (dB)

Fluent

Eulerian

k-⍵ SST

Grace

Tomiyama

Hosokawa


Burns

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47

Page 56 of 58

Table 2. Flow conditions for the wall peak, flattened wall peak and core peak case. Letters "MT" denotes the MT Loop test of (Prasser et al., 2003; Lucas et al., 2005) while letter "L" denotes the test of Liu (1998). Test cases

Wall peak

Flattened wall peak Core peak

Hibiki Sun L21B MT063 Lucas L11A MT041 Fu Run1

Pipe diameter (mm) 50.8 50.8 57.2 51.2 51.2 57.2 51.2 48.3

JL (m/s)

JG (m/s)

0.491 0.615 1.000 1.016 0.405 0.500 1.016 2.607

0.055 0.049 0.140 0.031 0.011 0.120 0.011 1.275

Reported αG (%) dB (mm) 2.50 1.90 3.03 4.50 4.50 2.94 4.50 9.30

10.00 6.04 10.60 2.64 1.89 15.20 1.00 25.70

Reb

Eo

29.25031 333.5765 686.2777 658.8234 879.6546 658.2015 371.7109 15130.62

0.847999 0.489804 1.245663 2.747516 2.747516 1.172762 2.747516 11.73495

Table 3. Different grids used in the simulation of two phase pipe flow. Grids

Number of cells

ncross section

nz

Grid 1

480000

1200

400

Grid 2

800000

1600

500

Grid 3

1200000

2000

600


Measurement location (L/D) 60 75 60 60 66 60 60 60

Page 57 of 58 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47


Table 4. Flow conditions for the validation of correlations for wall peak and core peak cases. Test cases

Andreu Liu and Bankoff Hosokawa and Tomiyama

Pipe diameter (mm) 52.00 38.00 25.00

JL (m/s)

JG (m/s)

Predicted dB (mm)

αG (%)

Reb

Eo

0.50 0.50 1.00 0.376 1.087 1.00 0.50

0.05 0.10 0.30 0.347 0.027 0.02 0.018

3.59 4.62 6.84 4.80 4.10 3.70 3.00

6.00 13.00 22.00 30.09 2.42 1.46 2.31

1212.186 1006.717 625.1324 3308.907 8.013275 1471.609 898.6364

1.748655 2.896004 6.347861 3.126062 2.280777 1.857456 1.221118


Measurement location (L/D) 98 36 68


TOC Graphic

Radial distribution of αg for wall, flattened wall and core peak case Gas volume fraction (-)

0.15 Hibiki et al.

(a)

(b)

0.1

0.05

Comparison of measured and predicted pressure drop

0 0.2

0.4 0.6 0.8 Radial Position r/R (-)

1

0.06

(d) Lucas et al.

16000

(c) + 15%

0.04 Predicted Pressure drop (Pa/m)


0

12000

0.02

0 0 Gas volume fraction (-)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 58 of 58

0.2

0.4 0.6 Radial Position r/R (-)

0.8

1

0.6

(f) Fu et al.

- 15% 8000

4000

0.4

(e) 0.2

0 0

0 0

0.2

0.4 0.6 0.8 Radial Position r/R (-)

1


4000

8000 12000 Measured Pressure drop (Pa/m)

16000

Computational Fluid Dynamics Simulations of the AirâWater Two

Recommend Documents

Computational Fluid Dynamics Simulations of the AirâWater Two