Bubble motion and reaction in different viscous liquids

Reactive bubble columns are omnipresent in the chemical industry. The layout of these columns is still limited by correlations and therefore improved simulation techniques are required to describe the complex hydrodynamics/reaction interaction. In this work, we focus on the numerical and experimental study of the viscosity influence on bubble motion and reaction using an Euler-Lagrange framework with an added oscillation and reaction model to bring the column layout base closer to a predictive level. For comparison and validation, experimental data in various water-glycerol solutions was obtained in a cylindrical bubble column at low gas hold-up, where the main parameters such as bubble size, motion, and velocities were detected. Glycerol leads thereby to a change in viscosity and surface tension. Further, the surface tension was modified by addition of a surfactant. The bubble oscillating motion in low to higher viscosity could be described using an Euler-Lagrange framework and enables a description of industrial bubble flows. In addition, the simulations were in good agreement concerning reactive mass transfer investigations at higher viscosity of the liquid which led to an overall lower mass transfer compared to the cases with lower viscosity.


Introduction
Reactive bubble columns are widespread in chemical process engineering. The overall performance of these reactors relies on the interaction between the hydrodynamics, interfacial mass transfer, and chemical reactions. The layout of these columns is commonly based on simplified integral models that are not able to track the complex interactions between the local hydrodynamics and the reactions. Hence, a detailed knowledge about the ongoing interactions is required. For the simulation of bubbly flows, various approaches are used, ranging from the detailed resolution of individual bubbles to the prediction of industrial reactors with millions of bubbles inside. Single bubble investigations mainly include the resolution of the interface, shape deformation, and bubble oscillations that are directly simulated (Lörstad and Fuchs, 2004;Dijkhuizen et al., 2010;Pesci et al., 2018). The coarse scale Euler-Euler (EE) approach treats both phases as interpenetrating continua. The particle size distribution can be accounted for population balance modelling (Wang et al., 2006;Bhole et al., 2008;Petitti et al., 2010;Renze et al., 2014;Hlawitschka et al., 2016Hlawitschka et al., , 2017Schäfer et al., 2019). The interactions between the phases are modelled in this case. Further implementations to account for reactions were made, that enabled a prediction of concentration profiles in bubble columns (Renze et al., 2014;Gruber et al., 2015;Hlawitschka et al., 2017). The following summary highlights the research based on the Euler-Langrange (EL) approach. The EL model describes bubbles as point volumes (particles) acting under Newtonian dynamics. The forces of the surrounding fluid and neighboring bubbles determine the path of these Lagrangian particles through the domain. The continuous phase itself is calculated using the Navier-Stokes equations and is coupled to the interaction forces by source terms. Each bubble follows a different path where all forces are calculated for each bubble individually. This is one of the main advantages compared to EE methods. Sokolichin et al. (1997) investigated different numerical solution techniques combined with EL and EE simulations and could predict the flow patterns in bubble columns for a long period of time. Lain et al. (1999) included bubble wake-generated Vol. 4, No. 1, 2022, 26-38 Experimental andComputational Multiphase Flow https://doi.org/10.1007/s42757-020-0072-4 Nomenclature a bubble short half axis (m)   (Deen et al, 2004). Darmana et al. (2007) enabled a first detailed modelling of the hydrodynamics and reactive mass transfer in bubble columns and reached a good trend regarding mass transfer compared to experimental results. Gong et al. (2007) applied the EL mass transfer approach for the case of ozone dissolution in a 3-D water filled tank in regard to obtain an optimum bubble size. The transition of flow regime was studied by Simonnet et al. (2008). By including the gas hold-up in the drag correlation, the transition from homogeneous to heterogeneous regime could be well described. Later, Lau et al. (2011) applied a multiscale approach to determine an optimum drag correlation for EL bubble swarm simulations. Shams et al. (2011) focused on the cavitation of bubbles and the volume change of bubbles.
The bubble volume oscillation effects on the liquid phase could be described based on the Rayleigh-Plesset equation. Gruber et al. (2013Gruber et al. ( , 2015 stated that monodisperse bubbles are sufficient to correctly predict the flow features, while breakage and coalescence are required to predict the correct specific interfacial area. The effects of coalescence and breakage decreased compared to the shrinkage of bubbles due to fast reactions. Lau et al. (2014) Sommerfeld et al. (2018) highlighted the influence of tumpling bubble motion and concluded its relevance to mass transfer, while a direct comparision to mass transfer experiments was not presented. However, in the simulation of bubble columns, many of these models do not describe in detail the motion of individual bubbles, e.g., oscillating or tumbling motion of the bubbles (Fig. 1), which is described by Aybers and Tapucu (1969) based on experiments. The motion of the bubbles is dependent on the Reynolds number and the individual bubble motion ranges from rectilinear motion, motion on a helical path, first plane then helical motion, plane motion to rectilinear motion with rocking. A first numerical description for the complex rising behavior was found by Mougin and Magnaudet (2006). For this reason, a rotation equation has been introduced to the EL framework-namely in the premade solver from OpenFOAM "sprayFOAM" (Weber, 2018). The solver was recently further extended to account for the reactions taking place in the liquid phase, as well as to account for bubble coalescence and breakage. In particular, by modelling the shape of the bubbles and the resulting spiraling motion, a more detailed mapping of the bubble flow was made possible. From an experimental point of view, some data about oscillative motion of single bubbles exists in literature. The values presented in the following extracted from the pictures of oscillative motion from literature, while the diameter of the bubble was directly given by the authors. It may be mentioned here that the oscillative motion may depend on the specific experimental set-up, such as flow rate through the distributor, distributor type, and its geometry. Aybers and Tapucu (1969) were among the first to describe the bubble motion in liquids itself, focusing on the bubble velocity and general motion. In their specific case, the wavelength of a single rising bubble with a diameter of 1.76 mm was approximately 70 mm. The path and oscillations of single bubbles were investigated for millimeter sized bubbles in water and polymer solutions by Shew and Pinton (2006). The oscillation wavelength was around 70 mm and the amplitude could be approximated to 2.2 mm for a 2.24 mm bubble. Liu et al. (2015) analyzed the bubble dynamics in stagnant water and glycerol solutions. The respective measured amplitude in horizontal direction was determined to 2.4 mm for a bubble of 2.01 mm diameter. The corresponding wavelength could be approximated to 50 mm. Sharaf et al. (2017) investigated the shapes and paths of single bubbles rising in different glycerol solutions. For a 25% water-glycerol solution, the bubble amplitude could be approximated from the experiments based on the given bubble diameter to approximately 40 mm. Experiments in bubble columns using high viscous liquids were either made in plug flow setups (e.g., de Azevedo et al., 2017) or at higher gas hold-up, where in the second case, bubble interactions such as coalescence may have a dominant effect on bubble size and with that on the hydrodynamics (Wilkinson and Dierendonck, 1990;Urseanu et al., 2003;Kulkarni and Joshi, 2005;Kajero et al., 2012;Besagni et al., 2017;Laupsien et al., 2017). Kulkarni and Joshi (2005) described no effect or only a slight effect of viscosity on bubble size. Better measurement techniques and algorithms however improved the sensitivity of the analysis and led to better results in recent years (Lichti and Bart, 2017).
In this work, we extended previous numerical work on bubble oscillation using an Euler-Lagrange code (Weber, 2018) to study the motion in higher viscous liquids, at different surface tensions, and to investigate the applicability of the aforementioned model to new conditions. To the authors' knowledge, the effect of viscosity on bubble oscillative motion has not been studied before using Euler-Langrange modelling. In addition, mass transfer at higher viscosity is analyzed in this study. Therefore, experimental data from a cylindrical bubble column using different glycerol/water and sodiumdodecylsulfat (SDS) concentrations is used. The relevant parameters, such as bubble oscillation amplitude and wavelength, bubble velocity, and CO 2 concentration fields are analysed.

Numerical methods
In the following description, the subscripts c and b indicate the continuous phase and the bubbles, respectively. The implementations are based on the solver "sprayFoam" in OpenFOAM version 2.4.0.

Liquid phase hydrodynamics
The continuous phase is treated as incompressible fluid. The basis for the calculation is a modified Navier-Stokes equation: The continuous velocity vector is named c u . Further, the density is given by , ρ the pressure by , p the viscosity by , μ and the source term, e.g., responsible to couple the phases, is named f. The turbulence in the liquid phase is modelled based on the standard RANS k- model (Launder and Spalding, 1972). In this model, the turbulence is based on the kinetic energy k and the energy dissipation . Turbulent eddies are not resolved directly by the model, but the total amount of kinetic energy is tracked instead. Both scalar fields (k and ) are obtained by the following transport equations: where S ij is the mean strain rate. The turbulent eddy viscosity, t , μ is based on the relation of k 2 to : The influence of bubbles to turbulence is accounted by a bubble induced turbulence model presented by Rzehak and Krepper (2013). The source term of bubble induced turbulent kinetic energy is with F D as the drag force and c b u u as the relative velocity between liquid and bubble. The turbulent dissipation source term ε S includes the time scale τ given as where d b is the bubble diameter and k the turbulent kinetic energy. Based on the model constant

Disperse phase hydrodynamics
The disperse phase is treated in the Lagrangian framework.
The simulated particles are supposed to belong to the cell, where the center of mass is located. The movement of particles to a new position x i+1 is straightforward. The calculation of the particles velocities is based on an implicit Euler algorithm.
The acceleration a i is calculated using the (effective) particle mass and the sum of all forces acting on the bubble.

Bubble hydrodynamics
The movement of bubbles is determined on a number of different acting forces: The acting forces include the buoyancy and weight force F B , the drag force F D , the lift force F L , and the virtual mass force F VM , as well as the wall lubrication force F W . Further models, such as the bubble dispersion force F BD can be included. The forces include the gravitational acceleration g the bubble's volume V, and the turbulent kinetic energy k and the phase fraction . The model parameters are indicated by C i : denotes the material derivative. The constants used in the Eqs. (12)-(16) are taken from literature as a basis. It will be shown later, that sufficiently good result could be obtained with the applied constants. For an industrial layout of a bubble column, the values could be further validated by single bubble experiments using the respective system (Hlawitschka et al, 2014) or Direct Numerical Simulations (DNS) (Mühlbauer et al., 2019). The lift and drag model constants C L and C D are based on the correlations of Tomiyama (2002Tomiyama ( , 2004. Tomiyama noted that the correlations for C L accidentaly yield to a good prediction in air-water systems while it was developed experimentally for higher viscous systems, such as glycerolwater (Tomiyama, 2002). The drag coefficient was also tested for a wide range of aspect ratios in distilled water, while it was stated that the aspect ratio is one of the dominant parameters on velocity and not the surfactant itself. Additionally, the swarm influence presented by Rusche (2002) is applied, but has a minor effect there the average hold-up is low (< 1%). The virtual mass coefficient C VM is set to 0.5 (Delnoij et al., 1997) and therefore assumes an independency of the void fraction. Overall, this value led to suffiecient good results in a variety of publications. The wall force coefficient C W is modelled based on the correlation of Antal et al. (1991) and the turbulent dispersion force coefficient C BD is 0.1 accordingly to Lahey Jr. et al. (1993). Additionally, the turbulent dispersion accounts for the collision of bubbles and turbulent eddies. The random dispersion model that is compared to the oscillation model in this work is based on the derivation of Smith and Milelli (1998). In the random dispersion model, assuming an isotropic turbulence, eddies are traveling through the liquid in a uniformly random direction with a specific lifetime. In the model, the turbulent eddy lifetime t E is given, after which a new eddy is calculated as The movement direction of the eddy is uniformly chosen. The velocity follows a normal distribution with a variance dependent on the turbulent kinetic energy k : This turbulent velocity is accounted by the underlying continuous phase velocity for each bubble individually and is applied to calculate the different bubble forces. The bubble interactions, such as coalescence, can be accounted for by stochastic algorithms, such as the stochastic collision equation by O'Rourke (1981): A collision of two bubbles leads to a further investigation of the coalescence probability coal P , that relates the film rupture time t r to the t c contact time (Coulaloglou and Tavlarides, 1977). The film rupture time calculation is based on the model of Chappelear (1961). The breakup model is based on the impact of turbulent eddies (Coulaloglou and Tavlarides, 1977). Both coalescence and breakage play a minor role in this investigation due to a low hold-up inside the bubble column and were barely observed in the experiments.

Ellipsoidal bubble model
The bubble changes its shape in pure water from spherical to ellipsoidal to spherical cap shape depending on the size of the bubble and liquid properties. In many applications, spherical bubbles and ellipsoidal bubbles are observed. Therefore, the bubble can be approximated by an oblate spheroid, representing an ellipsoid with two different axes for height and width (Fig. 2). The described implementation of motion is based on the Jeffery's equation (Jeffery, 1922;Junk and Illner, 2007).
The ratio of the axes describes the shape factor (eccentricity) of the bubble: A shape factor of unity describes a perfect sphere. The rotation of the bubbles is accounted for by an additional force, represented as a modified drag term. The oscillation of a real bubble only occurs in a certain range (Aybers and Tapucu, 1969), the force only acts in the range of 500 <Re b < 1300.
[ ] s D S , 500, 1300 The side force is based on the drag force and the direction vector p. The bubble path amplitude needs to be calibrated by the factor s β based on experiments, where higher values imply a larger amplitude. The bubble rotation is calculated based on an explicit Euler algorithm: The orientation  p of the bubble relaxes against the outer rotation G with the correction factor . γ The additional term p T includes a random rotation due to turbulence, where R s is applied to scale the effect of randomness. Both factors have to be applied based on experiments. The outer rotation is described as The first term in the outer rotation equation describes the effect of the surrounding liquid transferred to the bubble and can be adjusted by J 1 , the second term is the rotation by the shear stress in the surrounding fluid and is influenced by J 2 and the third term is responsible to describe the ambition of a bubble to orient along the gravitational direction and is dependent on J 3 . The constants were introduced by Weber (2018) to adjust the theoretical model of ellipsoidal motion to experimental data.

Mass transfer
The species mass transport in the continuous fluid described by an additional scalar transport equation: The transferred mass per cell volume between the phases is accounted by i R : The transferred mass is calculated based on the two film theory, where a homogeneous and constant species mass fraction is considered in the bubbles.
The equilibrium concentration ,interface i Y at the gas-liquid interface is described by Henry's law:

Chemical reaction
The investigated chemical reaction is the chemisorption of carbon dioxide CO 2 in sodium hydroxide solution NaOH (aq) . The first step is the absorption of 2 CO from the gas phase into the liquid phase: The equation was simplified to account for a constant reaction coefficient, by neglecting temperature and activation energy influence. The applied constant for the first forward reaction is A = 9162, and for the first backward reaction A = 2.4×10 -4 . The second reaction constants are 10 6 and 168 for the forward and backward reaction, respectively. The mass transfer enhancement describes the increased mass transfer due to reaction. The enhancement factor is close to unity for the experiments having a lower pH value of 10. The general equation to describe the enhancement factor also for higher pH values is used in this work for the simulation (Hlawitschka et al., 2017): For smaller pH values, the enhancement factor equals 1. The previous described source term R j in the species transport equation is hence extended by the mass transfer due to reaction: Herein, stoichiometric factors are given by ik ν and the reaction rate is described by . k ω The subscript i indicates a certain species and k denotes a certain reaction. The total number of reactions observed is labeled by K. Finally, the reaction rate depends on the forward and backward reactions and the species concentrations that were previously described.

Experimental setup
The bubble column used in the experiments has a diameter of 0.142 m and the dimensions are summarized in Table 1. The experimental bubble column is equipped with a four-needle sparger (Fig. 3). The needles are aligned in-line in the centre of the column with a spacing of 22 mm. The nozzles themselves are made of stainless steel with an inner diameter of 0.25 mm. Around the bubble column, a rectangular box is installed to enable refractive index matching for the optical measurement techniques, used in the experiments. The bubble column is filled to a level of 730 mm with the corresponding liquid. As liquids, deionized water is used with amounts of glycerol ranging from 0 to 50 wt%. In addition, the surfactant sodiumdodecylsulfat (SDS) was applied in concentrations up to 6 mg/g liquid, leading to 9 different solutions. The CO2 gas volumetric flow rate was set constantly to 7.5 h -1 air equivalent volume flow. All experiments were executed at atmospheric pressure and room temperature, approximately 20 °C. The hydrodynamics were analyzed using particle image velocimetry (PIV; liquid phase) at 3.3 Hz and particle tracking velocimetry (PTV; bubbles) at an acquisition frequency of 100 Hz. A 2-tracer laser induced fluorescence (2T-LIF) set-up was applied to study the concentration changes in the sodium hydroxide solution with an initial pH value of 9.5. The first tracer, uranine, is thereby used to detect the concentration changes; the second tracer is applied to reduce the effects of reflections of the laser at the bubble surface, inhomogeneities in the laser sheet and bubble shadows. For the experimental determination  of the pH development due to reaction, the measurement technique is limited in the range of pH 5.5 to 8.5, in which the tracer uranine is active. All details of the experimental setup, the measurement systems and the measurement results can be found in Kováts et al. (2020).

Simulation
The dimensions of the column, such as diameter and height, correspond to the experimental investigation. The inlets are represented by four point sources in the center of the column ( Table 2). The previously described numerical framework in Section 2 is applied as baseline model. However, e.g., the wall influence to the bubbles becomes negligible small in most parts of the column, since the bubbles are located mainly in the center of the bubble column. As boundary conditions, the walls are treated by a no-slip condition for the continuous phase. The top of the column is defined as slip condition for the continuous phase, there only a slight increase of the water level will be observed due to the low gas hold-up. At the outlet, atmospheric pressure condition is used, while the walls have a zero gradient pressure boundary condition. The simulation is started for 100 s without reaction to obtain a premixed condition in the column. This velocity field is used as starting point for the reaction simulations, while no bubbles, are initially inside the column for these simulations. For the calculations, the oscillation factors 1 2 , , J J and 3 J were taken from previous investigations (Weber, 2018), where they were adapted for the CO 2 /water system to 0.05, 0.05, and 2, respectively. The correction factors s , β , γ and R s were set to 0.15, 650, and 3, accordingly. The parameters were calibrated based on experimental data of gas bubbles in water using the same experimental setup (Weber, 2018). Therefore, the optimization software Dakota (https://dakota. sandia.gov/) was used to find the best parameter set on a variety of parameter combinations. As target function, the characteristic bubble motion parameters, frequency and amplitude were used.

Mesh study
Different mesh sizes were investigated to study convergence of the newly adapted bubble rotation scheme with respect to bubble oscillation and amplitude for the case of pure water. The average bubble amplitude in the experiment (Fig. 4) is 2.1 mm, and the wavelength is 45 mm at a gas volumetric flow of 7.5 h -1 air equivalent. The average amplitude and wavelength of the three investigated meshes are depicted in Table 3. It can be seen that only slight changes in these values can be observed and the implementation of the bubble oscillation is not directly dependent on the mesh resolution. The finally applied mesh is based on 140,100 hexahedra cells leading to an approximated minimum cell edge length of 3.9 mm, slightly larger than the expected maximum bubble size and is depicted in Fig. 5. This mesh enables a good representation of the velocity field close to the wall. In addition, the concentration fields in the liquid phase are averaged with the cell volume and require a fine mesh to be comparable to the experiments. The simulation time is acceptable with a simulation time of approximately 12 h for 200 s of real time simulation using a time step size of 0.001 s and a parallelization at eight cores.

Results and discussion
The glycerol concentration was varied in three steps: zero concentration, 25% glycerol, and 50% glycerol. The SDS concentration is 0 mg/g (T0), 2 mg/g (T1), and 6 mg/g (T2). The liquid properties, measured for these solutions, are listed in Table 4. The bubble size and shape is influenced by the system properties such as viscosity and surface tension. In this study, the boundary conditions for the simulations were set to the experimentally determined inlet values. During bubble rise, the size is mainly influenced by the reaction (mass transfer) and pressure change. In Table 5, the average spherical equivalent bubble sizes measured at a column height of 50 mm are depicted along with the resulting average vertical bubble velocity at steady state conditions (no reaction). More than 100,000 bubbles were analyzed for each case. The standard deviation for each case was  In the case of higher viscosity (50% glycerol), the individual bubbles started to follow a straight path close to the inlet in the experiments and mixed further upwards the column. This has to be accounted for an adequate modelling in the simulation. To show the effect of the bubble oscillation from different implementations, three different cases were analyzed: The first case does not include random dispersion or the newly implemented oscillation scheme. In the second case, random dispersion is enabled, while in the third case, the newly implemented scheme is enabled. The results are compared in Fig. 6 for the case of pure water. Without random dispersion or additional oscillation term (Fig. 6, left), the bubbles follow the liquid flow field, rising in bubble chains that do not interact further. The random dispersion approach (Fig. 6, centre) leads to an earlier mixing of the bubbles close to the nozzles at the bottom. The individual bubble chains do intermix after a short distance. In the case of the implemented oscillation scheme (Fig. 6, right), the bubbles start to oscillate close to the disperser (bottom) and start to intermix.
The resulting bubble movement from experiment and simulation is compared in Fig. 7 as an overlay. The figure depicts the lowest 30 cm of the column. The simulated bubble is shown as a sphere in green, and the bubble obtained from experiment is given in black. It can be seen that a higher glycerol concentration and a higher SDS concentration leads to a decrease in bubble size oscillation close to the inlet. In the case of 0% glycerol and 6 mg/g SDS (T2), the transition from straight bubble rise to the helical bubble movement is observed at a height of approximately 0.018 m (Fig. 7, red line).  At a glycerol concentration of 25% and the same amount of SDS, the transition starts at a height of 0.03 m and at a glycerol concentration of 50%, the transition is at 0.18 m. The simulation predicts a slightly higher amplitude of the oscillations at high glycerol concentrations compared to the experiment, while the transition to helical motion is still described rather correctly. Overall, a good agreement could be found concerning bubble oscillation.
The time-averaged vertical velocities of the liquid phase from the PIV experiment and simulation in the case of 25% glycerol are depicted in Fig. 8 and Fig. 9 at a column height  of 210 and 410 mm, respectively. The experimental vertical velocity reaches a value of 0.56 m/s and in the simulation, the velocity reaches an average value of 0.57 m/s in the bubble region at a height of 210 mm. In total, the velocity peaks in the simulation (Fig. 8) are more pronounced than the velocity profile obtained from the experiments and may be a result of bubble passage. Around the bubbles, approximately one millimeter, the velocity cannot be detected by the measurement technique; but here, the highest velocities are expected. At a column height of 410 mm, the maximum velocity is 0.41 m/s in the experiment and 0.40 m/s in the simulation.
The influence of liquid viscosity on the vertical velocity of the liquid phase is shown in Fig. 10. At a column height of 410 mm, the liquid velocity increases with increasing glycerol concentration (0 mg/g surfactant), which was also found in the experiments (symbols in Fig. 10).
The influence of the glycerol concentration on the reaction development is investigated in Fig. 11. It can be seen, that the lowest viscosity (0% glycerol, 0 mg/g SDS) leads to the fastest decrease of the pH value. Higher glycerol concentration leads to a slower mass transfer and pH value decrease. A good agreement between the simulation (lines in Fig. 11) and experiment (symbols in Fig. 11) could be obtained.

Conclusions
The layout of reactive bubble columns is still challenging due to the complex interaction of hydrodynamics and mass transfer. Thereby, the motion of the bubble, especially in Fig. 10 Comparison of the liquid vertical velocity at a glycerol concentration of 0% (blue), 25% (green), and 50% (red), with 0 mg/g surfactant and at a liquid height of 410 mm from simulation (lines) and experiment (dots).

Fig. 11
pH-value development for the cases of 0% glycerol (blue), 25% glycerol (green), and 50% glycerol (orange) and 0 mg/g SDS at a column height of 370 mm. Simulation results are presented as lines, experiments as symbols.
higher viscosity liquids was mainly neglected in the numerics. A numerical framework to describe individual bubble motion, including oscillation, is used to our knowledge for the first time to study a bubble column with an Euler-Lagrange code at low to higher viscosity including reactive mass transfer. In a first step, the extension of bubble rotation is validated with experimental results at different liquid viscosities. It could be shown that the extended solver is able to describe the complex transition from straight bubble movement to helical bubble movement and therefore improves the description of the bubbles' spatial distribution in columns. At a flow rate of 7.5 h -1 air equivalent and in pure water, the measured wavelength of bubble oscillation is in average 45 mm and the amplitude 2.1 mm. The distance of transition from straight motion at bubble generation, to oscillating bubble motion, is longer with higher glycerol and SDS concentration. With the new CFD framework, these measured parameters, as well as velocity fields and bubble positions can be predicted with good agreement for the different viscosity conditions. An adaption of the main models, such as drag, lift, and oscillation model was not necessary to account for the higher viscosity and could be calculated straight forward. In a second step, reactive mass transfer investigations at higher viscosity led to an overall lower mass transfer compared to the cases with lower viscosity. The developed framework enables therefore a good prediction of the hydrodynamics and mass transfer and enables an industrial application for the layout of bubble columns.