Research Article  Open Access
CFD Simulation of Polydispersed Bubbly TwoPhase Flow around an Obstacle
Abstract
This paper concerns the model of a polydispersed bubble population in the frame of an ensemble averaged twophase flow formulation. The ability of the moment density approach to represent bubble population size distribution within a multidimensional CFD code based on the twofluid model is studied. Two different methods describing the polydispersion are presented: (i) a moment density method, developed at IRSN, to model the bubble size distribution function and (ii) a population balance method considering several different velocity fields of the gaseous phase. The first method is implemented in the Neptune_CFD code, whereas the second method is implemented in the CFD code ANSYS/CFX. Both methods consider coalescence and breakup phenomena and momentum interphase transfers related to drag and lift forces. Airwater bubbly flows in a vertical pipe with obstacle of the TOPFLOW experiments series performed at FZD are then used as simulations test cases. The numerical results, obtained with Neptune_CFD and with ANSYS/CFX, allow attesting the validity of the approaches. Perspectives concerning the improvement of the models, their validation, as well as the extension of their applicability range are discussed.
1. Introduction
Many flow regimes in Nuclear Reactor Safety Research are characterized by multiphase flows, with one phase being a liquid and the other phase consisting of gas or vapor of the liquid phase. The flow regimes found in vertical pipes are dependent on the void fraction of the gaseous phase, which varies, as void fraction increases, from bubbly flow to slug flow, churn turbulent flow, annular flow, and finally to droplet flow at highest void fractions. In the regimes of bubbly and slug flows, a spectrum of different bubble sizes is observed. While dispersed bubbly flows with low gas volume fraction are mostly monodispersed, an increase of the gas volume fraction leads to a broader bubble size distribution due to breakup and coalescence of bubbles. The exchange area for mass, momentum, or heat between continuous and dispersed phases thus cannot be simply modeled based on the knowledge of just the void fraction and a mean diameter. Moreover, the forces acting on the bubbles may depend on their individual size which is the case not only for drag but also for nondrag forces. Among the forces leading to lateral migration of the bubbles, that is, acting in normal direction with respect to the main drag force, bubble lift force was found to change the sign as the bubble size varies. Consequently, in the context of pipe flows, this leads to a radial separation between small and large bubbles and to further coalescence of large bubbles migrating toward the pipe center into even larger Taylor bubbles or slugs.
An adequate modeling approach must consider all these phenomena. The paper presents two different approaches both being based on the Eulerian modeling framework. On one hand, a generalized inhomogeneous multiple size group (MUSIG) model was applied, for which the dispersed gaseous phase is divided into N inhomogeneous velocity groups (phases), each of these groups being subdivided into bubble size classes. On the other hand, a moment density method is used to model the evolution of the bubble size distribution function along the flow, where main statistics of the distribution being described with a reduced set of transport equations. For both models, bubble breakup and coalescence processes are taken into account by appropriate models.
2. Model for the Local Size Distribution of a Bubble Population
2.1. Closure Models for Momentum Exchange and for Bubble Coalescence and Breakup
While modeling a twophase flow using the Euler/Eulerian approach, the momentum exchange between the phases has to be considered. Apart from the drag acting in flow direction, the socalled nondrag forces acting mainly perpendicularly to the flow direction must be considered. Namely, the lift force, the turbulence dispersion force, the virtual mass force, and the wall force play an important role.
The turbulent dispersion force acts on smoothing the gas volume fraction distribution and can be evaluated either from a single expression related to drag turbulent contribution, [1] or deduced from a model for dispersed phase agitation contribution to the momentum balance, for example, the Tchen model [2]. To avoid the maximum gas fraction at the wall, Tomiyama et al. [3, 4] propose a wall force.
The lift force considers the interaction of the bubble with the shear field of the liquid. For a single bubble, it reads where is the liquid density, the bubble diameter, and and are the bubble and liquid velocities, respectively. The classical lift force, which has a positive coefficient , acts in the direction of decreasing liquid velocity. In case of cocurrent upwards pipe flow, this is the direction towards the pipe wall. Numerical [5] and experimental [3] investigations showed that the direction of the lift force changes its sign if a substantial deformation of the bubble occurs. Tomiyama [4] investigated single bubble motion and derived the following correlation for the coefficient of the lift force from these experiments: This coefficient depends on the modified Eötvös number which is formulated for the maximum horizontal bubble size in flow direction [6] and on the Reynolds number based on bubble size. The sign of the lift coefficient, and, consequently, the direction of the lift force, depends, thus, on the bubble size diameter. This behavior, originally found for single bubbles of air in glycerol, was also established by different experiments for gaswater polydispersed flows (e.g., [7]) and for different fluids too.
For several flow configurations, this bubble size dependency of the lift force direction can lead to the separation between small and large bubbles. This effect has been shown to be a key phenomenon for the development of the flow regime.
The bubble coalescence model takes into account the random collision processes between two bubbles. The applied model is based on the work of Prince and Blanch [8]. The bubble breakup model considers the collision between a liquid turbulent eddy of a certain characteristic size and a bubble (see [9]). Coalescence and breakup mechanistic models are common to both methods for taking into account the polydispersion in size of the bubble population. We introduce corresponding “breakup” and “coalescence” numerical coefficients and used to scale the original correlations.
2.2. Population Balance Approach
2.2.1. The MUSIG Model by Lo
In principle, the Eulerian twofluid approach can be extended to simulate a continuous liquid phase and several gaseous dispersed phases solving the complete set of balance equations for each phase. The investigations, however, showed that for an adequate description of the gas volume fraction profile including a population balance model decades of bubble size classes would be necessary. In a CFD code, such a procedure is limited by the increased computational effort to obtain converged flow solutions. To solve this problem, the multiple size group model first implemented by the code developers in CFX4 solves only one common momentum equation for all bubble size classes (homogeneous MUSIG model, see [10], Figure 1). Mathematically, the multiple size group model (MUSIG) is based on the population balance method and the twofluid modeling approach. The dispersed phase is divided into size fractions. The population balance equation is applied to describe the mass conservation of the size fractions taking into account the interfraction mass transfer resulting from bubble coalescence and breakup. This model approach allows a sufficient number of size fraction groups required for the coalescence and breakup calculation to be used and has found a number of successful applications to largescale industrial multiphase flow problems.
Nevertheless, the assumption also restricts its applicability to homogeneous dispersed flows, where the slip velocities of particles are almost independent of particle size and the particle relaxation time is sufficiently small with respect to inertial time scales. Thus, the asymptotic slip velocity can be considered to be attained almost instantaneously. The homogeneous MUSIG model described above fails to predict the correct phase distribution when heterogeneous particle motion becomes important. One example is the bubbly flow in vertical pipes, where the nondrag forces play an essential role on the bubble motion. The lift force was found to change its sign, when applied for large deformed bubbles, which are dominated by the asymmetrical wake ([3], see Section 2.1). The lift force in this case has a direction opposite to the shearinduced lift force on a small bubble. For this reason, large bubbles tend to move to the pipe core region resulting in a core void maximum, whereas a nearwall void peak is measured for small bubbles. The radial separation of small and large bubbles cannot be predicted by the homogeneous MUSIG model. This has been shown to be a key mechanism for the establishment of a certain flow regime.
2.2.2. New Strategy: the Inhomogeneous MUSIG Model
A combination of the consideration of different dispersed phases and the algebraic multiple size group model was proposed to combine both the adequate number of bubble size classes for the simulation of coalescence and breakup and a limited number of dispersed gaseous phases to limit the computational effort [11]. The inhomogeneous MUSIG model was developed in cooperation with ANSYS CFX and is implemented in CFX since the version CFX10 ([12–14] see Figure 2).
In the inhomogeneous MUSIG model, the gaseous dispersed phase is divided into a number of the socalled velocity groups (or phases), where each of the velocity groups is characterized by its own velocity field. Further, the overall bubble size distribution is represented by dividing the bubble diameter range within each of the velocity groups in a number , , bubble subsize fractions. The population balance model, considering bubble coalescence or bubble breakup is applied to the subsize groups. Hence the mass exchange between the subsize groups can exceed the size ranges assigned to the velocity groups resulting in mass transfer terms between the different phases or velocity groups.
The lower and upper boundaries of bubble diameter intervals for the bubble size fractions can be controlled by either an equal bubble diameter distribution, an equal bubble mass distribution, or can be based on user definition of the bubble diameter ranges for each distinct bubble diameter fraction. The subdivision should be based on the physics of bubble motion for bubbles of different size, for example, different behavior of differently sized bubbles with respect to lift force or turbulent dispersion. Extensive model validation calculations have shown that in most cases, or , velocity groups are sufficient in order to capture the main phenomena in bubbly or slug flows [15, 16].
2.3. The Moment Density Method Approach
The moment density method proposes an alternative way to model polydispersion. This method allows to model the time and space evolution of a realistic distribution with the help of a very reduced number of transport equations, for example, Kamp et al. [17], Hill [18]. This is an interesting property, especially from a numerical point of view, with regard to alternative methods like population balance methods. In the present study, we restrict our attention to adiabatic bubbly flows of an incompressible gas inside a continuous liquid phase, and only consider a dispersion in bubbles size. We here present the basic formalism of the moment density method.
2.3.1. A Distribution Function Parameterized by Its Statistical Moments
The moment density method requires the use of an approximate representation of the bubble population thanks to a presumed shape continuous function for the bubble size distribution function , the one being thus totally determined by a finite number of parameters. The moment densities are intensive physical quantities characterizing mean properties of the population, like the density of bubble area or volume, the socalled volumetric interfacial area and the volumetric fraction of the dispersed phase. Bubbles being assumed as spherical and of diameter , it yields
In this work, we introduce as being simply quadratic in the bubble size , where is the local density number of bubbles and is the mean diameter of the distribution. It is easy to show that , , , and are related through: The knowledge of the two moment densities and is thus sufficient to reconstruct the whole distribution function . The choice for these two parameters is dictated by their importance in the description of the twophase flow, especially concerning the transfers between the phases.
2.3.2. Transport Equations for the Size Distribution
The transport of the distribution function obeys a LiouvilleBoltzmann equation from which one can derive transport equations for the main moment densities as well as for their transport velocities [19]. The main idea of the method consists in solving a finite number of moment densities transport equations in the framework of the twofluid model. This solving allows to determine the distribution parameters and thus the whole population evolution. Thus we consider the solving of the following transport equations: where is the source term of the LiouvilleBoltzmann equation. The transport velocities and of and , respectively, are assumed to be equal to the mean transport velocity of the bubbles (no correlation between dispersions in velocity and in size is taken into account at the current stage of development of the model). is defined by where is the velocity of a bubble of size . It satisfies the transport equation:
2.3.3. Closure Relations
The closure of the system of (3)–(8) concerns models (i) for coalescence and breakup events through , (ii) for the bubble size evolution along its trajectory , (iii) for the hydrodynamic forces acting on individual bubbles, through , as well as (iv) for the kinetic stress tensor . The closure issue is based on the mechanistic model of all these phenomena at the scale of a single bubble. Source terms of the equations are then evaluated from the integral of individual contributions over the whole bubble population, whose size distribution function is known. The simplicity of the analytical expression (4) for allows to easily derive expressions for the integrals from the classical models described in Section 2.1. In this study, since the dispersed phase is considered as incompressible and incondensable, the bubble size evolution along its trajectory is zero. Detailed expressions for the source terms can be found in [20] .
It is worth pointing out that solving the system of (3)–(8), coupled with the system of balance equations for the continuous phase, is sufficient to describe the evolution of a broad spectrum of bubble sizes. This simplicity of the formulation is advantageous especially in terms of computational cost. Further extension of the model, taking into account the dispersion in velocity, can be considered using a similar formalism. For example, sizerelated velocity distribution can be introduced to model sizedependent bubble migration due to lift force, like it is done for the inhomogeneous MUSIG model. Moreover, this method has already been successful to consider more general particles velocity dispersion, for example, the work of Fox [21] concerning crossing trajectories.
3. The Experiment–Bubbly Flow Around an Obstacle
In the presented experiment performed at the TOPFLOW facility of Research center DresdenRossendorf (FZD), the large test section with a nominal diameter of DN200 was used to study the flow field around an asymmetric obstacle (see Figure 3). This is an ideal test case for CFD code validation, since the obstacle creates a pronounced threedimensional twophase flow field. Curved stream lines, which form significant angles with the gravity vector, a recirculation zone in the wake and a flow separation at the edge of the obstacle are common in industrial components and installations.
The wiremesh technology was applied to measure the gas volume fraction and the gas velocity in different distances up and downstream the obstacle [22]. The sensor provides detailed data on the instantaneous flow structure with a high resolution in space and time. In particular, they allow visualizing the structure of the gasliquid interface [23].
The tests were performed both for air/water and for steam/water. In the current paper, only adiabatic air/water tests were considered. The parameters are summarized in Table 1.

4. Numerical Settings
4.1. CFX and Population Balance Method
Pretest calculations using ANSYS/CFX and applying a monodispersed bubble size approach were performed for the conditions of test run 074 ( m/s, m/s) (see [24–26]). In the calculation, a fluid domain was modeled 1.5 m upstream and downstream the obstacle. Half of the tube including a symmetry boundary condition set at the plane of the geometry was simulated. In the present paper, the inhomogeneous MUSIG model approach was applied to air/water obstacle experiments run 096 ( m/s, m/s) and run 097 ( m/s, m/s). In the presented calculations for run 096 and run 097, 25, and 20, respectively, subsize gas fractions representing equidistant bubble sizes up to 25 mm and 20 mm, respectively, were simulated, assigned to 2 dispersed gaseous phases. The first 6 size groups were assigned to the first gaseous phase (or velocity group) and the remaining size groups were assigned to the second gaseous phase. The bubble size distribution measured at the largest upstream position was set as an inlet boundary condition for the calculation.
4.2. Neptune_CFD and Moment Density Method
Neptune_CFD code [27], which is based on the twofluid model formulation, has been used to perform numerical simulations using the moment density method. The test run 074 has been simulated using a calculation domain corresponding to a 1 m long half pipe centered on the obstacle of around 150 000 cells that takes benefit of the symmetry of the experimental setup. Neptune_CFD calculations of test 074 are based on a  model for the liquid flow, the Tchen correlation being used to model the fluctuations related to bubbles turbulent motion. Bubble hydrodynamics model include both laminar and turbulent evaluations of drag, lift, and virtual mass forces, thanks to a drift model for the relative velocity. Inlet boundary conditions correspond to flat profiles for and the velocities, allowing to recover the experimental surface averaged values for superficial velocities and mean bubble diameter. The surface averaged volumetric fraction and mean diameter are used to reconstruct the inlet bubble size distribution function. The numerical method is based on a pressurebased method. Mass, momentum, and energy equations are coupled by an iterative procedure within a time step. Spatial discretization is based on finite volume framework on unstructured meshes, all the variables being computed at the center of cells.
5. Comparison of Measured and Calculated Results
5.1. The Main Observed Phenomena
Both the steadystate ANSYS CFX calculations applying the inhomogeneous MUSIG model, and the Neptune_CFD calculations applying the moment density method, could reproduce all qualitative details of the flow structure of the twophase flow field around the diaphragm. The structure of the flow for the here considered test cases 074, 096, and 097 are essentially similar. These different tests have been selected for the purpose of investigating certain phenomena which are more or less pronounced.
The numerical results have been compared to threedimensional wiremesh sensor data in Figure 4 (ANSYS CFX for the run 096) and Figure 5 (Neptune_CFD for the run 074). The water velocity and the total gaseous void fraction are represented. All qualitative details of the structure of the twophase flow field around the obstacle could be reproduced.
(a)
(b)
Shortly, behind the obstacle a strong vortex of the liquid combined with the accumulation of gas is observed. The measured and calculated shape and extension of the recirculation area agree very well. Upstream the obstacle, a stagnation point with lower gas content is seen in experiment and calculation. Details, like the velocity and void fraction maxima above the gap between the circular edge of the obstacle and the inner wall of the pipe, are also found in a good agreement between experiments and calculations. In the unobstructed cross sectional part of the tube a strong jet is established. Main discrepancies between experimental and Neptune_CFD results concern the volumetric fraction upward (below) the obstacle and can be associated to the flat profiles used as inlet bottom boundary conditions.
The structure of the flow is studied in more detail in the following sections.
5.2. Phenomena in the Wake of the Obstacle
5.2.1. Size Distribution of the Bubble Population
More detailed understanding of the flow situation can be gained, considering the bubble size distribution. According to the applied bubble breakup model of Luo and Svendsen [9], bubble breakup can be expected in regions showing highturbulent eddy dissipation. Figure 6 presents maximum values of the numerically evaluated turbulent eddy dissipation at the edges of the obstacle. At the same time, the applied bubble coalescence model of Prince and Blanch [8] indicates strong importance of coalescence in regions of bubble accumulation, that is, in the wake behind the obstacle. Both bubble coalescence (see gas accumulation shown in Figures 4 and 5) and bubble breakup (see distribution of turbulence dissipation Figure 6), which might partially compensate each other, are expected shortly behind the obstacle.
Figure 7 presents measured crosssectional averaged bubble size distributions upstream ( m), shortly behind ( m) and downstream the obstacle ( m) for the run 096. The measurements show in the bubble accumulation zone at m the crosssectional average shows a shift toward larger bubbles. The experimentally measured bubble size distribution of run 074 is essentially similar to the run 096. It can be seen on the volumetric fraction maps, Figure 8. Above the obstacle, where turbulent intensity is strong, breakup is not experimentally attested, contrarily to the Luo and Svedsen [9] model prediction, but coalescence occurs, leading to the formation of big bubbles (over 7 mm) that concentrate thereafter in the central part of the pipe.
(a)
(b)
(c)
(d)
Both ANSYS CFX and Neptune_CFD calculated bubble size distributions, however, show a shift of the mean bubble diameter toward smaller bubbles shortly behind the obstacle when both coalescence and breakup are taken into account. In the calculations, the bubble breakup is overestimated. The corresponding results are provided by Figure 9 for the run 096 and on Figure 10 for run 097. This disagreement was found not solvable by simple changes of breakup or coalescence coefficients, which were set here to . Similar deviations would arise at other locations of the flow domain. Neptune_CFD results for the run 074 with standard coefficients are shown in Figure 11 (the details can be found in [20]).
This suggests to perform computations for which breakup is neglected, but that still consider the Prince and Blanch [8] coalescence model with standard coefficient. The corresponding Neptune_CFD results show that both bubble size repartition and order of magnitude are consistent with experimental results: big bubbles of more than 7 mm are created and accumulate in the recirculation zone just above the obstacle (see Figure 14). Big bubbles stay in the LHS of the pipe, the unobstructed crosssectional part of the pipe is free of big bubbles, that is consistent with the experimental observation.
As partial conclusions, (i) the Luo and Svendsen breakup model tends to overestimate the breakup of TOPFLOW experimental tests, whereas (ii) the use of the Prince and Blanch model for coalescence within the framework of the moment density method allows to predict satisfying evolution of bubble size across the flow.
5.2.2. Bubbles Streamlines
More detailed effects of lateral motion of small and large bubbles can be revealed by studying bubble streamlines and by analyzing lift forces acting on bubbles of different size. On one hand, the liquid velocity flow carries the small bubbles into the region behind the obstacle (see Figures 12 and 14, righthand side for the bubble streamlines). Lateral deviation due to lift force is illustrated on Figure 13 for the lift force arrows with the population balance method and Figure 14, righthand side for the lift coefficient with the moment density method. On the other hand, the air accumulation in the wake region leads to bubble coalescence and the generation of large bubbles as revealed by the analysis of experimental results. This phenomenon is underestimated in the calculations taking breakup into account.
(a)
(b)
Caused by the lift force, large bubbles are redirected into the downstream jet (see Figure 13) once they can be formed in the wake by coalescence. The streamline representation (see Figure 12) clearly shows this phenomenon for large bubbles already present in the upstream flow.
In the actual version of the moment density method used in Neptune_CFD, dynamics of the bubble population is estimated using a single transport velocity. The averaged lift contribution takes into account both the bubble size distribution and the Tomiyama correlation. When a majority of bubbles are locally above the critical Eötvös number, the lift coefficient changes its sign. In this case, the direction of the lift force is changed, as it can be seen on Figure 14. As a consequence, the bubble streamlines above the obstacle deviate to the center of the pipe. This explains the low value of the dispersed phase volumetric fraction near the lefthand side part of the pipe wall above the obstacle for Neptune_CFD numerical results (see Figure 5). This result is fully consistent with experimental results.
As a partial conclusion concerning the bubble streamlines calculations, both methods showed their ability to consider the effect of bubble size on the lateral deviation of bubble streamlines due to lift force. This provides a more precise understanding and a more accurate prediction of bubbles repartition across the flow.
5.3. Phenomena in the Jet
In the crosssectional area beside the obstacle, a strong jet is established creating strong shear flow. The resulting phenomena are more pronounced with increasing water velocity, like in run 097, where the liquid velocity was increased to m/s. Figure 16 represents measured and ANSYS CFX calculated crosssectional gas fraction distributions for this run. In the most downstream cross section of the measurements an almost gas bubble free region in the centre of the jet is found. Bubbles are collected at the edges of the jet in regions of largest water velocity gradient.
The streamline representation of the ANSYS CFX calculations, however, (Figure 12 for run 096, which is fairly similar to run 097) indicates large bubbles being directed into the jet caused by the lift force. This discrepancy between experiment and ANSYS CFX calculations can possibly be explained by the strong water velocity gradient near the jet. This strong shear flow induces bubble breakup which is not yet considered in the model of Luo and Svendsen [9]. In the tests, the big bubbles could migrate toward the jet, but be fragmented at the relatively sharp boarder of this jet. Only a small fraction of the small bubbles created by this breakup process can enter the jet by action of the turbulent dispersion force.
The gas distribution resolved by bubble size classes (see Figure 15) shows clearly the deviations between measurements and calculations. In the experiment, the gas accumulation behind the obstacle leads to a strong coalescence and the creation of large bubbles. In the CFXcalculations, this effect is underestimated.
(a)
(b)
In Neptune_CFD calculations based on the moment density method the bubble breakup was neglected. As attested by Figure 14, lefthand side, in this case the bubbles present in the jet are small bubbles. Big bubbles created are concentrated in the center of the pipe where the largest value of volumetric fraction can be found (Figure 5). The structure of the gas repartition in the jet for both calculation and experiments is represented on Figure 17. In the upper part, the structure of the flow of the numerical calculation is consistent with experimental results. The most important deviation is observed in the region just above the obstacle . A gain in accuracy of these results could be achieved by considering sizedependent velocity dispersion within the moment density method, in particular for bubble sizedependent lateral migration in the regions where there exists a strong mixing of different bubble sizes (mainly just above the obstacle as far as run 074 is concerned, see Figure 8).
Nevertheless, the main trends of the bubble dynamics downward the obstacle are in agreement with experimental results (see Figure 16 for CFX and Figure 17 for Neptune_CFD).
6. Summary and Perspectives
In this study, we focused on the model of the bubble size distribution in the numerical simulation of bubbly flows. More deep understanding of the flow structure is possible when considering a more accurate characterization of the polydispersion. For upward twophase flow in vertical pipes the core peak in the crosssectional gas fraction distribution could be reproduced very well both by the moment density and by the MUSIG population balance methods. For complex flows, the general threedimensional structure of the flow could be well reproduced in the simulations.
These test cases of pipe flow with internal obstacle demonstrate the complicated relationship and interference between sizedependent bubble migration, bubble coalescence, and breakup effects for real flows. With an appropriate given distribution function, the numerical effort of the moment density method is lower compared to the multiple bubble size group method (MUSIG). On the other hand, applying the MUSIG method the simulation of a flow situation allows to deal with more general shapes for the distribution function. Both methods enable to consider the effect of polydispersion in size on the bubble population dynamics, in particular on the evaluation of interphase momentum transfers associated with lift and drag. The inhomogeneous population balance model, using several velocity fields for the bubbly phase, is able to deal with size separation of a locally polydispersed in size population, whereas the moment density method accounts for local diversity in bubble hydrodynamics thanks to a singleaveraged contribution of interphase transfers. These two promising refinements of the twofluid model for bubbly flows have shown their ability to recover consistent description of the lift force in the upper part of the flow, where complex flow structures are observed.
While the closure models on bubble forces, which are responsible for the simulation of bubble migration, allow to explain the bubbles hydrodynamic behavior observed experimentally, clear deviations occur for bubble coalescence and breakup. Both methods based on similar coalescence and breakup models lead to the same conclusion: in the simulations of the TOPFLOW experiments, the Luo and Svendsen model leads to an overestimation of the breakup that appears as negligible in the experiments. On the other hand, the coalescence model of Prince and Blanch seems able to recover the correct bubble size, if used solely, as attested by corresponding Neptune_CFD calculations. For both methods, the presently applied models describing bubble breakup and coalescence could be proven as weak points in numerous CFD analyses. These bubble breakup and coalescence models depend to a large extent on the turbulence properties of the twophase flow, which were not measured and could not be validated in the pipe flow test cases. Therefore, further investigations are necessary to determine whether the currently used multiphase flow turbulence models deliver appropriate and verifiable quantities that can be used for the description of bubble dynamics processes.
Extensions of both the moment density method and the MUSIG method to nucleate boiling regime numerical simulation are in progress. This includes the phenomena of compressibility, phasechange, and wall nucleation. To model the bubble sizedependent lateral migration phenomenon, the moment density method should also include a model for the bubble velocity distribution. This can be done using a similar formalism to the present model for bubble size distribution function.
Nomenclature:  Specific interfacial area [m^{−1}] 
:  Lift force coefficient [] 
:  Bubble diameter [m] 
:  Lagrangian derivative 
:  Eötvös number [] 
:  Size distribution function [m^{−4}] 
:  Breakup and coalescence related variations of [m^{−4}s^{−1}] 
Breakup, coalescence coefficients []  
:  Lift force [kg m s^{−2}] 
:  Superficial velocity [m s^{−1}] 
:  Bubbles number density [m^{−3}] 
:  Number of velocity groups [] 
:  Number of subsize groups [] 
:  Reynolds number [] 
:  Velocity [m s^{−1}] 
:  Bubble velocity [m s^{−1}] 
:  Axial coordinate [m] 
:  Liquid 
:  Gas 
:  Volumetric fraction [] 
:  Density [kg m^{−3}]. 
Acknowledgments
Part of this study is carried out at the Institute of Safety Research of the FZD as a part of current research projects funded by the German Federal Ministry of Economics and Labour, Project nos. 150 1265 and 150 1271. The other part of this study carried out at the Institut de Radioprotection et de Sûreté Nucléaire (IRSN) has been achieved in the framework of the Neptune project, financially supported by CEA (Commissariat à l'Énergie Atomique), EDF (Électricité de France), IRSN, and AREVANP. The authors express their gratitude to the technical TOPFLOW (FZD) team.
References
 A. D. Burns, T. Frank, I. Hamill, and J.M. Shi, “The favre averaged drag model for turbulent dispersion in Eulerian multiphase flows,” in Proceedings of the 5th International Conference on Multiphase Flow (ICMF '04), pp. 1–17, Yokohama, Japan, MayJune 2004, paper no. 392. View at: Google Scholar
 E. Deutsch and O. Simonin, “Large eddy simulation applied to the motion of particles in stationary homogeneous fluid turbulence,” in Turbulence Modification in Multiphase Flows, vol. 110, pp. 35–42, ASME FED, Portland, Ore, USA, 1991. View at: Google Scholar
 E. A. Ervin and G. Tryggvason, “The rise of bubbles in a vertical shear flow,” Journal of Fluids Engineering, vol. 119, no. 2, pp. 443–449, 1997. View at: Publisher Site  Google Scholar
 R. O. Fox, “Introduction and fundamentals of modelling approaches for polydisperse multiphase flows,” in Multiphase Reacting Flows: Modelling and Simulation, D. L. Marchisio and R. O. Fox, Eds., vol. 492 of CISM Courses and Lectures, pp. 1–40, Springer, New York, NY, USA, 2007. View at: Publisher Site  Google Scholar
 T. Frank, P. J. Zwart, J.M. Shi, E. Krepper, D. Lucas, and U. Rohde, “Inhomogeneous MUSIG model—a population balance approach for polydispersed bubbly flows,” in Proceedings of the International Conference on Nuclear Energy for New Europe, Bled, Slovenia, September 2005. View at: Google Scholar
 Th. Frank, P. J. Zwart, E. Krepper, H.M. Prasser, and D. Lucas, “Validation of CFD models for mono and polydisperse airwater twophase flows in pipes,” in OECD/NEA International Workshop on The Benchmarking of CFD Codes for Application to Nuclear Reactor Safety (CFD4NRS '06), Garching, Germany, September 2006, paper no. B632. View at: Google Scholar
 Th. Frank, H.M. Prasser, M. Beyer, and S. Al Issa, “Gasliquid flow around an obstacle in a vertical pipe—CFD simulation and comparison to experimental data,” in Proceedings of the 6th International Conference on Multiphase Flow (ICMF '06), Leipzig, Germany, July 2007, paper no. 135. View at: Google Scholar
 A. Guelfi, D. Bestion, M. Boucker et al., “NEPTUNE: a new software platform for advanced nuclear thermal hydraulics,” Nuclear Science and Engineering, vol. 156, no. 3, pp. 281–324, 2007. View at: Google Scholar
 D. P. Hill, The computer simulation of dispersed twophase flows, Ph.D. thesis, Imperial College of Science, Technology and Medecine, University of London, London, UK, 1998.
 A. M. Kamp, A. K. Chesters, C. Colin, and J. Fabre, “Bubble coalescence in turbulent flows: a mechanistic model for turbulenceinduced coalescene applied to microgravity bubbly pipe flow,” International Journal of Multiphase Flow, vol. 27, no. 8, pp. 1363–1396, 2001. View at: Publisher Site  Google Scholar
 E. Krepper, Th. Frank, D. Lucas, H.M. Prasser, and P. J. Zwart, “Inhomogeneous MUSIG modela population balance approach for polydispersed bubbly flows,” in Proceedings of the 6th International Conference on Multiphase Flow (ICMF '06), Leipzig, Germany, July 2007, paper no. 375. View at: Google Scholar
 E. Krepper, M. Beyer, Th. Frank, D. Lucas, and H.M. Prasser, “Application of a population balance approach for polydispersed bubbly flows,” in Proceedings of the 6th International Conference on Multiphase Flow (ICMF '06), Leipzig, Germany, July 2007, paper no. 378. View at: Google Scholar
 E. Krepper, D. Lucas, and H.M. Prasser, “On the modelling of bubbly flow in vertical pipes,” Nuclear Engineering and Design, vol. 235, no. 5, pp. 597–611, 2005. View at: Publisher Site  Google Scholar
 S. Lo, June 1996, Application of the MUSIG model to bubbly flows, AEAT1096, AEA Technology.
 H. Luo and H. F. Svendsen, “Theoretical model for drop and bubble breakup in turbulent flows,” AIChE Journal, vol. 42, no. 5, pp. 1225–1233, 1996. View at: Publisher Site  Google Scholar
 H.M. Prasser, A. Böttger, and J. Zschau, “A new electrodemesh tomograph for gasliquid flows,” Flow Measurement and Instrumentation, vol. 9, no. 2, pp. 111–119, 1998. View at: Publisher Site  Google Scholar
 H.M. Prasser, E. Krepper, and D. Lucas, “Evolution of the twophase flow in a vertical tube—decomposition of gas fraction profiles according to bubble size classes using wiremesh sensors,” International Journal of Thermal Sciences, vol. 41, no. 1, pp. 17–28, 2002. View at: Publisher Site  Google Scholar
 H.M. Prasser, T. Frank, M. Beyer, H. Carl, H. Pietruske, and P. Schütz, “Gasliquid flow around an obstacle in a vertical pipe experiments and CFD simulation,” in Proceedings of the Annual Meeting on Nuclear Technology, vol. 16, Aachen, Germany, May 2005. View at: Google Scholar
 H.M. Prasser, M. Beyer, H. Carl et al., “Evolution of the structure of a gasliquid twophase flow in a large vertical pipe,” Nuclear Engineering and Design, vol. 237, no. 15–17, pp. 1848–1861, 2007. View at: Publisher Site  Google Scholar
 M. J. Prince and H. W. Blanch, “Bubble coalescence and breakup in airsparged bubble columns,” AIChE Journal, vol. 36, no. 10, pp. 1485–1499, 1990. View at: Publisher Site  Google Scholar
 P. Ruyer, N. Seiler, M. Beyer, and F. P. Weiss, “A bubble size distribution model for the numerical simulation of bubbly flows,” in Proceedings of the 6th International Conference on Multiphase Flow (ICMF '06), Leipzig, Germany, July 2007, paper no. 484. View at: Google Scholar
 J.M. Shi, P.J. Zwart, Th. Frank, U. Rohde, and H.M. Prasser, “Development of a multiple velocity multiple size group model for polydispersed multiphase flows,” 2004, Annual Report of Institute of Safety Research. Forschungszentrum Rossendorf, Dresden, Germany. View at: Google Scholar
 O. Simonin, “Continuum modelling of dispersed twophase flows,” in Combustion and Turbulence in Two Phase Flows, VKI Lecture Series 199602, Von Karman Institute for Fluid Dynamics, Brussels, Belgium, 1996. View at: Google Scholar
 A. Tomiyama, I. Sou, I. Zun, N. Kanami, and T. Sakaguchi, “Effects of Eötvös number and dimensionless liquid volumetric flux on lateral motion of a bubble in a laminar duct flow,” in Advances in Multiphase Flow, pp. 3–15, Elsevier Science, Amersterdam, The Netherlands, 1995. View at: Google Scholar
 A. Tomiyama, “Struggle with computational bubble dynamics,” in Proceedings of the 3rd International Conference on Multiphase Flow (ICMF '98), pp. 1–18, Lyon, France, June1998. View at: Google Scholar
 R. M. Wellek, A. K. Agrawal, and A. H. P. Skelland, “Shapes of liquid drops moving in liquid media,” AIChE Journal, vol. 12, no. 5, pp. 854–862, 1966. View at: Publisher Site  Google Scholar
 P. A. Zwart, A. Burns, and C. Montavon, “Multiple size group models,” Tech. Rep. CFX5.7, AEA Technology plc, Oxfordshire, UK, November 2003. View at: Google Scholar
Copyright
Copyright © 2009 E. Krepper et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.