#### Abstract

Bubble Sauter Mean Diameter (SMD) in gas-liquid multiphase system is of particular interest and the quantification of gas characteristics is still a challenge today. In this contribution, multiphase Computational Fluid Dynamic (CFD) simulations are combined with Population Balance Model (PBM) to investigate the bubble SMD in baffled stirred tank reactor (STR). Hereby, special attention is given to the phenomenon known as the fluid macro-instability (MI), which is a large-scale low-frequency fluid velocity variation in baffled STRs, since the fluid MIs have a dominating influence on the bubble breakage and coalescence processes. The simulations, regarding the fluid velocity, are validated with Laser Doppler Anemometry (LDA) experiments, in which the instant radial velocity is analyzed through Fast Fourier Transform (FFT) spectrum. The frequency peaks of the fluid MIs are found both in the simulation and in the experiment with a high degree of accuracy. After the validation, quantitative predictions of overall bubble SMD with and without MIs are carried out. Due to the accurate prediction of the fluid field, the influence of the fluid MI to bubble SMD is presented. This result provides more adequate information for engineers working in the field of estimating bubble SMDs in baffled STRs.

#### 1. Introduction

In multiphase mixing processes of the chemical and biotechnological industry, massive means of agitation are required. Among the various types of mixing tanks, the stirred tank reactors (STRs) with different kinds of impellers are the most widely process units. Therefore, many investigations have been carried out in the past to predict the fluid hydrodynamics which have a fundamental influence on the overall mixing and multiphase performance of the STRs.

One of these particular phenomena of interest is the fluid macro-instability in baffled STRs. Several experimental results have shown the existence of a large-scale low frequency fluid velocity with a periodical variation [1, 2]. In currently available literature it is assumed that two types of the fluid MI exist. The first type of the MI is established by the periodical variation of the flow circulation pattern [3–5]. The oscillations of fluid velocity magnitude are driven by turbulent vortex shedding from the interaction between the flow discharge stream and the baffle-body, which leads to pressure reflections of the wall jets. The second kind of MI is generated by the so-called “whirlpool” vortex (or processing vortex) rotating around the impeller shaft with a uniform frequency [6, 7]. Details on conditions of the occurrence of two fundamental types of the MI are given in the comments [8] from Kresta and Roussinova, as well as the comments reply [9] from Yianneskis to the editors. Nevertheless, in such gas-liquid multiphase STRs, less attention has been paid to the influence of fluid MIs on the bubble size distribution, particularly the bubble Sauter Mean Diameter (SMD) in gas-liquid STRs. The fluid MIs resulting in the fluid turbulence condition variation are believed to have the dominating impact on the bubble breakage and coalescence during the bubble SMD prediction. In a baffled gas-liquid STR, the bubble size distribution plays an important role in the hydrodynamics and mass transfer behaviors. It determines the bubble rising velocity and gas residence time and governs the gas holdup and gas-liquid interfacial area. There is insufficient literature providing either the quantitative results or qualitative estimations of the MI influence on the bubbles SMD. The numerical estimation of the influence of the MIs on the bubble size distribution is presented.

In this contribution, in order to obtain a comprehensive understanding of the transient flow character of the fluid MIs and the bubble SMD under the MIs, Population Balance Model (PBM) based Computational Fluid Dynamic (CFD) simulations of a baffled STR with a Rushton turbine are performed. The fluid side of these simulation results is validated with Laser Doppler Anemometry (LDA) experiments. Afterwards, numerical estimations of the influence of the fluid MI type one, which is caused by the vortex shedding of the baffle, on the bubble SMD are carried out. The simulation results with the PBM-CFD modeling could reveal the dynamics behavior of bubbles and might contribute to achieve the optimal design of reactors.

#### 2. Theoretical

##### 2.1. Modeling the Fluid Dynamics in G/L STR

The STR, which is being investigated in this article, operates with a gas injection at the bottom of the impeller shaft. Hence, a multiphase model is required in the CFD simulation, which consists of two phases: a continuous phase (liquid) and a dispersed phase (gas) in the tank. The blades of the Rushton turbine agitate the liquid near the impeller, forming the flow convection circulations in the STR. Thus, a fluid turbulent model and an agitation model need to be taken into account. Additionally, the PBM is integrated with the CFD simulation to describe the bubble breakage/coalescence behavior in the baffled STR. Therefore, during the development of the model, these interactions and requirements must be considered. The Eulerian model is suitable for this case when the gas-liquid slip velocity and the bubble SMD distribution are of interest. Considering the Euler multiphase model, the RANS (realizable* k*-*ε* model) is used to describe the turbulent behavior of our system, since Gimbun et al. [10] reported that the fluid velocity obtained by realizable* k*-*ε* model illustrated the minimum deviation with experiments. So as to estimate the bubble SMD, the PBM is used to predict the bubble size distributions along with fluid dynamic investigations in the multiphase system in STRs [11, 12]. In order to investigate the influence of the MI on SMD, the Moving Reference Frame (MRF) method is applied for the impeller agitation mechanism, since this method can be implemented in both the steady and transient iteration.

There are many breakage and coalescence kernels available for bubbly flows in present literatures. Several researchers [13, 14] made a comparison between different combinations of kernels [15–18]. Their findings suggest that there is no great difference between the mean SMD using the different kernels when the parameters are modified. Laakkonen employed a modified version of Prince and Blanch’s kernel [15] for coalescence kernel and his own breakage kernel in gas-liquid stirred tanks. Their simulation results show a good agreement with experimental measurements [19]. Due to the success in Laakkonen’s results, his breakage and coalescence kernels are used in this contribution. The PBM-CFD modeling structure flow diagram is shown in Figure 1.

##### 2.2. Population Balance Modeling

The general transport equation for the particle number density function in PBM is given as

Herein, is the number density function for certain bubble diameter.is the gas phase velocity in the mesh cell. is defined as the growth rate, which is neglected for our purposes. andare the birth and death rate due to coalescence, whileandare the birth and death rate due to breakage.is the source term.

The birth rate of bubbles with volumedue to coalescence is given bywhereis the coalescence rate and. Here,is the collision frequency andis the coalescence efficiency.

The death rate of bubbles with volumedue to aggregation is displayed by

The birth rate of the bubble with volumedue to breakage is shown by

Here, the breakage rate is expressed as.is the bubble breakage frequency, andis the probability density function (PDF) of the breaking bubble.

The death rate of the bubble of volumedue to breakage is given by

In this work, the bubble coalescence is mainly induced due to fluid turbulence assuming binary coalescence. The modified Prince and Blanch’s coalescence kernel [15] is implemented to express the collision frequencyand the coalescence efficiency:whereis the turbulent dissipation rate.andare the diameters of the two bubbles.is the model constant, which is set to 0.88 for stirred reactors. is the maximum possible gas holdup and equal to 0.8 [20]. is the ratio of the distance between bubbles and the bubble turbulent path length [21].

The coalescence efficiencyis regarded as the liquid drainage as trigger step and the bubbles deforming with immobile surfaces. The constantis usually taken equal to , as proposed by Laakkonen et al. [19]. The bubble breakage is assumed to occur through bubble interactions with turbulent eddies.

The Laakkonen breakage kernel is applied to describe the breakage frequencyand the breakage PDF.

The constants , and are set to be 2.52, 0.04, and 0.01 [19]. With these equations, the population balance of the bubbles in the STR is carried out.

##### 2.3. Computation Environment and Numerical Strategy

Three-dimensional transient and steady simulations are carried out with the commercial CFD software Ansys FLUENT 16.2. The iteration time step is set to be 0.001 s. The residual errors are set to be 10^{−3} for transient case and 10^{−5} for steady case. The initial injected bubble diameter is set to be 1.6 mm, which is calculated with the correlation proposed by Kazakis et al. [22]. The Bin information of the solution method (Class Method, CM) in PBM is listed in Table 1. Approximately grids are applied for the meshing process, which has been found to be adequate regarding the balance between computation accuracy and assumption, since Petitti et al. reported that grids are more than sufficient to predict the mean velocity with the RANS method in an STR with a similar geometry in this study [23]. The PBM-CFD simulations are carried out on a Lenovo ThinkStation with two Xeon E5-2680 processors and NVidia GPU acceleration.

#### 3. Experimental

LDA has proven to be adequate for the measurement of flow fields in stirred tank reactors, since its noncontact measurement has no influence on the original flows in the reactor as well as its high measurement sensitivity due to the optical interference methodology [24]. LDA investigations have identified the formation of the trailing vortex pair behind each impeller blade and provided information concerning the flow periodicity in the impeller vicinity and its associated increase of the turbulence levels [25]. For the case in this contribution, the LDA experiment is capable of validating the CFD simulation extreme case, where no gas injection and thus no bubble breakage/coalescence occur. In other words, only a part of the model can be validated with this experiment. This, nevertheless, is useful, as it shows that the model is valid in a certain area. In latter section of this article, the details on the experimental set-up as well as the experimental results are given.

The schematic diagram of the baffled STR is shown in Figures 2 and 3 and the cylindrical coordinate system is implemented in the radial, axial, and azimuthal coordinates with the origin at the centre of the impeller. This experiment is performed in a cylindrical flat-bottom agitated vessel with a diameter of = 286 mm. The vessel is placed in a glass cube trough with 500 mm in side length. Four symmetrically spaced glass baffles are vertically mounted and separately installed to the tank. Dimethyl sulfoxide is filled into the tank until the height is equal to its diameter. For the stirring, a six-bladed Rushton turbine is used in the tank. The experiments are performed under atmospheric pressure and 20°. In order to minimize the measurement error due to the different reflection ratio between water and normal glass, special DURAN glass with the same reflection ratio 1.473 as the dimethyl sulfoxide is applied for the vessel wall. The fluid physical properties are listed in Table 2. DANTEC LDA Probes are set up for the validation experiment. A burst spectrum analyzer is employed to process the Doppler signals and the time series radial velocity component is measured. The sampling time is over 300 s to establish sufficient data. Hereby, the sampling rate is set to 1 kHz and approximately 300,000 signal points are taken during one experimental run. The anticipation errors of the mean and RMS radial velocity component are around 1.9–6.5% and 8.3–16.7%. The tracer sampling redundancy validation ratio is over 80%.

#### 4. Results and Discussion

The PBM-CFD simulations are carried out for two cases. In the first case, the injected gas flow rate is zero where no bubble coalescence and breakage behaviors are involved and the CFD simulation results are validated by the LDA experiment. In the second case, the injected gas flow rate is 0.065 vvm. Here, on the basis of verified simulation models, the influence of the fluid MI on the bubble SMD is estimated among the transient steady PBM-CFD simulation and multiblock model [19] with the same geometry and boundary conditions.

##### 4.1. CFD Simulation and LDA Validation for MI

For the results of the first case, the fluid velocity vector diagram on the plane (*r/* = −0.5~0,* z/* = −0.5~0.5, and = 45°) with transient simulation is shown in Figure 4, where two flow circulations can be observed above and below the impeller, respectively. In order to investigate the fluid MI, the liquid velocity vector diagrams are presented in time series (from 110 s to 114 s) in Figure 4, where the circulation length changes occurring in the simulations become apparent. The periodical flow circulation pattern exchange can be observed. In the time frame from 110 s to 112 s, the flow circulations above the impeller become smaller, while the other circulations below the impeller become larger. Towards 114 s the flow circulations above the impeller become larger again, while the ones below become smaller. This is shown by the drawn heights , , and in Figure 4, which qualitatively represent the change of the flow pattern circulation. This variation is observable periodically throughout the simulation, which is a clear indication of the first type fluid MI (flow circulation pattern periodical variation).

The fluid radial velocity is measured by the LDA system at the point (*r/* = 0.32,* z/* = 0.38, and °), which is near to the baffle where the vortex shedding effect is strong and remarkable. The raw experimental data is shown as scatter points in Figure 5(a). An obvious magnitude variation and an ambiguous periodical oscillation can be obtained from these signal points. In order to further investigate the fluid MIs, FFT analysis is carried out in the Matlab with numerically resampled frequency of 50 Hz, shown in Figure 5(b). A frequency peak can be captured from the FFT diagram.

**(a)**

**(b)**

Figure 5(a) shows that the mean radial velocity recorded from the LDA experiment is −0.105 m/s at the point (*r/* = 0.32,* z/* = 0.38, and °). Meanwhile, the mean velocity obtained from the CFD simulation is −0.973 m/s. This is a clear indication that the mean radial velocity from the simulation lies in the desired range. However, the amplitude variation in the LDA (from the maximum 0.151 m/s to the minimum −0.353 m/s) is much greater than the amplitude variation in the CFD simulation (from the maximum −0.034 m/s to the minimum −0.146 m/s). Details on this are shown in Figure 6(a). This huge difference is due to measurement noise involved in the LDA experiment. Therefore, Fourier fitting is implemented to wipe out the measurement noise, shown in Figure 6(a), which is the red profile in Figure 5(a). A periodical oscillation can be observed in the fitting diagram. The same procedure is also applied in the CFD simulation. The Fourier fitting parameters are shown in Table 3. Now, the two Fourier fitting profiles from LDA experiment and CFD simulation share a relatively good conformity. This result is displayed in Figures 5(a) and 6(a).

**(a)**

**(b)**

The MI frequency peaks in simulation (0.342 Hz) and in experiment (0.368 Hz) are both captured in the FFT spectrum analysis. In other publications, it has already been discovered that the frequency of the MI in baffled STRs is linearly related to the impeller rotational speed for fixed impeller types [26, 27]. For our case, these similar results are also found in CFD simulations and LDA experiment in Figure 7, which is to validate the CFD models. The slope of the fitting liner is defined as the dimensionless MI frequency, which is independent of the impeller rotational speed. This result can be also verified by the simplified correlation approach from Paglianti et al. [28]. The dimensionless MI frequency in a STR with / = 0.333 and the same Rushton turbine is predicted as 0.084. Our research shows almost identical results (our simulation: 0.088; experiment: 0.082), which is again a similar result to the experiment with the POD method from 0.08 [29]. The good agreement of simulation and experiment results illustrating the CFD models can predict the accurate fluid field, which provides the basis to estimate the influence of the type one MI to the bubble SMD in the baffled STR.

##### 4.2. Influence of MI on Bubble SMD

Moving on, the second part of the CFD simulations is looked at. This part considers gas injection into the STR. Here, no experimental validation is carried out. Nevertheless, the PBM-CFD simulation has been proven to be a reliable method in predicting the bubble SMD in baffled STR [10, 23, 30]. From the results of the previous section we can securely say that the existence of the fluid MIs in baffled STRs has been proven, in both experiments and CFD simulations. Since the fluid velocity is one of the most significant parameters influencing bubble breakage and coalescence kernels, the MI phenomenon must be taken into consideration in the PBM-CFD simulation. Therefore, of interest is the influence of fluid MI on the overall SMD at the 250 rpm rotational speed and gas flow rate = 0.065 vvm in the baffled STR. The comparison is carried out among transient, steady cases and multiblock model [19], so that the influence of the MI on SMD in PBM-CFD simulations is investigated.

The contour of the SMD on plane ( = −0.5~0.5,* z/* = −0.5~0.5,= 45°) in transient and steady cases is demonstrated in Figure 8. The initial bubbles are injected from the porous sparge. The bubbles swirling upwards carry on the breakage and coalescence process due to the fluid turbulence generated by the impeller. The small bubbles flow downwards with the fluid circulation under the impeller and the big bubbles rise to the top with the other circulation above the impeller. The larger bubbles accumulate around the shaft, where the low dissipation rate leads to more coalescence. The smaller bubbles appear near the rotational zone, where the breakage process is dominant due to the high turbulence shear rate. The SMD in various subregion among transient, steady PBM-CFD simulation and multiblock model is shown in Figure 8. The SMD difference is more noticeable among the transient case, the steady case, and the multiblock model in zone 3, zone 8, and zone 10, where the influence of the fluid MI is apparent in transient simulation. The histograms of the overall SMD distribution are shown in Figure 9. The SMD distribution varies slightly at the different time steps, since the bubble breakage and coalescence kernels result in the bubble SMD change with the fluid velocity and turbulence dissipation which is caused by the fluid MIs. Nevertheless, the value of the overall SMD diameter becomes steadily more constant as it is shown in Figure 10. This demonstrates that the breakage and coalescence behavior happening in the overall STR has reached an equilibrium. This also means that the fluid MI does not lead to the overall bubble SMD size oscillation and the equilibrium of the bubble breakage and coalescence in the STR is relatively stable in the final state.

Next, in order to investigate the influence of the fluid MI on the overall SMD, the comparisons among transient PBM-CFD (with MIs), steady PBM-CFD (without MIs), and multiblock model (without MIs) are carried out. The main advantage of the multiblocks model is that the calculation time is quite smaller and it can be implemented in the online operation and optimization. In the multiblock approach, the liquid velocity, the liquid turbulence kinetics, and the liquid dissipation energy rate of the subregions are obtained from the CFD approach via the volume average. Moreover, the same breakage kernel and coalescence kernel are applied in both approaches. Therefore, the three simulations are performed with consistent liquid fluid field and bubble breakage/coalescence kernels. The overall mean SMD in the transient PBM-CFD simulation is 3.12 mm, which is smaller than that in the steady PBM-CFD simulation (3.30 mm) involving the fluid MI. This result is displayed in Figure 10. It is demonstrated here that the equilibrium of bubble breakage and coalescence lies more on the breakage side in the transient PBM-CFD simulation when the fluid MI is considered. Therefore, the fluid MI acts the role of slightly promoting the bubble breakage probability in the STR. The overall SMD in the multiblocks model is 3.41 mm, which is greater than that in the transient PBM-CFD simulation involving the fluid MI. The quantitative result of the influence of the MI on the bubble SMD in the baffled STR is provided in the CFD simulation with the same breakage and coalescence kernels as in the multiblocks model. Thus, for cases in which the multiblocks model is currently used and thus the MI influence is not incorporated, our PBM-CFD method can yield an additional factor to increase the realistic representation of systems.

#### 5. Conclusions

In this contribution a multiphase model is created and combined with a PBM to investigate the MI phenomenon, as one of the most important fluid hydrodynamic characteristics, and its influence on the bubble SMD in baffled STRs. The calculated fluid MI from the CFD simulations is validated through an LDA experiment, which underlines the validity of the designed CFD model to a certain degree. This agreement of both CFD and LDA results confirms that transient iterations in CFD simulations are suitable for the prediction of the fluid MI in baffled STRs. The former established linear dependence of the MI frequency to the impeller rotational speed is also confirmed herein.

Looking at the multiphase simulations, two conclusions can be drawn. Our PBM-CFD simulations show that the overall SMD reaches a steady value after a certain time frame, even if the periodical oscillations of the fluid MI are considered. Compared to the steady PBM-CFD case, the transient PBM-CFD case results yield a smaller overall SMD prediction with the same mesh scheme under steady PBM-CFD in one fixed rotational speed. This difference shows that fluid MI has an influence on the final value of the SMD. An experimental validation is yet be carried out, so that an actual assessment can be eventually performed. Nevertheless, the MI obviously has an influence, which should definitely be considered. The smaller the SMD, the greater the mass transfer for the same gas holdup. The influence of type one macro-instability on the bubble size distribution is numerically investigated through PBM-CFD modeling, while other types of macro-instability were not considered in the PBM-CFD modeling. More effort should be made to couple multi-influence of the various types of macro-instability on the SMD investigation in the future work.

The method developed in this contribution can help improve results obtained from less defined models, such as the multiblocks model. Due to the experimental validation performed on the liquid side, the fluid field estimation of our PBM-CFD simulation is quite precise. The overall SMD in the multiblocks model is greater than that in the transient PBM-CFD simulation involving the fluid MI. With this article, the authors hope to assist engineers or scientists working in the field of SMD prediction—be it to improve reactions, mixing processes, or similar—to obtain more realistic results for their multiphase systems. The model developed herein may provide more adequate information for the estimation of breakage and coalescence kernels of bubbles, because of its precise prediction regarding the fluid field. Moreover, As the mechanism modeling of the hydrodynamics is more and more improved through the experimental validation, more advanced process control methodologies as multiple predictive control, cybernetics, and fuzzy logic will be implemented for the gas-liquid STR online operation and control in the future.

#### Symbols

*Latin Letters*

: | Geometrical parameter of the tank (mm) |

: | Geometrical parameter of the impeller (mm) |

: | Geometrical parameter of the blade (mm) |

: | Cylindrical coordinates |

: | Velocity (m/s) |

: | Birth rate due to breakage/coalescence |

: | Death rate due to breakage/coalescence |

: | Breakage frequency () |

: | Breakage probability density function |

: | Collision frequency () |

: | Coalescence efficiency |

: | Macro-instability frequency () |

: | Nondimensional macro-instability |

: | Impeller rotational speed (rpm or r/s) |

: | Bubble size number function |

: | Gas flow rate (vvm) |

: | Diameter of bubble (m) |

: | Turbulence kinetic (m^{2}/s^{2}) |

: | Volume of bubble (m^{3}) |

: | Source term |

: | Weber number. |

*Greek Letters*

: | Density (kg/m^{3}) |

: | Viscosity (kg/m-s) |

: | Surface tension (N/m) |

: | Turbulence energy dissipation (m^{2}/s^{3}) |

: | Gas holdup. |

*Subscripts*

: | Tank |

: | Baffle |

: | Impeller |

br: | Breakage |

co: | Coalescence |

: | Gas |

: | Liquid. |

#### Competing Interests

The authors declare that there is no conflict of interests regarding the publication of this paper.

#### Acknowledgments

This work is part of the Collaborative Research Centre “Integrated Chemical Processes in Liquid Multiphase Systems” coordinated by the Technische Universität Berlin. Financial support by the German Research Foundation (Deutsche Forschungsgemeinschaft, DFG) is acknowledged (TRR 63).