CET Volume 86 DOI: 10.3303/CET2186190 Paper Received: 7 November 2020; Revised: 13 February 2021; Accepted: 14 April 2021 Please cite this article as: Maluta F., Paglianti A., Montante G., 2021, CFD Prediction of Gas-filled Cavity Structures in Aerated Vessels Stirred with Multiple Rushton Turbines, Chemical Engineering Transactions, 86, 1135-1140 DOI:10.3303/CET2186190 CHEMICAL ENGINEERING TRANSACTIONS VOL. 86, 2021 A publication of The Italian Association of Chemical Engineering Online at www.cetjournal.it Guest Editors: Sauro Pierucci, Jiří Jaromír Klemeš Copyright © 2021, AIDIC Servizi S.r.l. ISBN 978-88-95608-84-6; ISSN 2283-9216 CFD Prediction of Gas-filled Cavity Structures in Aerated Vessels stirred with multiple Rushton Turbines Francesco Malutaa,*, Alessandro Pagliantib, Giuseppina Montantea a Dipartimento di Chimica Industriale ‘Toso Montanari’, Alma Mater Studiorum – Università di Bologna, viale Risorgimento 4, 40136, Bologna, Italy b Dipartimento di Ingegneria Civile, Chimica, Ambientale e dei Materiali, Alma Mater Studiorum – Università di Bologna, via Terracini 34, 40131, Bologna, Italy francesco.maluta@unibo.it In this work, a turbulent gas-liquid stirred vessel is investigated by a fully predictive computational method, based on the Two-Fluid Model equations in the Reynolds Averaged Navier-Stokes formulation. The typical geometrical characteristics and operating conditions of aerated industrial fermenters are selected, which pose specific challenges due to the complex distribution of the gas phase and the cavity formation on the rear of the flat blades of the Rushton turbines. The results show that successful predictions of the gas cavities and the subsequent power consumption reduction in different operating conditions are obtained. The differences between the impellers positioned at different axial heights are assessed both in terms of gassed to ungassed power ratio and of gas flow rate treated by each turbine, showing that very a different behavior is obtained on the lowest impeller, with respect to the upper ones, thus providing a tool to quantify the individual impeller performances. The spatial distribution of the gas phase dispersed in the system and the resulting local mass transfer coefficient distribution are also discussed, observing a local decrease of the kLa due to the cavity formation. The presented model may contribute to the characterization of gas distribution, power requirements and local kLa in industrial fermenters, advancing the state of the art of Computational Fluid Dynamics tools in the context of fermentation intensification. 1. Introduction Fermentation Intensification (FI) closely follows in the footsteps of Process Intensification to drive the chemical and biotechnological processes towards a more environmentally sustainable path (Noorman et al., 2018). On the other hand, fermentative processes have the critical disadvantages tied to the limited operating conditions range bearable by the microorganisms. For instance, large substrate concentration fluctuations may trigger different metabolic pathways, causing loss of productivity (Maluta et al., 2020a), and heterogeneous distribution of oxygen may result in slower microbic growth rate (Pigou & Morchain, 2015) with a consequent sub-optimal use of the fermenter volume (Morchain, 2017). For these reasons, being able to reliably predict the gas-phase distribution inside of an aerated fermenter is paramount in the context of fermentation intensification. Computational fluid dynamics (CFD) is a versatile tool to obtain local information on the phase distribution inside a fermenter, therefore, it can be used to design real size as well as scale-down equipment (Noorman & Heijnen, 2017). CFD models have been used extensively to simulate fermenters (Montante et al., 2013) and gas-liquid stirred tanks in general (Shi & Rzehak, 2018). These models are often based on the two- fluid formulation of the Reynolds Averaged Navier-Stokes (RANS-TFM) equations, since this approach is rather computationally cheap, and it proved successful in describing a large range of multiphase disperse systems (Shi & Rzehak, 2018). Despite their widespread use, models based on RANS-TFM still require validation with experimental results and applications in different operating conditions to test their robustness, before being reliably adopted to industrially relevant situations. Recently, Maluta et al. (2020b) proposed a fully predictive computational method for the simulation of industrial fermenters, which caught the gas cavity formation behind the flat blades of Rusthon turbines (RT). The issue related to the CFD prediction of the gas- cavities in stirred tanks is well known in literature (Lane et al., 2005), and it is especially pivotal since the 1135 formation of aerated cavities reduces the power transferred from the impeller to the fermenter volume and consequently the oxygen transfer rate (OTR). In this study, the above mentioned fully predictive method is applied to different operating conditions to assess the gas cavity formation, the gassed to ungassed power ratio, the volume fraction distribution and the mass transfer coefficient in a gas-liquid vessel stirred by multiple RT. Based on correlations from the literature for single impeller stirred tanks (Nienow et al., 1977), at four selected combinations of impeller speed, N, and gas flow rate, QG, the loading flow regime establishes and gas cavities are expected behind the blades of the lowest impeller. Thus, the four cases provide a benchmark for assessing if the simulation strategy allows to follow the expected variations of the fluid dynamics and the mass transfer characteristics of the equipment with the operating conditions variation. Attention is devoted to the different behaviour of each of the four impellers, originating from the different gas flow rates treated by each of them and the consequent cavity structure. 2. Computational Model The simulations are based on the numerical solution of the Two-Fluid model (TFM) in the context of the steady state RANS equations, as implemented in ANSYS Fluent 19.3, that read as: ∇ ⋅ ( ) = 0 (1) ∇ ⋅ ( ) = − ∇ + + ∇ ⋅ ( + ) + + (2) Where is the − ℎ phase volume fraction, and its density and mean velocity vector respectively. is the pressure shared by the two phases, is the gravitational acceleration vector, and are the laminar viscous stress tensor and the Reynolds stress tensor of the phase . This last was modelled with the standard k-ε model extended to multiphase flows by considering phase-averaged properties (namely the mixture k-ε model). and are the drag and turbulent dispersion interphase momentum exchange forces, respectively. As in previous studies concerning similar systems (Maluta et al., 2019), other interphase momentum exchange terms are neglected. The drag force was implemented via a User defined function (UDF), and it is expressed in terms of bubble terminal velocity, , as: = ( − )‖ − ‖( − ) (3) The turbulent dispersion force was modelled with the Burns et al. (2004) model: = ( − )‖ − ‖ , ∇ − ∇ (4) In which is the turbulent viscosity and , , is the turbulent Schmidt number, fixed equal to 0.9. The correlation developed by Grace as reported in Clift et al. (2005) was adopted to calculate the bubble terminal velocity: = . 43 . − 0.857 (5) With = ( − ) ( )⁄ being the Morton number, that is function of the water surface tension, , assumed equal to 0.072 N/m, = ( − ) / is the Eötvös number and is the Sauter bubble diameter. When the term in square brackets is between 2 and 59.3, and are equal to 0.94 and 0.757 respectively, while when it is bigger than 59.3, and are equal to 3.42 and 0.441, respectively. In the four operating conditions considered in this work, listed in Table 1, Eq (5) can be adopted for the terminal velocity calculation of both ellipsoidal and spherical air bubbles (density ρG=1.2 kg/m 3) in water (density ρL=998 kg/m 3 and viscosity μL=0.001 Pa·s). The constant dB value adopted in the simulations was estimated from the correlations by Alves et al. (2002) for the impeller zone, as a function of the gassed power consumption per unit volume, as reported in Table 1 together with the corresponding terminal velocity, the gas Flow numbers, = / and the Froude numbers, = / for each of the four investigated conditions, referred to as Cases in the following. The gassed power consumption, Pg, required for the bubble size estimation was referred to the lower impeller only, adopting the values available for single impeller stirred tanks. In particular, it was obtained by multiplying the ungassed power drawn by a single RT, based on the impeller power number, by the Relative Power Demand (RPD) defined as the ratio between the gassed and the ungassed power consumption. The RPD values reported in Table 1 were estimated by the correlations developed by Smith et al. (1987). 1136 Table 1: Main operating variables and dimensionless numbers of the Cases studied in this work Case Name N (rpm) QG (L/h) Fl Fr RPD dB (mm) Ut (cm/s) Case I 250 200 0.030 0.14 0.75 2.06 20.2 Case II 250 275 0.041 0.14 0.69 2.15 20.8 Case III 337 275 0.030 0.25 0.69 1.35 15.2 Case IV 337 375 0.041 0.25 0.63 1.23 14.1 A local mass transfer coefficient, kLa, model from the literature (Maluta et al., 2019) was used to obtain the spatial distribution of both the liquid side mass transfer coefficient, kL, and the specific interfacial area in each grid cell, a, in the context of the two-fluid model, and it reads as: = 0.4 . . (6) = 6α ⁄ , α ≤ 0.34 34 α / / , 0.3 < α ≤ 0.54 34 α / / , α > 0.5 (7) Where is the air diffusion coefficient in water equal to 2×10-9 m2/s, is the local turbulent dissipation rate, is the liquid kinematic viscosity and is the volume of the computational grid cell. This model allows the local calculation of the liquid side mass transfer coefficient based on the turbulent dissipation rate distribution and introduces a simplified closure model to obtain the specific interfacial area in the TFM frame when phase segregation is present. Further information on the interphase mass transfer model adopted in this investigation can be found in Maluta et al. (2019). 3. Numerical domain and solution procedure The set of model equations described in Section 2 was solved in a computational domain representing the typical geometrical characteristics of a laboratory scale aerobic fermenter. The tank diameter, T, was equal to 0.23 m and the vessel was equipped with four equally spaced baffles of width equal to T/10. The impeller diameters, D, were equal to T/3, the lowest impeller had an off-bottom clearance equal to T/2 and the distance between the impellers was equal to T. The total height of the liquid volume, H, was of 4T, thus resulting in four equal standard geometry stirred tanks of height T, with the corresponding RT positioned at a relative height of T/2. A ring sparger of diameter 0.4D was positioned below the lowest impeller, at an off-bottom distance of T/5. The fermenter geometry was discretized in a 1 million cells structured mesh. The gas injection was modelled with a velocity inlet boundary condition and the gas exit from the top of the fermenter with a degassing boundary condition. The no-slip condition was assumed for each phase on the solid walls of the tank. The rotation of the impellers was modelled with the multiple reference frame approach, and the solution was iteratively solved with a pseudo-transient solver with pseudo time steps of 0.001. Convergence was assumed when the scaled residuals reached a value lower than 10 -5, and the gas volume flowrate exiting the system equaled the inlet gas flow rate, ensuring the integral mass balance. 4. Results In this section, the predictions of the gas cavities behind the flat blades of the RT, the subsequent reduction in the power drawn by the impeller, the gas volume fraction distribution, and the liquid side mass transfer coefficient, kLa are presented. The gas cavities, identified as those cells where the gas volume fraction exceeds the selected threshold value of 0.95, as predicted for Case II are shown in Figure 1. The cavities developed in the lowest impeller (Figure 1d) are significantly larger than those obtained on the other three impellers (Figure 1a-c). This result is expected, since the air sparger is positioned below the lowest impeller, whereas in the upper impellers, gas enters the impeller zone after being distributed by the impellers below. Similar results are obtained for Case I (not shown for brevity), with smaller cavities especially in the three upper impellers, consistently with the lower gas flow rate with respect to Case II. The gas cavities obtained from the simulation of Case III are shown in Figure 2 for the two lower impellers. In these conditions, gas cavities develop just behind the blades of the lowest impeller (Figure 2b), whereas the gas volume fraction behind the upper impeller blades (Figure 2a), does not reach 0.95. This result is consistent with the expected result, since in both Case II and Case III the gas flow rate from the air sparger is 1137 equal to 275 L/h, while the impeller rotational speed for Case III is higher than Case II. Moreover, the gas cavities in Figure 2b are smaller than those observed in Figure 1d. The model, thus, manages to predict the different gas accumulation observed in different operating conditions, as well as the differences between the impellers positioned at different heights. The results of Case IV show a very similar behavior. The area of the projection of the cavities on the lowest impeller on a plane parallel to the tank bottom is 17% of the area swept by the impeller, equal to πD2/4, for Case I, 19% for Case II, and 7% for Case III and Case IV. Figure 1: Iso-surface of gas phase volume fraction equal to 0.95 for Case II on the impeller blades at z/T=3.5 (a), z/T=2.5 (b), z/T=1.5 (c) and z/T=0.5 (d). The origin of the axial coordinate, z, is on the tank bottom. The mesh on the RT is also shown. Figure 2: Iso-surface of gas phase volume fraction equal to 0.95 for Case III on the impeller blades at z/T=1.5 (a), z/T=0.5 (b). The origin of the axial coordinate, z, is on the tank bottom. The mesh on the RT is also shown. The gas accumulation behind the impeller blades results in a decreased power drawn by the impeller. To quantify this effect, the relative power demand, RPD, is shown in Figure 3a for the four impellers of Case I-IV. (a) (b) Figure 3: Relative power demand (a) and gas pumping number (b) for the four impellers in the different operating conditions studied in this work. Figure 3a shows that the lowest impeller is subject to the highest power reduction, due to the presence of the largest cavities. Case I shows a steady reduction of the RTD moving from the lowest impeller to the highest, whereas in the other Cases the RPD of the lowest impeller is smaller than the other impellers, in which the RPD is rather constant. The behavior of the upper impellers shown in Figure 1 and 2 is consistent with this 1138 finding, since very small differences are noticeable in terms of gas accumulation behind the blades. In Case I, the gas cavities have a constant decreasing size, moving from the lowest to the highest impeller. As also observed in Figure 1 and Figure 2, higher gas volume fractions accumulate behind the impeller blades of Case I and Case II, with respect to Case III and Case IV, thus leading to a higher reduction of the gassed power drawn by the impeller. The predictions of RPD on the lowest impeller compare reasonably well with the values obtained from the available correlations of Smith et al. (1987), with deviation from the predicted RPD for the lowest impeller of Cases I and II around 5 %, while the deviations for Cases III and IV are around 12 %. It must be emphasized that the model employed to obtain these results is completely predictive and the RPD is a consequence of the gas volume fraction accumulation behind the impeller blades. Figure 3b shows the impeller gas pumping number, NQ,G, obtained dividing the gas flow rate radially pumped by each impeller by ND3. As can be observed, the lowest impeller treats the highest gas flow rate, being the sparger positioned just below it. For each Case, it can be observed that NQ,G sharply decreases in the two upper impellers and in all cases they radially discharge an almost equal amount of gas. In fact for Cases I and II, the average gas pumped by the three upper impellers with respect to the gas pumped by the lowest impeller is 34 % for the impeller at z/T = 1.5, 13 % for the impeller at z/T = 2.5 and 9 % for the impeller at z/T = 3.5. For Cases III and IV, the average gas pumped by the three upper impellers with respect to the gas pumped by the lowest impeller is 47 % for the impeller at z/T = 1.5, 38 % for the impeller at z/T = 2.5 and 38 % for the impeller at z/T = 3.5. This result shows that the two upper impellers work with a very similar gas flow rate that is considerably smaller than the gas flow rate handled by the lowest impeller, while the impeller at z/T = 1.5 works in an intermediate condition. Consistently with the air flow rate injected in the system, the average pumping numbers increase moving from Case I to IV. The gas hold up distribution in the tank is shown in Figure 4, on a vertical plane midway consecutive baffles. Figure 4: Gas volume fraction distribution for Case I (a), Case II (b), Case III (c) and Case IV (d). Figure 5: kLa distribution in s -1 for Case II (a) and Case III (b). Figure 4 shows that the lowest impeller works in the loading regimes, as also predicted by correlations from the literature. Cases I and II (Figures 4a and 4b), highlight that the upper impellers also work in limited recirculation regimes, with very little gas recirculated below the impellers and a rather large portion of the space above each impeller devoid of gas phase. Figure 4a also shows that the amount of gas around the impellers decreases moving from the lowest to the highest impeller. Cases III and IV (Figure 4c and Figure 4d), show that also the upper impellers work in the limited circulation regime, as can be observed by the small zone devoid of gas above the impellers, but in these Cases a portion of the gas is circulated below the impellers, thus achieving a better distribution in the tank volume. In fact in these Cases, the absence of gas cavities in the upper impellers allows a better recirculation of the gas, with respect to Case I and II and the consequent more homogeneous distribution of dispersed phase. The spatial distribution of the liquid side mass transfer coefficient is shown in Figure 5 on a vertical plane midway two consecutive baffles, for Case II (Figure 5a) and Case III (Figure 5b). A more homogeneous distribution of kLa is apparent in Case III, due both to the higher turbulent dissipation rate and the higher hold up in the system. The presented model also consistently describes the reduction of the mass transfer coefficient due to the gas cavities, as can be observed in the impeller zone. The volume average mass transfer coefficient for Case II is 13.7 h-1, while for 1139 Case III it is 30.1 h-1. Previous studies (Maluta et al., 2020b) highlighted that failing to predict the gas cavities can lead to an overall overprediction of the mass transfer coefficient up to 28 %, therefore this computational model is a step forward in increasing the understanding of gas-liquid stirred tanks, and a predictive tool aiding the design of fermentative processes. The model can be integrated with species mass transfer to obtain the local oxygen transfer rate, as done in Maluta et al. (2020b), with biochemical reaction rates by introducing a sink term in the mass balance equation and a suitable kinetic equation, allowing the calculation of the oxygen uptake rate distribution. The resulting computational approach based on the current work will allow to obtain local information on the fermentation evolution and spatial distribution of oxygen consumption, thus allowing for a comprehensive modelling approach of realistic fermentation processes. 5. Conclusions The RANS-TFM simulations of a multiple-impeller aerated stirred tank where the cavity formation behind the flat impeller blades takes place provide a consistent picture of the complex fluid dynamics of the equipment under different operating conditions with a fully predictive methodology. The power reduction associated with the gas cavities is correctly predicted by the model, with deviations from the RPD derived by literature correlations of the order of 10 %. The different gas flow rate treated by each turbine, that depending on the conditions could reach up to 91 % between the lowest and the highest turbine, generates differences in the power drawn by each impeller and this reflects on the impeller regime that transitions from loading to conditions close to complete recirculation. The formation of the cavities reflects on the volumetric mass transfer coefficient, with local reductions clearly observable in the impeller blade portion of the fermenter volume. References Alves, S.S., Maia, C.I., Vasconcelos, J.M.T., Serralheiro, A.J., 2002, Bubble size in aerated stirred tanks, Chemical Engineering Journal, 89, 109-117. Burns, A.D., Frank, T., Hamill, I., Shi, J.M., 2004, The Favre averaged drag model for turbulent dispersion in Eulerian multi-phase flows, In 5th international conference on multiphase flow, ICMF (Vol. 4, pp. 1-17). Clift R., Grace J.R., Weber M.E., 2005, Bubbles, drops, and particles, Dover Publications, Inc., Mineola, New York. Lane, G.L., Schwarz, M.P., Evans, G.M., 2005, Numerical modelling of gas–liquid flow in stirred tanks, Chemical Engineering Science, 60, 2203-2214. Maluta, F., Paglianti, A., Montante, G., 2019, Modelling of biohydrogen production in stirred fermenters by Computational Fluid Dynamics, Process Safety and Environmental Protection, 125, 342-357. Maluta, F., Paglianti, A., Montante, G., 2020b, Two-Fluids RANS predictions of gas cavities, power consumption, mixing time and oxygen transfer rate in an aerated fermenter scale-down stirred with multiple impellers, Biochemical Engineering Journal, 107867. Maluta, F., Pigou, M., Montante, G., Morchain, J., 2020a, Modeling the effects of substrate fluctuations on the maintenance rate in bioreactors with a probabilistic approach, Biochemical Engineering Journal, 107536. Montante, G., Coroneo, M., Francesconi, J.A., Paglianti, A., Magelli, F., 2013, Computational analysis of a vortex ingesting bioreactor for hydrogen production. Chemical Engineering Transactions, 32, 721-726. Morchain, J., 2017, Bioreactor Modeling: Interactions Between Hydrodynamics and Biology, ISTE-Wiley edition, Elsevier, London, UK. Nienow A.W., Wisdom D.J., Middleton J.C., 1977, The effect of scale and geometry on flooding, recirculation, and power in gassed stirred vessels, Proc. 2ndEur. Conf. Mix., BHRA Fluid Eng., Mons, BE, pp. F1-1-F1- 16 Noorman, H.J., Heijnen, J.J., 2017, Biochemical engineering’s grand adventure, Chemical Engineering Science, 170, 677-693. Noorman H.J., Van Winden W., Heijnen J.J., Van Der Lans R.G.J.M., 2018, Intensified Fermentation Processes and Equipment, RSC Green Chemistry: Intensification of Biobased Processes, Royal Society of Chemistry, Croydon, UK. Pigou, M., Morchain, J., 2015, Investigating the interactions between physical and biological heterogeneities in bioreactors using compartment, population balance and metabolic models, Chemical Engineering Science, 126, 267-282. Shi, P., Rzehak, R., 2018, Bubbly flow in stirred tanks: Euler-Euler/RANS modelling, Chemical Engineering Science, 190, 419-435. Smith, J.M., Warmoeskerken, M.M.C.G., Zeef, E., 1987, Flow conditions in vessels dispersing gases in liquids with multiple impellers, Biotechnology Processes: Scale-up and Mixing, Ho, C.S. and Oldshue, J.Y., AIChe, New York, pp 107–115. 1140