Investigation of continuous stirred tank reactors for improving the mixing in anaerobic digestion: A numerical study

.


Introduction
The increasing population and growing economy resulted in increased energy demand at a higher price which is a prime challenge in today's world.This has put the thrust on the exploitation of fossilized fuels and necessitated the usage of renewable energy resources for human survival [1].The usage of fossil fuels enhances greenhouse gas emissions and global warming and results in overall climate change.One of the key solutions to alleviate these challenges is fulfilment of energy needs from treatment of waste using its biochemical conversion to useful energy.Anaerobic digestion (AD) is one such technology which is being practiced worldwide to produce biogas from organic fraction of substrate and regarded as a reliable source of energy [2].The produced biogas can be used to generate heat and energy.AD undergoes various biochemical reactions in an oxygen-free environment and under controlled temperature and pH conditions [3].For efficient conversion of biodegradable waste into biogas, it is mandatory to have biological, chemical, and physical uniformity within the digester which can be fulfilled only if, the reactant is thoroughly mixed [4].The economics of AD is majorly influenced by the mixing within the digester and accounts up to 20% of total energy input to digester for maintaining uniformity of controlling parameters [5].The mixing system are solicited to avoid the settling down the substrate in bottom of digester.Proper mixing arrangement can contribute to reduction in substrate size and separation of gas from substrate whereas inadequate mixing can damage the microbial species and can result in formation of dead zones in digester [5].The ways of mixing mainly involve are compressed gas, pumped recirculation, and mechanical mixing-based methods.The modes of biogas recirculation, impeller mixing, and slurry recirculation have been explored experimentally by feeding different total solids (TS) concentration in reactor.It has been reported that the biogas production remains unaffected if fed with less than 5% TS content whereas properly mixed digesters loaded with more than 10% TS can produce 10-30% more biogas than unmixed digester [6].Recently, researchers have started to anticipate the mixing phenomenon using the computational fluid dynamics (CFD) tool which is one of the prominent ways to optimize the mixing and temperature uniformity in the AD's [7].The impact of changes in geometry, blade designs, and operating parameters like speed, TS, etc. can be evaluated through CFD without running the field tests.The work on optimizing the energy yield for reducing the operational cost of the AD reactors is inevitable for making the AD technology more competitive.The recent progress in using computational techniques is an efficient way to predict the performance of AD reactors and can be used to optimize the shapes, design, and other operational parameters for improving the performance of AD reactors.The study on egg shape reactors using mechanical draft tube-based mixing reported that the mixing energy level increases exponentially while average velocity in the digester increases linearly with increases of rotational speed of the propeller [8].Same study observed that draft tube-based mixing indicates better mixing as compared to pump recirculation and egg-shaped digesters have efficient mixing than cylindrical shape.The mixing study of anaerobic digester using impeller and pump recirculation in egg shaped digester concluded that combination of pumped recirculation and impeller-based mixing yields the highest level of mixing [9].Whereas, in absence of using impeller, pumped recirculation efficiently mixes the plane parallel to feeding pipe and possess poorly mixed situation in plane perpendicular to feeding pipe.The reduction in poor mixed zone from 33.6% to 29.6% were observed using 45 • conical bottom in the lieu of flat bottom [10].According to Bridgeman [4], biogas yield is independent of speed of impeller in a digester for having 2.5% TDS sludge [5].The study including non-Newtonian fluid flow in both lab-scale and pilot-scale anaerobic digesters concluded that the increase in TS content, fluid shows non-Newtonian pseudo plastic behaviour more strongly [11].It was also observed that an increase in the inlet flow velocity of the fluid does not affect the overall flow pattern and insignificantly reduces the dead zones in the digester.In another study, the impeller blade angle of 30 • and mixing time of 5 min was found optimum in terms of methane yield and TS degradation as compared to impeller blade angles of 45 • and 60 • [12].The existing literature does not provide the details of the effect of different impeller sizes on temperature distribution in greater than 10% TS feed reactors.
Therefore, the current study was conducted to analyse the effect of size and using number of impellers on power requirement in 5 lab scale anaerobic reactors.Four cubical and one cylindrical shaped reactor with variable impeller size were chosen to perform the investigation.The velocity and temperature profiles were drawn for each case for analysing the mixing within digester.The current comparative study was focused on the development of an optimized anaerobic digester based on defined mixing quality measurement variables.Additionally, the correlation between the mixing energy level and other decisive parameters such as dead volume, velocity gradient and uniformity index are also discussed.

Methodology
Impeller-based mechanical mixing is a challenging phenomenon to simulate in 3D since it involves the creation of a precise threedimensional mesh and extensive computations.In this work, we have used the numerical models in steady-state and transient mode using Ansys Fluent 17.2.

Governing equations 2.1.1. Non-Newtonian fluid rheology
The sewage sludge behaves as a shear-thinning non-Newtonian fluid [5].Viscosity is independent of the shear rate for Newtonian fluids whereas it is not so for non-Newtonian fluids.The flow patterns for non-Newtonian fluids vary from the Newtonian fluids' flow [11].The behaviour of these fluids can be understood using models like the power-law model and the Bird-Carreau model.The Bingham model and Herschel-Buckley model are used in case of the existence of yield stress which must be exceeded for the flow of the fluid [13].For this study, the power-law model was used for modelling the sewage sludge as it is a computationally effective as well as the simple representation that relates the dynamic viscosity (η) and shear rate ( γ) as depicted in eq.(1&2) [13].
Where k and n are known as the consistency coefficient and the powerlaw index, respectively.For n equal to 1, it represents the Newtonian fluid.The sludge exhibits the shear-thinning behaviour for n less than 1, and shear-thickening for n larger than 1.The density relates to the TS as given by eq. ( 3) [5].
TS 5.4% mimics wastewater sludge very much.In current study, the mixing of a concentrated (TS 12.1%) sludge was investigated.The rheological properties of the non-Newtonian fluid for different values of TS can be fed to the solver according to the data presented in Table 1.

Flow equations
The mass conservation equation for the single-phase fluid flow in a three-dimensional system is represented by eq. ( 4), [15].For this simulation study, the sludge has been assumed to be an incompressible fluid.
The sludge particles experience two types of forces, namely surface forces, and body forces.The surface forces are there because of the pressure and the viscous stress.Body forces exist mainly in the form of gravitational force.Therefore, the rate of change of momentum is given by the net forces on the particle which is shown below in eq. ( 5), [15,16].

∂ ∂t
where p, τ, and g → represent fluid's pressure, viscous stress, and gravitational force acting on the body.As given by Metzner and Otto (1957), the Reynolds number for non-Newtonian fluids is calculated with the help of algebraic expression given by eq. ( 6), [17].
Where D m is the diameter of the impeller, N denotes its speed in rotations per second (rps), k, n, and ρ are presented in Table 1.
In Ansys fluent, turbulent heat transport is modelled using Reynold's analogy to turbulent momentum transfer.The energy equation is therefore given as; ∂ ∂t Here, E is the internal energy, and H is total enthalpy which is defined as; Different turbulence models for the mixing of non-Newtonian fluids have been studied and used by many researchers.Standard k-ε, Realizable k-ε (RKE), Standard k-ω, and Shear Stress Turbulence k-ω models

Table 1
The rheological properties of the sludge for TS 5.4% and 12.1% [10], (Original data of Achkari Begdouri [14].have been investigated [18][19][20].The standard k-ε model is used in the cases of large Reynold's number [21].In general, the standard k-ε model lacks performance in cases related to swirling and mixing [22].The RKE turbulence model performs well for rotational models [23].Also, RKE offers similar results as the Reynolds Stress model with less computation time [19].The recommended model for turbulence is the standard k-ω model or the realizable k-ε model.The RKE method as the turbulence model was used for this study because this model shows the least amount of error for the fluids at TS 12.1% [19].The transport equations for calculation in the RKE model are depicted below in eq. ( 9)&10 [24].The transport equations for calculation in the RKE model are presented in the Ansys Fluent Theory guide [25] which are originally taken from Ref. [24].
∂ ∂t ∂ ∂t where, μ t is computed as follows, It follows some constraints of stresses which shows consistency with the turbulent flow physics.

Characterisation of mixing
The dead volume or stagnant zone is a widely accepted criterion to quantify the mixing.The stagnant zone is defined as the part of the reactor that has a velocity less than 5% of the maximum velocity in the digester [26].The fruitful understanding has provided by this parameter while doing comparative studies.There is no established value to conclude that good mixing is achieved however a dead volume lower than a predefined benchmark, can be understood as better mixing.Hence, for an improvement in mixing our objective is to reduce the dead volume and keep it as low as possible.In present study, optimized dead volume which is the volume lower than 5% of a predefined value of velocity was used.The predefined value was chosen as 0.64 m/s which is the maximum velocity possible in case of digesters D1 and D3.Other digesters have higher values of velocities as well.Since these higher values are not representing a major portion of the digesters instead just the miniscule volume around the blades of the impeller.
The water industry utilizes the velocity gradient as an important criterion to analyse the mixing process [5].Also, it helps in analysing the dead zones and mixed zones.For larger digesters, having a turnover time of 20-30 min, the critical velocity gradient lies in the range of 50-80 s − 1 [26].However, the critical value of velocity gradient for lab-scale digesters lies in the range of 7.2 s − 1 to 14.3 s − 1 [27].In a Cartesian coordinate system, it is calculated by the expression presented in eq.(11) where, u x , u y , and u z are the velocity of the cell in x, y, and z directions respectively [9,28].
The uniformity index is a very well-established tool to study the homogeneity of mixing in a vessel [29,30].For this work, we are using uniformity in velocity (UI v ) instead of the concentration of the fluid.UI v can be calculated using the expression given in eq. ( 12) [31], where v i and V i are the velocity and the volume of the ith cell.v denotes the volume-averaged velocity in the digester presented in eq. ( 13).UI v ranges from 0 (homogeneously mixed) to 1 (no homogeneity).The previous literature chosen the values ranges between 0.05 and 0.1 for good mixing but for our comparative study, a lower value of UI v than a predefined criterion can be said to achieve better homogeneity at a steady state of mixing [29].
The mixing energy level (MEL) is defined as the power consumed per unit volume of the mixing fluid.It can be evaluated with the help of the expression depicted in eq. ( 14) [32].Power is denoted by P and V represents the volume of the anaerobic digester.Any fixed value of MEL is debatable for achieving proper mixing but still, we can argue that if MEL is lower than a benchmark value then the mixing is efficient.MEL is 5-8 W/m 3 for those digesters which have a 20-30 min turnover time [33].
The power required by the impeller for mixing can be estimated with the help of the speed of the impeller (N) and the torque (T) experienced by it.Power can be estimated by the algebraic expression presented below in eq. ( 15), [8].Torque was observed as the moment of pressure and viscous forces experienced by the impeller blades during the simulation.

Model description
In present investigation, the substrate was assumed to be homogenous using physical treatment (less than 2 mm particle size) before adding to the digester as reported by various experimental studies [34] and subsequently, the fluid was taken as non-Newtonian by assuming it single phase which is also used by researchers in earlier studies [16,29].In current work, Five different digesters namely D1, D2, D3, D4 and D5 were created using cubical or cylindrical parts and inverted pyramid or   rotational speeds accordingly.Standard wall functions were chosen with no-slip, boundary conditions at the walls.SIMPLE i.e., the semi-implicit method for pressure-linked equations method was used for the calculations.For x, y, z, k, and ε the convergence condition was set at 10 − 5 .
The calculation was started with setting k = 0.001, and ε = 0.01 at a reference temperature of 295 K.

Grid independence test
A grid convergence index (GCI) is a tool that provides uniformity in the presented results.It provides an error band for the solution on a grid concerning a reference grid.It can be calculated using two grids but the recommendation is for using three or more grids [36].Five grids numbered from A to E were created using tetrahedral elements.E is the finest grid among all listed in Table 3.The maximum skewness of the element was used as the parameter for element quality and a skewness lower than 0.9 was targeted.GCI and time required for 1500 iterations along with the number of elements and maximum skewness is presented in Table 3. Root mean square error (e rms ) in the result can be calculated using eq.( 17) and the refinement ratio (r) was calculated using eq.( 18).Then, GCI was calculated by algebraic expression presented in eq. ( 16), [5]; here, e rms is the root mean square value of error in velocity between

Table 3
The number of elements, maximum skewness, refinement ratio, GCI, and the time required to complete 1500 iterations for the respective mesh.the coarse and fine grids observed at 1000 points, h 1 and h 2 represent the number of elements, u m,1 and u m,2 was the velocity in coarse and fine grids, respectively.Recommended refinement ratio is more than 1.1 between coarse and fine grids [36].The safety factor (F s ) is taken as 3 or 1.25 for two or more grids, respectively [37].GCI for grid D performed noticeably better than grids B and C.Although grid E outperformed grid D in terms of GCI, grid D was selected for simulation in this study due to its maximum skewness and the time needed for calculation as shown in Table 3.

Validation
Model validation was conducted using the CFD results for a cylindrical digester reported in literature [5].The dead volume was extracted for TS levels presented in Table 1 at the 100-rpm speed of the flat blade impeller and the comparison is depicted in Fig. 2(a).Also, the dead volume was obtained for mixing at 30, 50, 100, and 200 rpm operation while keeping the TS level fixed at 5.4%, and the comparison of the same is presented in Fig. 2(b).CFD results were within the band of 2% variation in dead volume from the literature values in Ref. [5].

Effect of stirrer speed on mixings
The previous literature study recommended the operating speed of 80 rpm in a cylindrical digester which was using a six-blade impeller [5].For maintaining the physical, chemical, and biological homogeneity due to cubical shape of digester, the mixing was performed using a two-blade impeller at 240 rpm for TS 12.1% in actual digester.The simulation for the same design was performed at 50, 100, 150, 200, 250, 400, and 600 rpm for TS 5.4%.The dead volume was optimized by keeping the maximum speed at 0.133 ms − 1 for other cases.At 250 rpm, the dead volume was below our target value of 10% but the velocity gradient was slightly lower than 7.2 s − 1 .At 400 rpm, the dead volume comes down to be lower than 5% and the velocity gradient improved to 8.54 s − 1 .As depicted in Fig. 3(a), keeping the mixing speed at 300 rpm was optimum operating point whereas Fig. 3(b) suggested the mixing speed of about 500 rpm.Increasing speed from 300 rpm to 500 rpm provided the improvement in dead volume and velocity gradient but nearly 4-fold increase in power consumption.Further study was done at a mixing speed of 250 rpm to maintain a low power consumption while leaving scope for improvement through the use of additional impellers or enhanced vessel designs.

Effect of impeller diameter and number of impellers on mixing
The mixing in the digester can be improved by increasing the diameter of the impeller and increasing the number of impellers [38,39].The first digester (D1) is the current digester with an impeller diameter of 5 cm, second digester (D2) is the improved digester with an impeller diameter of 7 cm as described in Table 2. Further digesters D3 and D4 were prepared using an impeller assembly of two impellers (centre to centre separation 8 cm) in place of the single impeller used in D1 and D2, respectively.Velocity profiles for D1, D2, D3, and D4 are presented in Fig. 4 for TS 12.1%.
Increasing the impeller diameter from 5 cm to 7 cm improved the dead volume from 66.64% to 45.85% for the single impeller configuration.Similarly, the dead volume was improved from 47.6% to 23.4% for the assembly configuration.The velocity gradient advanced to 3.6 s − 1 and 4.5 s − 1 but remain far from 7.2 s − 1 .A huge increment in power consumption of 4.7 times for a single impeller and 2 times for the assembly configuration was observed.

Vessel design
A digester with a large impeller assembly, discussed in section 3.2 as D4, could be further optimized by changing the vessel design.To improve the large dead volume at the vessel's edges for D4 shown in Fig. 4, the digester, D5, was constructed with an outer body that is composed of a conical bottom body and an upper body that is cylindrical as described in Table 2. Total height of the digester is 27 cm which is 1.3 times the major diameter of the vessel which makes it close to the desired value of 1.5 times [40].Velocity profiles, for D4 and D5 on a horizontal plane, passing through the impeller and other locations, are D2 shows a 30% improvement in DV with 5% better uniformity for 5 times more power consumption than D1.Digester D4 shows 50% improvement in DV with 11% better uniformity for 3 times more power consumption than D3.Digester D3 shows a 28% improvement in DV with 18% better uniformity than D1 but with 1.5 times more power consumption.As depicted in Fig. 5 (a)(b) &(c), digester D5 is better than D4 since DV, VG and UI v showed improvement with nearly the same or lower power consumption.
The substrate was heated as shown in Fig. 1 using a 250-W heater.
The vessel's temperature distribution was also shown, much like the velocity distribution, in order to track how the temperature profiles varied over a period of 10 min.Temperature profiles are only shown and discussed for D1 in order to prevent showing too many images.The temperature profiles from t = 1 to t = 10 min, at intervals of 1 min, are shown in Fig. 6 on a vertical plane which passes through the centre of the digester.Fig. 7 shows temperature profiles at t = 2 min and t = 10 min on five distinct horizontal planes.The temperature was 295K in the beginning.But, as soon as the heater was turned on and the impeller began to disperse the heat, the temperature rose.It took about 10 min to raise the temperature from 295 K to 310 K, however, Figs. 6 and 7 show that the substrate in the upper portion of the digester and close to the borders is not heating up.This analysis shed light on the mixer's  function in the digester's heat distribution.Without the mixer, the temperature within the digester would probably have been uneven, with some regions warmer than others depending on how near they were to the heater.

Conclusions
CFD is an effective tool to analyze the performance of a mixer.The optimum operating speed for a two blades impeller was found to be 250 rpm which is quite similar to 80 rpm using a 6 blades impeller as found in literature [5].Increasing the diameter from 5 cm to 7 cm brought improved performance based on dead volume, velocity gradient, and uniformity index.This improvement came at a cost of 2-4.7 times increase in power consumption.The cylindrical vessel provides better performance than the cubical vessel based on dead volume and uniformity index while the velocity gradient did not show any significant variation with a change in the vessel shape.The information obtained during the thermal analysis can also assist make such systems' design and operation better by ensuring that heat is dispersed uniformly and optimizing the positioning of the heater and mixer within the digester.The temperature is increased from 295 K to 310 K in approximately 10 min.Finally, a cylindrical digester (D5) is recommended with a 7 cm impeller and operating speed of 250 rpm with a two-blade impeller.The cylindrical equivalent of digester D3 can also be used for saving energy in the mixing process.

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.
inverted cone parts.Table2displays the height and side length of the pyramidal and cubical parts for D1, D2, D3, and D4.Similar to that, D5's cylindrical and conical part's radius and height are likewise shown.The diameter of the impeller and the number of impellers utilised in each digester are shown in Table2& Fig.1.The height of impeller blades were 20 mm in each digester.A rotation zone of 35 mm height, and has a diameter of 80 mm for D1, D3, and 120 mm for the remaining three digesters i.e., more than 1.5 times the size of the impeller, was built[11].The study was performed using (2 blades − 1 impeller) in D1, D2 and (4 blades-2impellers) D3, D4 and D5, respectively.A grid refinement study was conducted, and quality 4 mesh was chosen for this study as discussed in section 2.4.The vertical plane view and mesh for the geometry of D3 are displayed in Fig.1.The mesh was developed with the same quality for all the other digesters.Pressure-based steady-state solver was used.For the simulation of rotating fluid, Multiple Reference Frame (MRF) model or Sliding Mesh Model (SMM) could be used.SMM model produces slightly better results, but it requires higher computational effort than the MRF model.The difference in accuracy obtained from the SMM model and MRF model was not very significant[20].The multiple reference frame (MRF) method was used for modelling the impeller-based mixing[35].It divides the working volume into zones having different rotational speeds so different zones were assigned the

Fig. 1 .
Fig. 1.(a) Vertical 3-D view for geometry of D3 and (b) Mesh used for D3 (C) Top view of mesh (d) Mesh structure near heater wall (e) Shape and dimensions of impellers used in D1 (f) D2 (g) D3 (h) D4 &D5.

Fig. 4 .
Fig. 4. Velocity profiles on the horizontal plane at which the impeller is lying.

Fig. 6 .Fig. 7 .
Fig. 6.Temperature profile on the vertical plane for D1 at different time scale.

H
. Thakur et al.

Table 2
Details of reactors used in current investigation.