Next Article in Journal
Modified Transboundary Water Interaction Nexus (TWINS): Xayaburi Dam Case Study
Previous Article in Journal
1-D Vertical Flux Dynamics in a Low-Gradient Stream: An Assessment of Stage as a Control of Vertical Hyporheic Exchange
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Modeling of the Free-Surface Vortex-Driven Bubble Entrainment into Water

1
Helmholtz-Zentrum Dresden-Rossendorf (HZDR), Institute of Fluid Dynamics, Bautzner Landstr. 400, 01328 Dresden, Germany
2
Department of Mechanical & Industrial Engineering, Faculty of Engineering, Universitas Gadjah Mada, Jalan Grafika No. 2, Yogyakarta 55281, Indonesia
*
Author to whom correspondence should be addressed.
Water 2020, 12(3), 709; https://doi.org/10.3390/w12030709
Submission received: 20 December 2019 / Revised: 13 February 2020 / Accepted: 3 March 2020 / Published: 5 March 2020
(This article belongs to the Section Hydraulics and Hydrodynamics)

Abstract

:
The recently developed GENTOP (Generalized Two Phase Flow) concept, which is based on the multifield Euler‒Euler approach, was applied to model a free-surface vortex—a flow situation that is relevant for hydraulic intake. A new bubble entrainment model has been developed and implemented in the concept. In general, satisfactory agreement with the experimental data can be achieved. However, the gas entrainment can be significantly affected by several parameters or models used in the CFD (Computational Fluid Dynamics) simulation. The scale of curvature correction C s c a l e in the turbulence model, the coefficient in the entrainment model C e n t , and the assigned bubble size to be entrained have a significant influence on the gas entrainment rate. The gas entrainment increases with higher C s c a l e values, which can be attributed to the stronger rotation captured by the simulation. A smaller bubble size gives higher gas entrainment, while a larger bubble size leads to a smaller entrainment. The results also show that the gas entrainment can be controlled by adjusting the entrainment coefficient C e n t . Based on the modeling framework presented in this paper, further improvement of the physical modeling of the entrainment process should be done.

1. Introduction

A free-surface vortex (see Figure 1) may exist in a wide range of scales; it can be as small as a “bathtub vortex” [1,2,3] or can be as big as an ocean whirlpool [4]. The topic of a free-surface vortex is often found in the discussion of hydropower plants, nuclear reactors and other applications using pumps. The supply of water for irrigation, domestic, industry, and power generation is usually taken from rivers or reservoirs through an intake that is located near the surface [5]. Insufficient submergence (a short distance between the water surface and an intake) may lead to the formation of a free-surface vortex that can induce gas entrainment into the intake [5].
A free-surface vortex and its associated gas entrainment may lead to several operational and safety problems [5,6,7,8,9]. They may cause mechanical damage and loss of performance in fluid machinery such as turbines and pumps [6,7]. A swirl in a sump leads to rotational flow in a pipe, which may reduce the performance of the pump [8,9]. If such flow is unsteady, it may also cause fluctuating loads on pump bearings [8]. The gas entrainment induced by a free-surface vortex will reduce the delivery of a pump (1% air reduces the efficiency of a centrifugal pump by 5–15%) [6]. This reduction may cause a severe problem such as the overtopping of a dam [3], which may lead to a safety hazard and cause loss of life [7].
Vortex-induced gas entrainment is also an important issue for more specific applications such as nuclear reactors. In sodium-cooled fast reactors (SFRs), an inert cover gas such as argon is used and maintained above the sodium surface to accommodate the volume changes of sodium and prevent the contact of sodium with air [11,12]. Several mechanisms could lead this cover gas to entrain the sodium region, including the entrainment driven by a free-surface vortex [12]. The gas entrainment may cause safety and operational problems and, for this reason, becomes an essential issue in SFR safety analyses [13,14,15,16]. This issue has been intensively investigated by the Japan Atomic Energy Agency (JAEA) and its institutional partners [16,17,18,19,20,21,22,23,24,25,26,27,28,29,30]. The problems associated with gas entrainments in SFRs are: changes in reactivity when the gas reaches the core [11,14,18,31,32,33], burnout of the fuel pin due to the trap of large bubbles [31], thermal stresses in the reactor structure [32,33], pump cavitation and fluctuations in pump discharge [14,33], reduction in the heat transfer efficiency [11,33], disturbance in the prompt detection of fission products leakage from failed fuel pins [34], disturbance of electromagnetic sensors used for shutdown systems [34], and trouble with acoustic or ultrasonic instrumentation such as boiling noise detector [31,32,34]. The gas entrainment issue is not only applicable to SFRs, but also to Pressurized Water Reactors (PWRs) and Boiling Water Reactors (BWRs). During mid-loop operation, gas may entrain into the Reactor Coolant System (RCS) of PWRs due to a free surface vortex and then be sucked into the Decay Heat Removal System (DHRS), which finally leads to a disturbance of the instrumentation [35]. A total failure of DHRS may occur if the void fraction exceeds 15% [35]. In the BWRs, the gas entrainment due to a free surface vortex may occur at the suction inlet from a condensation chamber/wet-well [36,37].
Computational fluid dynamics (CFD) may help to design a safer process, minimizing the aforementioned risks associated with gas entrainment due to a free-surface vortex. Generally, the previous CFD works available in the literature can be divided into two parts: single-phase and two-phase computations. In a single-phase simulation, the deformation of the free surface is not considered and the free surface is defined as a free slip boundary [38,39]. When the mesh resolution is sufficient and the appropriate turbulence model is used, the velocity fields can be calculated using single-phase simulation, as reported by [17,40]. However, to the best of our knowledge, the estimation of the gas entrainment rate has never been performed by a single-phase CFD. In addition, the direct observation of a free-surface vortex from the single-phase simulation is not possible. A post-processing method is required to judge the occurrence and the location of the vortex, e.g., Q criterion [38]:
Q = 1 2 ( Ω 2 S 2 ) > 0 ,
where Q is the second invariant of velocity gradient tensor, Ω is the vorticity tensor, and S is the strain rate tensor. The above equation states that a free surface vortex exists when the strength of rotation is bigger than the local strain rate [38].
In the case of two-phase simulation, usually, the volume of fluid (VOF) model is employed [20,40,41]. The computation is performed in a fixed grid solving only one momentum equation, which is shared by both fluids [42]. Generally, a very fine mesh is required to resolve the interface, e.g., [40]. In the case of bubble entrainment driven by the free-surface vortex, the entrainment process needs to be resolved. This requires a very fine mesh, especially at the tip of the gas core. To reduce the number of computational cells used in the simulation, one may refine the mesh only in the tip region. However, to apply the local refinement, the precise location of the tip should be known before the simulation. This will be difficult to realize in practice since the tip location is transient in nature [43]. Therefore, a fine mesh should be used in a larger predicted region or over the whole computational domain instead of only at the tip region, which leads to very high computational costs. When there is a large computational domain such as a hydropower plant, this becomes a significant disadvantage of the method. An option to overcome this problem is using the Euler‒Euler model, wherein the gas and liquid have their own sets of momentum equations. The possible advantage is that a coarser mesh may be used to model this gas entrainment. However, appropriate closure models need to be used, including the gas entrainment model.
The GENTOP (Generalized Two Phase) concept, based on the Euler‒Euler model, which has been developed at Helmholtz-Zentrum Dresden-Rossendorf (HZDR), aims to handle various flow conditions where multiscale interfacial structures exist [44]. The potency of this concept has been demonstrated in several types of flows such as impinging jet [44], bubble column [44], a dam-break case with obstacle [45], churn turbulent flow in the vertical pipe [46], and boiling case [47,48]. However, to increase the effectiveness, continuous improvement and validation of various flow conditions are required. For that purpose, the applicability of the GENTOP concept to the flow case of bubble entrainment driven by a free-surface vortex was investigated in this work. A new bubble entrainment model considering the physics of the flow has been developed and implemented in the GENTOP concept in this study.

2. Model Description

2.1. The GENTOP Concept

The Euler‒Euler model was used in this work, whereby the following continuity and momentum equations were solved across the fixed computational mesh [49]:
t ( α j ρ j ) + . ( α j ρ j u j ) = S j
t ( α j ρ j u j ) + . ( α j ( ρ α u j × u j ) ) = α j p + . ( α j μ j ( u j + ( u j ) T ) ) + M j + S M j
M j = F d r a g + F l i f t + F w a l l + F T D + F V M .
In the above equations, the liquid and gas phases are assumed to be adiabatic and incompressible. The subscript j represents the phase, while α j , ρ j , u j , t, S j , and p are the volume fraction, density, velocity vector, time, mass source, and pressure, respectively.   M j and S M j denote the interfacial forces and the momentum sources (i.e., due to external body forces), respectively [49].
The flow is represented by a continuous liquid phase, one velocity group for the polydispersed gas phase (dg) and one for the potentially continuous gas (cg) in the frame of the GENTOP concept [44]. The closure models are assigned specifically for each of the gas fields (i.e., dg and cg). For the dispersed gas phase, dg, the closure models based on the baseline model concept [50] are used as given in Table 1. These closure models influence the simulation results obtained in the complex rotating gas‒liquid flow case, as discussed, e.g., in [51,52]. To avoid numerical instability problems, the virtual mass force (Equation (17)) was not used in this study.
The detection of the statistically resolved large interface between gas and liquid is done in the GENTOP concept by the free surface function φ f s , which is formulated as follows [44]:
φ f s = 0.5   t a n h ( a f s Δ x ( | α c g | | α c g | c r i t ) ) + 0.5 ,
where | α c g | is the volume fraction gradient of the potentially continuous gas and | α c g | c r i t is its critical value, which is given by [44]:
| α c g | c r i t = 1 n Δ x .
The blending coefficient a f s and n were set to 100 and 4, respectively [44].
To support the transition from dg to cg by the agglomeration of cg when it reaches a critical amount, a so-called clustering force is used in the GENTOP concept [44]:
M L c l u s t = M c g c l u s t = c c l u s t ( 1 φ f s ) φ c l u s t ρ L α L ,
where c c l u s t is the coefficient that defines the strength of the force [44]. The coefficient was set to c c l u s t = 0.1, similar to the value used by [44] for the impinging jet case. The blending function φ c l u s t is used to define the location where the force should be activated [44]:
φ c l u s t = ( 0.5   t a n h [ a B ( α c g α c l u s t , m i n ) ] + 0.5 ) . ( 0.5 × t a n h [ a B ( α c l u s t , m a x α c g ) ] + 0.5 ) ,
where a B is the blending coefficient, set to 20 [44]. The force is blended in when α c g is equal to 0.2 (i.e., α c l u s t , m i n = 0.2) and blended out when it exceeds 0.8 (i.e., α c l u s t , m a x = 0.8) [44].
To account for the surface tension, the following surface tension model of [64] that was implemented by [46] in the GENTOP concept was used:
M S T n = σ . ( α c g | α c g | ) α c g .
The interfacial area density and drag coefficient for cg in the GENTOP concept are calculated, respectively, as [44]:
A D , c g = ( 1 φ m o r p h ) A D , b u b b + φ m o r p h A D , f s
C D , c g = ( 1 φ m o r p h ) C D , b u b b + φ m o r p h C D , f s ,
where φ m o r p h is a blending function, formulated as follows [44]:
φ m o r p h = 0.5   t a n h ( a B ( α c g α c g , c r i t ) ) + 0.5 ,
where α c g , c r i t is the critical volume fraction of cg, set to 0.3 [44]. The interfacial area density of bubble A D , b u b b and free surface A D , f s are calculated as follows [44]:
A D , b u b b = 6 α c g d c g
A D , f s = ( 2 | α L | | α c g | ) ( | α L | + | α c g | ) .
The drag coefficient of bubble C D , b u b b is calculated using Equation (6), while that of the droplet and free surface are given by [44]:
C D , f s = m a x ( 0.01 , ( 2 ( α L τ W , L + α c g τ W , c g ) ) ρ L ( u L u c g ) 2 ) .
Since the focus of this present work was to investigate the capability of the newly proposed entrainment model (described in the next section) to model the bubble entrainment driven by a free-surface vortex, bubble break-up and coalescence models as well as the complete coalescence model in the GENTOP concept were deactivated in this study. It is also important to note that liquid droplets were not considered in this work since they were not observed or measured in the experiment.

2.2. The New Entrainment Model

The pressure balance for a stable gas core is presented in Figure 2. The hydraulic pressure in the liquid phase is higher than in the gas core due to the higher density of the liquid. In the stable condition, the vortex gas core in the rotational flow does not collapse due to the contribution of centrifugal force, which counteracts the pressure from the liquid side. The pressure balance can be formulated as follows:
ρ G g Δ z + f c   Δ z = ρ L g Δ z ,
where Δ z is the distance to the liquid level and f c   is the centrifugal force density, which can be calculated based on the forced vortex model as follows:
f c   = ρ L ω 2 r e q 4 .
where ω is the vorticity magnitude and r e q is the radius of the gas core at the equilibrium state. Substituting Equation (31) to Equation (30) yields:
ρ G g Δ z + ρ L ω 2 r e q Δ z 4 = ρ L g Δ z
ρ G g + ρ L ω 2 r e q 4 = ρ L g
r e q = 4 ( ρ L ρ G ) g ρ L ω 2 .
The underlying assumption of this entrainment model is that the detachment of the bubble from the gas core depends on the relation between the radius of the gas core related to the turbulent length scale L t :
L t = k 3 / 2 ε .
The entrainment Q e n t is assumed to be proportional to the ratio between turbulent length scale L t and the gas core radius r e q :
Q e n t ~ L t r e q .
The gas entrainment per unit volume and time then can be formulated as follows:
Q e n t = { 0     L t r e q   Q e n t = C e n t k 3 / 2 ε ω 2 ρ l g ( ρ l ρ g ) ( α l U l )   L t > r e q   .
The gas entrainment model is active only when L t > r e q . C e n t is a dimensionless coefficient and ω is the liquid vorticity magnitude. In the GENTOP concept, the entrainment model was then implemented as a source for the continuity equation as follows:
m ˙ e n t = φ f s . Q e n t   .   α c g . ρ G .
This entrainment model controls the conversion of continuous gas cg into dispersed gas dg in the region detected as “interface” in the GENTOP concept, which is defined by φ f s .
The gas entrainment model formulated in Equation (37) is different from the entrainment model proposed by [65], which has the formulation:
Q g ( x ) = C e n t g e n t k ( x ) u n n ( x ) .
In the proposed entrainment model described in Equation (37), the vorticity magnitude is used to take into account the rotation and to identify the location of the gas core tip, which is not included in the entrainment model of [65].

2.3. The Turbulence Model

The two-equation k‒ω-based shear stress transport (SST) model proposed by [66] and the dispersed zero equation model are usually used for the liquid and gas phase, respectively, in the GENTOP concept [44,45]. However, the eddy viscosity models are insensitive to streamline curvature and system rotation [49] due to their isotropic nature [39]. Therefore, the curvature correction implemented in the SST model by [67] was used for the liquid phase in this work. The correction function f r is used in the model to control the turbulence production term P k as follows [49]:
P k P k . f r
f r = m a x ( 0 , 1 + C s c a l e ( f ˜ r 1 ) )
f ˜ r = m a x { m i n ( f r o t a t i o n   ,     1.25 ) ,   0 } ,
where f r o t a t i o n is the empirical function proposed by [68]:
f r o t a t i o n = ( 1 + c r 1 ) 2 r 1 + r ( 1 c r 3 t a n 1 ( c r 2 r ˜ ) ) c r 1 .
C s c a l e in Equation (41) is the scale of curvature correction available in ANSYS CFX, which can be adjusted [49]. The arguments r and r ˜ in Equation (43) are given by [49]:
r = S Ω
r ˜ = 2 Ω i k S j k ( D S i j D t + ( ε i m n S j n + ε j m n S i n ) Ω m R o t ) 1 Ω D 3 ,
where S i j and Ω i j represent the strain rate and vorticity tensor, respectively [49]:
S i j = 1 2 ( U i x j + U j x i )
Ω i j = 1 2 ( U i x j U j x i ) + 2 ε m j i Ω m r o t ,
where
S 2 = 2 S i j S i j
Ω 2 = 2 Ω i j Ω i j
D 2 = m a x ( S 2 , 0.09 ω 2 ) .
The empirical constants c r 1 , c r 2 , and c r 3 in Equation (43) are set to 1.0, 2.0, and 1.0, respectively [49]. To consider bubble-induced turbulence, additional source terms based on [69,70] were used in the turbulence model.

3. CFD Set-Up

For the present work, an experimental dataset of bubble-type entrainment of the Ikoma experiment [71,72] was selected. The computational domain for this CFD study was adapted from the work of [72] with a slight modification to improve the boundary condition for the outlet. The computational domain, the mesh, and the defined boundaries are presented in Figure 3. The geometry consists of three major parts: the cylindrical test section, the connection pipe, and the bubble catcher or the buffer tank. In the Ikoma experiment, water enters the cylindrical test section through the tangential inlet and leads to the formation of a free-surface vortex [71]. In the case of bubble-type entrainment, bubbles are detached from the gas core and then entrained into the connection pipe [71]. Finally, the entrained bubbles are separated from the water in the buffer tank [71].
To investigate the influence of mesh size on the gas entrainment rate, two different meshes, namely Mesh A and Mesh B, which consist of 131,712 and 408,728 hexahedral cells, respectively, were used in this work. The meshes are refined in the central region and at the bottom of the cylindrical vessel, at the gas‒liquid interface and in the connection pipe (see Figure 3a). The liquid level in the cylindrical test section h l and the exit velocity from the vessel V e are 60 mm and 0.66 m/s, respectively. The definition of this exit velocity is the superficial velocity of water in the connection pipe (see Figure 3b). The water flow rate Q can be calculated as follows:
Q = V e π D 2 4 ,
where D is the diameter of the outlet pipe (D = 8 mm). From the calculation, it is known that the water flow rate is 2 L/min. The inlet velocity V i n can be determined by the following equation:
V i n = Q L i n h l ,
where L i n is the width of the inlet slit ( L i n = 10 mm). This results in V i n being equal to 0.055 m/s. V o u t can be calculated in a similar manner, resulting in a velocity of 0.33 m/s.
The boundaries that were applied to the computational domain are shown in Figure 3b. A uniform liquid velocity of 0.055 m/s was defined at the inlet. An opening boundary condition was applied at the top of the cylindrical vessel, while a degassing boundary was used at the top of the buffer tank. Uniform velocity boundary condition of V o u t = 0.33 m/s was applied at the outlet. All other boundaries were defined as no-slip wall and free-slip wall for the liquid phase and gas phase, respectively. Resolving the viscous sublayer is not required in the simulations since a wall function assuming a smooth wall is used [73].
In this study four bubble classes, G1‒G4, having a diameter from 1 mm to 7 mm, are defined within the dispersed gas velocity group dg, while the last group G5, having a diameter of 9 mm or larger, is classified as continuous velocity group cg (see Table 2). The investigations were performed in several steps. First, the investigation of the influence of the scale of curvature correction was carried out by performing simulations using various C s c a l e . For this investigation, a fixed entrainment coefficient C e n t = 2 × 10−4 was used to control the mass transfer from cg to G1.
Next, the influence of the entrained bubble size distribution was investigated by performing three separate simulations, namely Sim. A, Sim. B, and Sim. C (see Table 2). In Sim. A, only the smallest bubble class G1 is entrained from the potentially continuous gas. In Sim. B, G1 and G2 with equal entrainment fractions are assigned to be entrained. In Sim. C, all bubble classes are assigned to be entrained with uniform entrainment fraction. The scale of curvature correction C s c a l e = 1 and the entrainment coefficient C e n t = 2 × 10−4 were used in this step. The investigation of the influence of the entrainment coefficient was also carried out by performing simulations with two different C e n t . The mass transfer from cg to G1 with the scale of curvature correction C s c a l e = 1 was set in the simulations. Finally, the influence of the mesh size was also investigated and presented in the last section.

4. Results and Discussion

4.1. The Influence of the Turbulence Model

The interfacial deformation and bubble entrainment using various C s c a l e with a fixed value of C e n t = 2 × 10−4 are presented in Figure 4a,b, respectively. In the case of C s c a l e = 0, the free surface is almost flat. At this condition, almost no bubble entrainment can be observed. Increasing C s c a l e to 0.5 leads to dimple formation. A spot of the dispersed gas fraction is observed, indicating that the entrainment model is activated. However, no bubble entrainment (with α d g equal to or more than 0.02) into the connection pipe can be observed. Increasing C s c a l e to 1 leads to the formation of the gas core. At this condition, the bubble entrainment into the connection pipe is observed, confirming that the entrainment model successfully takes into account the bubble entrainment phenomenon.
To have a quantitative description regarding the influence of C s c a l e on the bubble entrainment, the transient plot of continuous, bubble, and total gas entrainment are shown in Figure 5a‒c, respectively. Those plots confirm that almost no gas entrainment is produced other than C s c a l e = 1. The averaged gas entrainment rate shown in Figure 5d shows that the gas entrainment using C s c a l e = 1 is mostly obtained from the entrainment of the dispersed gas phase. This is the expected condition since the bubble type entrainment is observed in the experiment. Figure 5e shows that the gas entrainment rate using C s c a l e = 1, which is the default value in ANSYS CFX, is the most suitable value since it leads to the entrainment rate closest to the experimental value.
To reveal the influence of C s c a l e on velocity fields, the velocity vectors represent the tangential projection of the liquid velocity are evaluated on horizontal planes H1, H2, and H3, which are located at an axial distance of 10, 20, and 30 mm from the bottom surface of the cylindrical vessel, respectively. These planes are in the liquid region to avoid the distraction of the velocity fields due to the presence of the gas core. Figure 6 presents the velocity vectors and the liquid vorticity magnitude as the background color. The velocity vector confirms that for all cases the flow is rotated and the center of rotation is not located precisely at the center of the cylindrical vessel. The magnitude of the rotation represented by the vorticity is largely affected by the selected C s c a l e used in the simulation. Increasing C s c a l e leads to a higher vorticity value that contributes to a higher gas entrainment. This can be understood since, according to Equation (37), the gas entrainment rate is proportional to the square of the vorticity magnitude. It can also be observed from the figure that the vorticity magnitude increases as the distance from the bottom of the vessel or the intake pipe decreases.

4.2. The Influence of Entrained Bubble Size Distribution

The contour of dispersed gas fractions at t = 6 s plotted on the central vertical plane is presented in Figure 7 for the different simulations listed in Table 2. The results of Sim. A clearly show bubbles being entrained into the connection pipe and then accumulating in the buffer tank (see Figure 7a). Although bubble entrainment can also be observed in Sim. B, it was not as significant as in Sim. A (see Figure 7b). Some of the dispersed gas still collected in the top of the gas core. In Sim. C, the dispersed gas mostly accumulated near the free surface close to the gas core instead of moving downward to the intake pipe (see Figure 7c). This may be because the terminal velocity for larger bubbles is bigger and the downward liquid velocity is not high enough to force bubbles to move into the connection pipe. The fact that some of the bubbles produced in Sim. B and Sim. C did not merge with the potentially continuous gas above the free surface is due to the inactivation of the complete coalescence model. This also shows that the complete coalescence is an important feature in GENTOP to bring the modeling of this flow behavior closer to the real condition, as we know from the literature (e.g., [13]) that there is a condition in which some of the entrained bubbles return to the free surface instead of flowing to the outlet. This inactivation in the simulations elucidates the role of the entrainment model alone.
The transient profile of the potentially continuous gas cg, the dispersed gas dg, and the total gas entrainment for all simulations are presented in Figure 8a‒c, respectively. The transient profile of cg entrainment of all simulations shows that cg entrainment is negligible. The transient profile of dg entrainment shows a relatively stable curve after 4 s of simulation in the case of Sim. A and Sim. B. The dg entrainment in Sim. A is always higher than in Sim. B. In contrast, dg entrainment in Sim. C is always negligible. Figure 8d shows the averaged entrainment of cg and dg. This plot confirms that the entrainment of the dispersed gas phase contributes most of the gas entrainment. The total gas entrainment obtained in Sim. A is around 2.6 times that obtained in Sim. B. Sim. A gives the best fit to the experimental value (see Figure 8e). The results show that the entrainment rate that is controlled by the entrainment model is highly influenced by the entrained bubble size distribution. In the present study, the smallest bubble class, G1, is suitable to obtain a higher entrainment.

4.3. The Influence of the Entrainment Coefficient

The interfacial deformation and bubble entrainment using two different entrainment coefficients C e n t are presented in Figure 9a,b, respectively. The isosurface of the potentially continuous gas shows no significant difference in the interfacial deformation between those two cases. In contrast, larger bubble entrainment is observed in the case of C e n t = 2 × 10−4.
To obtain a quantitative description, the transient profile of potentially continuous gas cg, dispersed gas dg, and total gas entrainment for simulations using two different values of C e n t are presented in Figure 10a‒c, respectively. The plot of cg entrainment shows almost no entrainment obtained from cg for both simulations. In contrast, dg entrainment shows a significant difference of gas entrainment behavior, where the entrainment obtained for C e n t = 2 × 10−4 is around 1.7 times that obtained for C e n t = 1 × 10−4. This shows that this coefficient can control the entrainment rate. The average cg and dg entrainment values show that dg entrainment is much higher than cg entrainment for both simulations (see Figure 10d). The result of the simulation using C e n t = 2 × 10−4, as shown in Figure 10d, is the same as in Figure 8e, so is in the best agreement with the experiment.
The velocity vectors represent the tangential projection of the liquid velocity on horizontal planes H1, H2, and H3, as presented in Figure 11. The background color represents the liquid vorticity magnitude. The influence of C e n t on the magnitude of the rotation represented by the vorticity magnitude seems to be minor. For all the evaluation planes, no significant difference can be observed from the simulation when using two different C e n t .

4.4. The Influence of the Computational Cell Size

The interfacial deformation for different mesh sizes is presented in Figure 12. The deformation and gas core length are significantly influenced by the cell size. For the same C s c a l e (i.e., C s c a l e = 1), the gas core length in the simulation using the finer mesh (i.e., Mesh B) is about 70% larger than the one using Mesh A (see Table 3). The differences in the deformation and gas core length affect the gas entrainment rate. The gas entrainment rate for Mesh B is three times that of Mesh A. As discussed earlier, the coefficient of the curvature correction, C s c a l e has an important influence on the shape of the free-surface vortex. For this reason, this parameter is modified. A similar vortex shape and gas core length are obtained for Mesh A with C s c a l e = 1 and Mesh B with C s c a l e = 0.85. As shown in Table 3, a similar entrainment rate is obtained for this case. This means that the entrainment model itself has the ability to be mesh-independent. The observed mesh dependency mainly results from the different curvature correction required for different mesh sizes.
Table 3 also shows the large dependence of the calculated vorticity magnitude on the mesh size. The dependence is caused by the inability of the curvature correction method to predict the same value of the vorticity magnitude for different cell sizes. As a consequence, the gas entrainment also differs largely when a different cell size is used. The gas entrainment model is expected to provide a mesh-independent solution when the turbulence model is able to estimate the mesh-independent magnitude vorticity.

5. Conclusions

Turbulence modeling is important for capturing the development of a free-surface vortex and gas entrainment. In the CFD simulation using the SST turbulence model, curvature correction is required to increase the sensitivity on the curvature or rotation. A lower value of C s c a l e leads to only a shallow dimple, while higher C s c a l e provides a better gas core formation. The bubble entrainment increases with the increase in C s c a l e , which can be attributed to the stronger rotation captured by the simulation. The results obtained by using a higher C s c a l e are characterized by the higher vorticity magnitude. The increase of vorticity magnitude closer to the connection pipe is observed for all CFD simulations presented in this paper.
A new entrainment model considering the physics of the phenomena was implemented in the GENTOP concept. In general, a bubble entrainment rate of the same order as in the experiment can be achieved. In the case presented in this work, the CFD simulation where the entrainment model was used to control the entrainment model from the potentially continuous gas into the smallest bubble class group gave the highest gas entrainment and the best fit with the experiment. Using a higher diameter of entrained bubbles led to a smaller gas entrainment. The entrainment coefficient can be adjusted to control the entrainment rate without significantly altering the local velocity fields or the vorticity magnitude. A higher coefficient leads to higher gas entrainment. The gas entrainment rate also significantly affects the gas entrainment rate, which mainly results from the different curvature correction required for different mesh sizes.

Author Contributions

Conceptualization, R.A.P. and D.L.; methodology, R.A.P. and D.L.; simulation, R.A.P.; investigation, R.A.P. and D.L.; writing—original draft preparation, R.A.P.; writing—review and editing, R.A.P. and D.L.; supervision, D.L. All authors have read and agreed to the published version of the manuscript.

Funding

R.A.P. received a scholarship from the Indonesian‒German Scholarship Programme (IGSP).

Acknowledgments

The authors would like to express their gratitude to Dr. Frank Blömeling from TÜV NORD SysTec GmbH & Co. KG; Dr. Kei Ito from Japan Atomic Energy Agency (JAEA); and Aljaž Skerlavaj from Kolektor Turboinstitut Slovenia for valuable discussions.

Conflicts of Interest

The authors declare no conflict of interest.

Nomenclature

Latin symbols
A D interfacial area density (m−1)
C D drag coefficient (dimensionless)
C e n t entrainment coefficient (m s−1)
C L lift coefficient (dimensionless)
C s c a l e scaling coefficient of curvature correction (dimensionless)
C V M virtual mass coefficient (dimensionless)
C W wall force coefficient (dimensionless)
C μ shear-induced turbulence coefficient (dimensionless)
Dpipe diameter (m)
d B bubble diameter (m)
d maximum horizontal dimension of a bubble (m)
E o Eötvös number (dimensionless)
E o modified Eötvös number (dimensionless)
F force (N m−3)
g gravitational acceleration (m2 s−2)
M momentum transfer term (kg m−2 s−2)
p pressure (Pa)
R e Reynolds number (dimensionless)
S M momentum source due to external body forces (kg m−2 s−2)
W mean velocity magnitude (m s−1)
u velocity vector (m s−1)
t time (s)
ydistance to the wall (m)
z axial distance (m)
Greek symbols
α volume fraction (dimensionless)
Δ x characteristic cell length scale (m)
μ dynamic viscosity (Pa s)
ρ density (kg m−3)
σ surface tension (N m−1)
σ T D turbulent Schmidt number (dimensionless)
φ blending function (dimensionless)
Subscripts and superscripts
f s free surface
j phase index
G gas
L liquid
V M virtual mass
W wall

References

  1. Shapiro, A.H. Bath-Tub Vortex. Nature 1962, 196, 1080–1081. [Google Scholar] [CrossRef]
  2. Forbes, L.K.; Hocking, G.C. The bath-plug vortex. J. Fluid Mech. 1995, 284, 43–62. [Google Scholar] [CrossRef]
  3. Andersen, A.; Bohr, T.; Stenum, B.; Rasmussen, J.J.; Lautrup, B. The bathtub vortex in a rotating container. J. Fluid Mech. 2006, 556, 121–146. [Google Scholar] [CrossRef] [Green Version]
  4. Gjevik, B.; Moe, H.; Ommundsen, A. Sources of the Maelstrom. Nature 1997, 388, 837–838. [Google Scholar] [CrossRef]
  5. Jain, A.K.; Raju, K.G.R.; Ramachandra, G.J. Vortex formation at vertical pipe intakes. Asce J. Hydraul. Div. 1978, 104, 1429–1445. [Google Scholar]
  6. Odgaard, A.J. Free-surface air core vortex. J. Hydraul. Eng. 1986, 112, 610–620. [Google Scholar] [CrossRef]
  7. Rindels, A.J.; Gulliver, J.S. An experimental study of critical submergence to avoid free-surface vortices at vertical intakes. In Project Report No. 224 University of Minnesota St. Anthony Falls Hydraulic Laboratory, Minneapolis; University of Minnesota: Minneapolis, MN, USA, 1983. [Google Scholar]
  8. Denny, D.F. An experimental study of air-entraining vortices in pump sumps. Proc. Inst. Mech. Eng. 1956, 106–125. [Google Scholar] [CrossRef]
  9. Chuang, W.L.; Hsiao, S.C.; Hwang, K.S. Numerical and experimental study of pump sump flows. Math. Probl. Eng. 2014, 2014. [Google Scholar] [CrossRef]
  10. Yang, J.; Andreasson, P.; Högström, C.-M.; Teng, P. The Tale of an Intake Vortex and Its Mitigation Countermeasure: A Case Study from Akkats Hydropower Station. Water 2018, 10, 881. [Google Scholar] [CrossRef] [Green Version]
  11. Banerjee, I.; Sundararajan, T.; Sangras, R.; Velusamy, K.; Padmakumar, G.; Rajan, K.K. Development of gas entrainment mitigation devices for PFBR hot pool. Nucl. Eng. Des. 2013, 258, 258–265. [Google Scholar] [CrossRef]
  12. Satpathy, K.; Velusamy, K.; Patnaik, B.S.V.; Chellapandi, P. Numerical simulation of liquid fall induced gas entrainment and its mitigation. International Journal of Heat and Mass Transfer 2013, 60, 392–405. [Google Scholar] [CrossRef]
  13. Baum, M.R. Gas entrainment at the free surface of a liquid: Entrainment inception at a laminar vortex. J. Br. Nucl. Energy Soc. 1974, 13, 203–209. [Google Scholar]
  14. Tenchine, D. Some thermal hydraulic challenges in sodium cooled fast reactors. Nucl. Eng. Des. 2010, 240, 1195–1217. [Google Scholar] [CrossRef]
  15. Takahashi, M.; Inoue, A.; Aritomi, M.; Takenaka, Y.; Suzuki, K. Gas entrainment at free surface of liquid, (I). J. Nucl. Sci. Technol. 1988, 25, 131–142. [Google Scholar] [CrossRef]
  16. Yamaguchi, A.; Tatsumi, E.; Takata, T.; Ito, K.; Hiroyuki Ohshima, H. Gas entrainment allowance level at free surface and gas dynamic behavior of sodium-cooled fast reactor. Nucl. Eng. Des. 2011, 241, 1627–1635. [Google Scholar] [CrossRef]
  17. Sakai, T.; Eguchi, Y.; Monji, H.; Iwasaki, T.; Ito, K.; Ohshima, H. Proposal of design criteria for gas entrainment from vortex dimples based on a computational fluid dynamics method. Heat Transfer Eng. 2010, 29, 731–739. [Google Scholar] [CrossRef]
  18. Kimura, N.; Ezure, T.; Tobita, A.; Kamide, H. Experimental study on gas entrainment at free surface in reactor vessel of a compact sodium-cooled fast reactor. J. Nucl. Sci. Technol. 2008, 45, 1053–1062. [Google Scholar] [CrossRef]
  19. Yamaguchi, A.; Takata, T.; Tatsumi, E.; Ito, K.; Ohshima, H.; Kamide, H. Rationalization of gas entrainment allowance level at free surface of Sodium-cooled fast reactor. In Proceedings of the NUTHOS-7: The 7th International Topical Meeting on Nuclear Reactor Thermal Hydraulics, Operation and Safety, Seoul, Korea, 5–9 October 2008. [Google Scholar]
  20. Ito, K.; Kunugi, T.; Ohshima, H.; Kawamura, T. Formulations and validations of a high-precision volume-of-fluid algorithm on nonorthogonal meshes for numerical simulations of gas entrainment phenomena. J. Nucl. Sci. Technol. 2009, 46, 366–373. [Google Scholar] [CrossRef]
  21. Ito, K.; Sakai, T.; Eguchi, Y.; Monji, H.; Ohshima, H.; Uchibori, A.; Xu, Y. Improvement of gas entrainment prediction method-introduction of surface tension effect. J. Nucl. Sci. Technol. 2010, 47, 771–778. [Google Scholar] [CrossRef]
  22. Kimura, N.; Ezure, T.; Miyakoshi, H.; Kamide, H.; Fukuda, T. Experimental study on gas entrainment due to nonstationary vortex in a sodium cooled fast reactor-comparison of onset conditions between sodium and water. J. Eng. Gas Turbines Power 2010, 132, 102908. [Google Scholar] [CrossRef]
  23. Koizumi, Y.; Ito, K.; Ohshima, H.; Ohtake, H. Study on gas entrainment rate into liquid from free surface by Vortex. In Proceedings of the ASME, International Mechanical Engineering Congress and Exposition, Vancouver, BC, Canada, 12–18 November 2010. [Google Scholar]
  24. Monji, H.; Shinozaki, T.; Kamide, H.; Sakai, T. Effect of experimental conditions on gas core length and downward velocity of free surface vortex in cylindrical vessel. J. Eng. Gas Turbines Power 2010, 132, 012901. [Google Scholar] [CrossRef]
  25. Ezure, T.; Kimura, N.; Miyakoshi, H.; Kamide, H. Experimental investigation on bubble characteristics entrained by surface vortex. Nucl. Eng. Des. 2011, 241, 4575–4584. [Google Scholar] [CrossRef]
  26. Ito, K.; Ohshima, H.; Nakamine, Y.; Imai, Y. Study on turbulent modeling in gas entrainment evaluation method. J. Power Energy Syst. 2012, 6, 151–164. [Google Scholar] [CrossRef] [Green Version]
  27. Koizumi, Y.; Ohte, N.; Hideki, K.; Ohno, S.; Ito, K. Measurement of gas entrainment rate from free surface by vortex. In Proceedings of the 2012 20th International Conference on Nuclear Engineering collocated with the ASME 2012 Power Conference, ICONE20-POWER2012, Anaheim, CA, USA, 30 July–3 August 2012. [Google Scholar]
  28. Ito, K.; Ezure, T.; Ohno, S.; Kamide, H. Evaluation model of bubble-type gas entrainment. In Proceedings of the the 15th International Topical Meeting on Nuclear Reactor Thermal-Hydraulics, NURETH15, Pisa, Italy, 12–17 May 2013. [Google Scholar]
  29. Ito, K.; Kunugi, T.; Ohshima, H. High-precision numerical scheme for vortical flow. Appl. Math. 2013, 4, 17–25. [Google Scholar] [CrossRef] [Green Version]
  30. Naosuke, O.; Koizumi, Y.; Hideki, K.; Ohno, S.; Ito, K. Effect of physical properties on gas entrainment rate from free surface by vortex. In Proceedings of the 2013 21st International Conference on Nuclear Engineering ICONE21, Chengdu, China, 29 July–2 August 2013. [Google Scholar]
  31. Winterton, R.H.S. Cover-gas bubbles in recirculating sodium primary coolant. Nucl. Eng. Des. 1972, 22, 262–271. [Google Scholar] [CrossRef]
  32. Eguchi, Y.; Yamamoto, K.; Funada, T.; Tanaka, N.; Moriya, S.; Tanimoto, K. Gas entrainment in the IHX vessel of top-entry. Nucl. Eng. Des. 1994, 146, 373–381. [Google Scholar] [CrossRef]
  33. Patwardhan, A.W.; Mali, R.G.; Jadhao, S.B.; Bhor, K.D.; Padmakumar, G.; Vaidyanathan, G. Argon entrainment into liquid sodium in fast breeder reactor. Nucl. Eng. Des. 2012, 249, 204–211. [Google Scholar] [CrossRef]
  34. Tenchine, D.; Fournier, C.; Dolias, Y. Gas entrainment issues in sodium cooled fast reactors. Nucl. Eng. Des. 2014, 270, 302–311. [Google Scholar] [CrossRef]
  35. Kim, S.N.; Jang, W.H. A study on the free surface vortex in the pipe system. J. Korean Nucl. Soc. 1992, 24, 311–319. [Google Scholar]
  36. Blömeling, F. Final Report: Generische numerische Untersuchungen zur Bestimmung der Mindestüberdeckung von Pumpenzuläufen zur Vermeidung von Luftmitriss; TÜV NORD SysTec GmbH & Co. KG: Hamburg, Germany, 2014. [Google Scholar]
  37. Pandazis, P.; Babcsány, B. Numerical and experimental investigation of surface vortex formation in coolant reservoirs of reactor safety systems. Kerntechnik 2016, 81, 477–483. [Google Scholar] [CrossRef]
  38. Škerlavaj, A.; Lipej, A.; Ravnik, J.; Škerget, L. Turbulence model comparison for a surface vortex simulation. Iop Conf. Ser. Earth Environ. Sci. 2010, 12, 012034. [Google Scholar] [CrossRef]
  39. Škerlavaj, A.; Škerget, L.; Ravnik, J.; Lipej, A. Predicting Free-Surface Vortices with Single-Phase Simulations. Eng. Appl. Comput. Fluid Mech. 2014, 8, 193–210. [Google Scholar] [CrossRef] [Green Version]
  40. Merzari, E.; Ninokata, H.; Wang, S.; Baglietto, E. Numerical Simulation of Free-Surface Vortices. Nucl. Technol. 2009, 165, 313–320. [Google Scholar] [CrossRef]
  41. Shi, X.M.; Yang, F.; Dai, R.; Chen, T.J.; Wu, Y.L. Simulation of free-surface vortex produced by a rotating cylindrical wall below a static barrel. Iop Conf. Ser. Earth Environ. Sci. 2012, 15, 052034. [Google Scholar] [CrossRef]
  42. Gauss, F.; Lucas, D.; Krepper, E. Grid studies for the simulation of resolved structures in an Eulerian two-fluid framework. Nucl. Eng. Des. 2016, 305, 371–377. [Google Scholar] [CrossRef]
  43. Cristofano, L.; Nobili, M.; Romano, G.P.; Caruso, G. Investigation on bathtub vortex flow field by Particle Image Velocimetry. Exp. Therm. Fluid Sci. 2016, 74, 130–142. [Google Scholar] [CrossRef]
  44. Hänsch, S.; Lucas, D.; Krepper, E.; Höhne, T. A multi-field two-fluid concept for transitions between different scales of interfacial structures. Int. J. Multiph. Flow 2012, 47, 171–182. [Google Scholar] [CrossRef]
  45. Hänsch, S.; Lucas, D.; Höhne, T.; Krepper, E. Application of a new concept for multi-scale interfacial structures to the dam-break case with an obstacle. Nucl. Eng. Des. 2014, 279, 171–181. [Google Scholar] [CrossRef]
  46. Montoya, G.A. Development and validation of advanced theoretical modeling for churn-turbulent flows and subsequent transitions. In Wissenschaftlich-Technische Berichte, Hemholtz-Zentrum Dresden-Rossendorf (HZDR); HZDR: Dresden, Germany, 2015. [Google Scholar]
  47. Höhne, T.; Krepper, E.; Lucas, D.; Montoya, G. A Multiscale Approach Simulating Boiling in a Heated Pipe Including Flow Pattern Transition. Nucl. Technol. 2019, 205, 48–56. [Google Scholar] [CrossRef]
  48. Höhne, T.; Krepper, E.; Montoya, G.; Lucas, D. CFD-simulation of boiling in a heated pipe including flow pattern transitions using the GENTOP concept. Nucl. Eng. Des. 2017, 322, 165–176. [Google Scholar] [CrossRef]
  49. ANSYS. ANSYS CFX-Solver Theory Guide, Release 19.2; ANSYS: Canonsburg, PA, USA, 2019. [Google Scholar]
  50. Lucas, D.; Rzehak, R.; Krepper, E.; Ziegenhein, T.; Liao, Y.; Kriebitzsch, S.; Apanasevich, P. A strategy for the qualification of multi-fluid approaches for nuclear reactor safety. Nucl. Eng. Des. 2016, 299, 2–11. [Google Scholar] [CrossRef]
  51. Putra, R.A.; Schäfer, T.; Neumann, M.; Lucas, D. CFD studies on the gas-liquid flow in the swirl generating device. Nucl. Eng. Des. 2018, 332, 213–225. [Google Scholar] [CrossRef]
  52. Putra, A.R.; Neumann-Kipping, M.; Schäfer, T.; Lucas, D. Comparison of Gas–Liquid Flow Characteristics in Geometrically Different Swirl Generating Devices. Energies 2019, 12, 4653. [Google Scholar] [CrossRef] [Green Version]
  53. Ishii, M.; Zuber, N. Drag coefficient and relative velocity in bubbly, droplet or particulate flows. Aiche J. 1979, 25, 843–855. [Google Scholar] [CrossRef]
  54. Ẑun, I. The transverse migration of bubbles influenced by walls in vertical bubbly flow. Int. J. Multiph. Flow 1980, 6, 583–588. [Google Scholar] [CrossRef]
  55. Tomiyama, A.; Tamai, H.; Zun, I.; Hosokawa, S. Transverse migration of single bubbles in simple shear flows. Chem. Eng. Sci. 2002, 57, 1849–1858. [Google Scholar] [CrossRef]
  56. Wellek, R.M.; Agrawal, A.K.; Skelland, A.H.P. Shape of liquid drops moving in liquid media. Aiche J. 1966, 12, 854–862. [Google Scholar] [CrossRef]
  57. 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, 635–652. [Google Scholar] [CrossRef]
  58. Tomiyama, A.; Sou, A.; Zun, I.; Kanami, N.; Sakaguchi, T. Effects of Eötvös number and dimensionless liquid volumetric flux on lateral motion of a bubble in a laminar duct flow. Multiphase flow 1995 1995, 1995, 3–15. [Google Scholar]
  59. Hosokawa, S.; Tomiyama, A.; Misaki, S.; Hamada, T. Lateral Migration of Single Bubbles Due to the Presence of Wall. In Proceedings of the ASME Joint U.S.-European Fluids Engineering Division Conference. FEDSM 2002, Montreal, QC, Canada, 14–18 July 2002. [Google Scholar]
  60. Burns, A.D.; Frank, T.; Hamill, I.; Shi, J.-M. The Favre Averaged Drag Model for Turbulence Dispersion in Eulerian Multi-Phase Flows; ICMF2004: Yokohama, Japan, 2004. [Google Scholar]
  61. Auton, T.R.; Hunt, J.C.R.; Prud’Homme, M. The force exerted on a body in inviscid unsteady non-uniform rotational flow. J. Fluid Mech. 2006, 197, 241–257. [Google Scholar] [CrossRef]
  62. Maxey, M.R.; Riley, J.J. Equation of motion for a small rigid sphere in a nonuniform flow. Phys. Fluids 1983, 26, 883–889. [Google Scholar] [CrossRef]
  63. Magnaudet, J.; Rivero, M.; Fabre, J. Accelerated flows past a rigid sphere or a spherical bubble. Steady straining flow. J. Fluid Mech. 2006, 284, 97–135. [Google Scholar] [CrossRef]
  64. Brackbill, J.U.; Kothe, D.B.; Zemach, C. A continuum method for modeling surface tension. J. Comput. Phys. 1992, 100, 335–354. [Google Scholar] [CrossRef]
  65. Ma, J.; Oberai, A.A.; Drew, D.A.; Lahey, R.T., Jr.; Hyman, M.C. A comprehensive sub-grid air entrainment model for RANS modeling of free surface bubbly flows. J. Comput. Multiph. Flows 2011, 3, 16. [Google Scholar] [CrossRef] [Green Version]
  66. Menter, F.R. Two-equation eddy-viscosity turbulence models for engineering applications. Aiaa J. 1994, 32, 1598–1605. [Google Scholar] [CrossRef] [Green Version]
  67. Smirnov, P.E.; Menter, F.R. Sensitization of the SST Turbulence Model to Rotation and Curvature by Applying the Spalart–Shur Correction Term. J. Turbomach. 2009, 131, 041010. [Google Scholar] [CrossRef]
  68. Spalart, P.R.; Shur, M. On the sensitization of turbulence models to rotation and curvature. Aerosp. Sci. Technol. 1997, 1, 297–302. [Google Scholar] [CrossRef]
  69. Rzehak, R.; Krepper, E. Bubble-induced turbulence: Comparison of CFD models. Nucl. Eng. Des. 2013, 258, 57–65. [Google Scholar] [CrossRef]
  70. Rzehak, R.; Krepper, E. CFD modeling of bubble-induced turbulence. Int. J. Multiph. Flow 2013, 55, 138–155. [Google Scholar] [CrossRef]
  71. Ikoma, Y.; Koizumi, Y.; Ito, K.; Ohshima, H. Bubble-Type Gas Entrainment into Liquid from Free Surface by Vortex. In Proceedings of the ICONE-19 the 19th International Conference on Nuclear Engineering, Osaka, Japan, 24–25 October 2011. [Google Scholar]
  72. Ito, K.; Koizumi, Y.; Ohshima, H.; Kawamura, T. Physics-basis simulation of bubble pinch-off. Mech. Eng. J. 2016. [Google Scholar] [CrossRef] [Green Version]
  73. Zidouni, F.; Krepper, E.; Rzehak, R.; Rabha, S.; Schubert, M.; Hampel, U. Simulation of gas–liquid flow in a helical static mixer. Chem. Eng. Sci. 2015, 137, 476–486. [Google Scholar] [CrossRef]
Figure 1. Free-surface vortex observed in Akkats hydropower station, Sweden. The image is taken from [10].
Figure 1. Free-surface vortex observed in Akkats hydropower station, Sweden. The image is taken from [10].
Water 12 00709 g001
Figure 2. Pressure balance for a stable gas core.
Figure 2. Pressure balance for a stable gas core.
Water 12 00709 g002
Figure 3. (a) Graphical representation of the computational domain and mesh; (b) CFD boundaries.
Figure 3. (a) Graphical representation of the computational domain and mesh; (b) CFD boundaries.
Water 12 00709 g003
Figure 4. (a) Average interfacial deformation represented by iso-surface of α c g = 0.5 and (b) average bubble entrainment represented by iso-surface of α d g = 0.02. The results are obtained from the simulations using C e n t = 2 × 10−4.
Figure 4. (a) Average interfacial deformation represented by iso-surface of α c g = 0.5 and (b) average bubble entrainment represented by iso-surface of α d g = 0.02. The results are obtained from the simulations using C e n t = 2 × 10−4.
Water 12 00709 g004
Figure 5. Transient profile of gas entrainment rate of: (a) potentially continuous gas cg, (b) dispersed gas dg, and (c) total gas. (d) Average gas entrainment rate of cg and dg and (e) total. The results are obtained from the simulations using C e n t = 2 × 10−4. The experimental gas entrainment rate is based on [71,72].
Figure 5. Transient profile of gas entrainment rate of: (a) potentially continuous gas cg, (b) dispersed gas dg, and (c) total gas. (d) Average gas entrainment rate of cg and dg and (e) total. The results are obtained from the simulations using C e n t = 2 × 10−4. The experimental gas entrainment rate is based on [71,72].
Water 12 00709 g005
Figure 6. Velocity vectors represent the tangential projection of the average liquid velocity on horizontal planes for cases: (a) C s c a l e = 0, (b) C s c a l e = 0.5, (c) C s c a l e = 1. The background color represents the average liquid vorticity magnitude. The results are obtained from the simulations using C e n t = 2 × 10−4. Horizontal planes H1, H2, and H3 are located 10, 20, and 30 mm above the bottom surface of the cylindrical vessel, respectively.
Figure 6. Velocity vectors represent the tangential projection of the average liquid velocity on horizontal planes for cases: (a) C s c a l e = 0, (b) C s c a l e = 0.5, (c) C s c a l e = 1. The background color represents the average liquid vorticity magnitude. The results are obtained from the simulations using C e n t = 2 × 10−4. Horizontal planes H1, H2, and H3 are located 10, 20, and 30 mm above the bottom surface of the cylindrical vessel, respectively.
Water 12 00709 g006
Figure 7. The contour of the dispersed gas fraction at t = 6 s plotted on the vertical central plane for the simulations of (a) Sim. A, (b) Sim. B, and (c) Sim. C. The results are obtained from the simulations using C e n t = 2× 10−4 and C s c a l e = 1.
Figure 7. The contour of the dispersed gas fraction at t = 6 s plotted on the vertical central plane for the simulations of (a) Sim. A, (b) Sim. B, and (c) Sim. C. The results are obtained from the simulations using C e n t = 2× 10−4 and C s c a l e = 1.
Water 12 00709 g007
Figure 8. Transient profile of gas entrainment rate of: (a) potentially continuous gas, (b) dispersed gas, and (c) total gas. Average gas entrainment rate of: (d) cg and dg and (e) total. The results are obtained from the simulations using C e n t = 2 × 10−4. The experimental gas entrainment rate is based on [71,72].
Figure 8. Transient profile of gas entrainment rate of: (a) potentially continuous gas, (b) dispersed gas, and (c) total gas. Average gas entrainment rate of: (d) cg and dg and (e) total. The results are obtained from the simulations using C e n t = 2 × 10−4. The experimental gas entrainment rate is based on [71,72].
Water 12 00709 g008
Figure 9. (a) Average interfacial deformation represented by isosurface of α c g = 0.5 and (b) average bubble entrainment represented by isosurface of α d g = 0.02. The results are obtained from the simulations using C s c a l e = 1.
Figure 9. (a) Average interfacial deformation represented by isosurface of α c g = 0.5 and (b) average bubble entrainment represented by isosurface of α d g = 0.02. The results are obtained from the simulations using C s c a l e = 1.
Water 12 00709 g009
Figure 10. Transient profile of gas entrainment rate of: (a) potentially continuous gas, (b) dispersed gas and (c) total gas. (d) Average gas entrainment rate of cg and dg and (e) total. The results are obtained from the simulations using C s c a l e = 1. The experimental gas entrainment rate is based on [71,72].
Figure 10. Transient profile of gas entrainment rate of: (a) potentially continuous gas, (b) dispersed gas and (c) total gas. (d) Average gas entrainment rate of cg and dg and (e) total. The results are obtained from the simulations using C s c a l e = 1. The experimental gas entrainment rate is based on [71,72].
Water 12 00709 g010
Figure 11. Velocity vectors represent the tangential projection of the average liquid velocity on horizontal planes for cases (a) C e n t = 1× 10−4 and (b) C e n t = 2 × 10−4. The background color represents the average liquid vorticity magnitude. The results are obtained from the simulations using C s c a l e = 1. Horizontal planes H1, H2, and H3 are located 10, 20, and 30 mm above the bottom surface of the cylindrical vessel, respectively.
Figure 11. Velocity vectors represent the tangential projection of the average liquid velocity on horizontal planes for cases (a) C e n t = 1× 10−4 and (b) C e n t = 2 × 10−4. The background color represents the average liquid vorticity magnitude. The results are obtained from the simulations using C s c a l e = 1. Horizontal planes H1, H2, and H3 are located 10, 20, and 30 mm above the bottom surface of the cylindrical vessel, respectively.
Water 12 00709 g011
Figure 12. Average interfacial deformation, represented by isosurface, of α c g = 0.5 for cases (a) Mesh A and C s c a l e = 1, (b) Mesh B and C s c a l e = 1 and (c) Mesh B and C s c a l e = 0.85.
Figure 12. Average interfacial deformation, represented by isosurface, of α c g = 0.5 for cases (a) Mesh A and C s c a l e = 1, (b) Mesh B and C s c a l e = 1 and (c) Mesh B and C s c a l e = 0.85.
Water 12 00709 g012
Table 1. Closure models for the polydispersed gas phase based on the baseline model concept [50]. The table is taken from [52].
Table 1. Closure models for the polydispersed gas phase based on the baseline model concept [50]. The table is taken from [52].
ForceFormulationRef.No.
Drag F d r a g = 3 4 d B C D ρ L α G | u G u L | ( u G u L ) [53](5)
C D , b u b b = m a x ( C D , s p h e r e , m i n ( C D , e l l i p s e , C D , c a p ) ) ,(6)
C D , s p h e r e = 24 R e ( 1 + 0.1 R e 0.75 ) , C D , e l l i p s e = 2 3 E o , C D , c a p = 8 3
Lift F l i f t = C L ρ L α G ( u G u L ) × ( × u L ) [54](7)
C L { m i n   [ 0.288   t a n h ( 0.121 R e ) ,   f ( E o ) ]              E o < 4 f ( E o )                                  f o r      4 < E o < 10 0.27                                  10 < E o                  [55] (8)
f ( E o ) = 0.00105 E o 3 0.0159 E o 2 0.0204 E o + 0.474 (9)
E o = g ( ρ L ρ G ) d 2 σ (10)
d = d B 1 + 0.163   E o 0.757 3 [56](11)
Wall lubrication F w a l l = 2 d B C W ρ L α G | u G u L | 2 y ^ [57](12)
C W ( y ) = f ( E o ) ( d B 2 y ) 2 [58](13)
f ( E o ) = 0.0217 E o [59](14)
Turbulent dispersion F T D = 3 4 C D α G d B | u G u L | μ L t u r b σ T D ( 1 α L + 1 α G ) α G [60](15)
σ T D = 0.9 [49](16)
Virtual mass F V M = C V M ρ L α G ( D G u G D t D L u L D t ) [61,62,63](17)
C V M = 0.5 (18)
Table 2. Setup for gas phase and entrainment fraction.
Table 2. Setup for gas phase and entrainment fraction.
Velocity Groupsdgcg
MorphologyPolydispersedContinuous
Bubble classesG1G2G3G4G5
Diameter [mm]1357≥9
Entrainment fraction (Sim. A)1000−1
Entrainment fraction (Sim. B)0.500.5000−1
Entrainment fraction (Sim. C)0.250.250.250.25−1
Table 3. Gas core length, gas entrainment rate, and maximum vorticity.
Table 3. Gas core length, gas entrainment rate, and maximum vorticity.
Method/Model L g
(mm)
Gas Entrainment Rate (m3/s) Max .   ω (1/s) around the Vortex Tip
Ikoma exp. 2.0 × 10−8 [71,72]
Mesh A, C scale = 1162.3 × 10−8369
Mesh B, C scale = 1276.9 × 10−8620
Mesh B, C scale = 0.85162.4 × 10−8384

Share and Cite

MDPI and ACS Style

Putra, R.A.; Lucas, D. Modeling of the Free-Surface Vortex-Driven Bubble Entrainment into Water. Water 2020, 12, 709. https://doi.org/10.3390/w12030709

AMA Style

Putra RA, Lucas D. Modeling of the Free-Surface Vortex-Driven Bubble Entrainment into Water. Water. 2020; 12(3):709. https://doi.org/10.3390/w12030709

Chicago/Turabian Style

Putra, Ryan Anugrah, and Dirk Lucas. 2020. "Modeling of the Free-Surface Vortex-Driven Bubble Entrainment into Water" Water 12, no. 3: 709. https://doi.org/10.3390/w12030709

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop