CFD simulation of a novel anaerobic-anoxic reactor for biological nutrient removal: Model construction, validation and hydrodynamic analysis based on OpenFOAM®

CFD simulation of a novel anaerobic-anoxic reactor for biological nutrient removal: model construction, validation and hydrodynamic analysis based on OpenFOAM® R. Blanco-Aguilera a,*, J.L. Lara b, G. Barajas b, I. Tejero a, R. Diez-Montero a,c a Group of Environmental Engineering, Department of Water and Environmental Sciences and Technologies, University of Cantabria, Avenida los Castros s/n, 39005 Santander, Spain b Environmental Hydraulics Institute (IHCantabria), University of Cantabria, Isabel Torres 15, 39011 Santander, Spain c Group of Environmental Engineering and Microbiology, Department of Civil and Environmental Engineering, Universitat Politècnica de Catalunya, c/ Jordi Girona 1-3, Building D1, E-08034 Barcelona, Spain * Corresponding author. E-mail address: [email protected]. (R. Blanco-Aguilera) Keywords: Computational Fluid Dynamics, Multi-environment, Turbulent flow, RTD analysis, Tracer tests, Biological Nutrient Removal Abstract AnoxAn is a novel multi-environment reactor for biological nutrient removal (BNR) from wastewater. Although its biological efficacy has been demonstrated on a pilot scale, hydrodynamics is observed to significantly affect the performance of AnoxAn. To study its complex hydraulic behaviour, a model based on Computational Fluid Dynamics 3D (CFD) is constructed using the OpenFOAM® open source toolbox and validated by experimental tests of Residence Time Distribution (RTD). Reactor elements represent a key factor in the modelling process. In this sense, the impeller of the anoxic zone is modelled as a flat disk, and the baffle after the anoxic zone as a porous media. According to CFD model simulations, stagnant, short-circuit zones and mixing quality are established and quantified. Finally, the influence on the hydrodynamics of reactor elements is also evaluated. The results of this detailed hydrodynamic analysis will form the basis for the design and optimization of scalable AnoxAn configurations. This is the accepted manuscript of the article that appeared in final form in Chemical Engineering Science 215 : (2020) // Aricle ID 115390, which has been published in final form at https://doi.org/10.1016/j.ces.2019.115390. © 2019 Elsevier under CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/) 1. Introduction For many years, the main objective of wastewater research has been to achieve the required efficiency of biological processes to meet regulations and preserve the ecological and healthy status of water bodies (rivers, lakes, reservoirs, oceans, etc.). Specifically, great efforts have been made to design and improve nutrient removal processes (i.e.: nitrogen (N) and phosphorus (P)) due to the increasing requirements of Wastewater Treatment Plants (WWTPs). However, conventional biological nutrient removal (BNR) processes in WWTPs require a complex treatment system and entail several environmental impacts. First, the aerobic reactor should be large enough to carry out ammonia oxidation (nitrification) and must be coupled with non-aerated compartments (anaerobic and anoxic). Therefore, a large volume is needed compared to organic matter removal processes. This issue becomes crucial in cases where land availability is limited and in cases of existing wastewater treatment plants that need to be upgraded to BNR plants. In addition, the high energy demand of the nitrification process, together with the need for recirculation pumping between the different compartments or reactors and the mixing of the non-aerated ones, result in a significant increase in energy consumption. In this context, multi-environmental biological reactors with high compaction and efficiency have been developed to reduce the energy consumption and land use of conventional BNR treatment trains (Kwon et al. 2005, Yerushalmi et al. 2011, Tejero et al. 2010, Tejero et al. 1991, Martin et al. 2012). Of special interest is the anaerobic-anoxic reactor AnoxAn, developed and patented by Tejero et al. (2010). AnoxAn is a continuous upward-flow sludge blanket reactor that unifies in a single reactor the anaerobic and anoxic zones necessary for the biological nutrient removal of conventional activated sludge from wastewater (Diez-Montero 2015, Diez-Montero et al. 2015, 2016). Due to its low energy consumption and minimal land use, the AnoxAn concept and technology can potentially be applied for the upgrade of a WWTP or in a new WWTP with limited space availability. In addition, primary settling tanks could be reused as anoxicanaerobic reactors to develop innovative and compact treatment systems (Diez-Montero et al. 2019). Finally, the anaerobic-anoxic biological functioning of AnoxAn is intended to be coupled to an aerobic reactor (for residual organic matter removal, phosphate uptake and nitrification) and a secondary settling unit (or final filtration step), to complete the BNR treatment train (Diez-Montero 2015). The AnoxAn reactor consists of an anaerobic zone at the bottom (receiving the influent wastewater), before an anoxic zone above (receiving a nitrate-rich current from a subsequent aerobic reactor). In addition, a clarification zone is achieved in the upper part, avoiding the escape of large quantities of biomass. One of the main objectives of the reactor setup is to establish an anoxic-anaerobic hydraulic separation, i.e. to maintain an insignificant concentration of nitrates in the anaerobic zone, achieving at the same time adequate mixing conditions in both zones and maintaining the continuous flow of the effluent through it. For this purpose, the reactor has independent mixing systems in each zone: a recirculation pump provides the mixture in the anaerobic zone and a mixing impeller in the anoxic zone. In addition, the deflectors and baffles improve hydraulic separation and retention of suspended solids within the reactor. Concretely, the upper baffle, BLAS® (Tejero et al. 1991) was originally conceived as a support media for biofilm, but in AnoxAn is used as a head loss generator reducing the velocity of the fluid flow. Specific elements with these characteristics and level of interference in the flow pattern increase the hydrodynamic complexity of the reactor, and could generate preferential flows and dead zones, reducing the overall performance of the system (Al-Sammarraee et al. 2009, Liu et al. 2018, Plascencia-Jatomea et al. 2015, Yan et al. 2015). The viability of anoxic-anaerobic hydraulic separation of AnoxAn was tested in a 48.4 L prototype by means of Residence Time Distribution (RTD) analysis (Diez-Montero et al. 2015). A hydraulic model based on compartments was constructed and validated with experimental traceability tests. This model was a combination of complete mixed compartments and plug flow with axial dispersion compartments, implemented to describe the non-ideal flow of the reactor. The model predicted with high accuracy the experimental records (local measurements without spatial resolution). Then, it was applied to evaluate hydraulic anoxic-anaerobic separation. However, this type of modelling cannot provide complete information on hydrodynamics within the reactor. In a later work, the biological behaviour of the reactor in the treatment of municipal wastewater was studied (Diez-Montero et al. 2016). The results demonstrated the feasibility of the reactor concept. However, it was observed that some physical characteristics of AnoxAn significantly affect its performance, hydrodynamics being clearly relevant. In order to optimise the reactor configuration and propose other reactor configurations applicable on a large scale, Diez-Montero (2015) pointed out a deeper and more complete hydrodynamic analysis. In fact, RTD experimental tests (Levenspiel, 1999), usually coupled to hydraulic models based on compartments such as the tank-in-series and the dispersion models, have been widely used for hydrodynamic analysis in wastewater treatment reactors. It was observed that mixing conditions (Hu et al. 2012, Olive et al. 2005, Yerushalmi et al. 2013), flow type and characteristics (Behzadian et al. 2013, Autumn et al. 2007, Gomez 2010, Ji et al. 2012, Saragai et al. 2010), dead volume (Autumn et al. 2007, Hu et al. 2012, Ji et al. 2012, Sarathai et al. 2010), channeling (Gómez 2010, Nemade et al. 2010, Zeng et al. 2005) and dispersal (Ji et al. 2012, Nemade et al. 2010, Yerushalmi et al. 2013, Zeng et al. 2005) were observed as the most important characteristics. However, experimental RTD analysis techniques require a lot of time and resources. In some cases, the complexity of experimental tests makes them impracticable in large-scale reactors (Fernandez, 2012). In addition, experimental RTD models and compartment-based hydraulic models do not contain any information on spatial flow and concentration camp resolution (Plascencia-Jatomea et al. 2015, Qi et al. 2013). The latter can be overcome by combining and developing advanced mathematical models and computational simulation. The application of numerical techniques to engineering has experienced great growth in recent decades, with Computational Fluid Dynamics (CFD) being one of the approaches with the greatest impact. The use of CFD in wastewater treatment processes is growing rapidly and is being applied in the resolution of complex problems (Angeloudis et al. 2016, Brannock et al. 2010, Liu et al. 2018, Klusener et al. 2007, Wicklein et al. 2016, Zhang et al. 2016). Regarding hydrodynamic analysis, combination of RTD and CFD for reactor analysis and optimization has been previously applied for different reactors (Brannock et al. 2010a-b, Climent et al. 2018, Le Moullec et al. 2008, Pereira et al. 2011, Plascencia-Jatomea et al. 2015, Terashima et al. 2009). The present work follows the validation performed by many of them regarding the hydrodynamic field, which is only based on RTD tests (Brannock et al. 2010a-b, Plascencia-Jatomea et al. 2015, Terashima et al. 2009). Moreover, all the studies aforementioned were carried out with commercial codes, being the present research an opensource approach that completes the already existing ones. Hydrodynamic in stirred reactors has also been modelled and analysed for several configurations (Bridgeman 2012, Bai et al. 2008, Choi et al. 2004, Qi et al. 2013). Different approaches can be undertaken for a correct representation of the impeller: being MRF (Multi Reference Frame) approach (Bai et al. 2008, Bridgeman 2012, Renade 2002, Wu et al. 2009) or momentum source approach common ones. The latter is a classical simplification that allows a significant computational resource saving and has been widely and successfully used in wastewater field, all of them carried out with commercial codes (Brannock 2003, Climent et al. 2018, Climent et al. 2019, Rehman 2016). Moreover, CFD modelling allows deeper hydrodynamic analysis including identification and location of dead zones, tracking of velocity profiles and flow patterns, mixing performance, and determination of the distribution of tracer concentration within the reactor (Climent et al. 2018, Dapelo et al. 2018, Michalopoulos et al. 2018, Plascencia-Jatomea et al. 2015, Terashima et al. 2009, Trab et al. 2015, Wei et al. 2019). In addition, this advanced knowledge will be essential for the operational optimization of reactors with complex hydrodynamic behaviour such as AnoxAn in terms of avoiding the presence of dead flow zones or preferred channelling in the design. Finally, a calibrated and validated CFD model is more efficient than other approaches in terms of testing other forms or combination of different elements in a reasonable time. Concretely, hydrodynamic evaluation is crucial in multi-environment reactors due to their specific shapes and combination of different elements, such as baffles and mixing devices, which interfere in the ideal hydraulic performance creating complex flow regimes (Kwon et al. 2005, Yerushalmi et al. 2011, DiezMontero, 2015). However, to our knowledge, only few hydrodynamic studies have been performed based on CFD for multi-environment reactors and its elements influence (Calder et al. 2013). Finally, it has been widely proved that an optimum hydraulic operation ensures an adequate biological performance efficiency (Climent et al. 2018, Arnaldos et al. 2018, Wei et al., 2019), being hydrodynamic understanding and analysis a critical step for designing process. The objective of this study is to develop a CFD model of the novel AnoxAn anaerobic-anoxic reactor and to analyse its hydrodynamic behaviour in order to identify the key features of the reactor configuration that cannot be achieved with conventional RTD experimental procedures. These results are crucial for the optimization and future design of large-scale reactor implementations. In addition, modelling and simulation of reactor elements such as baffles and deflectors would provide a deeper hydrodynamic understanding that could be applied for the development and optimization of other multi-environmental reactors and conventional water treatment reactors. The model is built in the OpenFOAM® (v18.12) open source toolbox (Weller at al. 1998), and calibrated and validated with RTD experimental tests and simulations of previous models. 2. Materials and methods 2.1. Description of experiments 2.1.1. Bench scale reactor setup To evaluate the model's ability to reproduce the hydrodynamics of AnoxAn, a series of experiments were conducted on a prototype reactor built on a bench scale (Diez-Montero et al. 2015). AnoxAn prototype, see Fig. 1a, is a 48.4L upflow reactor made of polymethyl methacrylate (PMMA). The reactor is divided into three different zones: a clarification zone at the top (4.0 L; 8% of the total volume), an anoxic zone below (32.0 L, 66% of the total volume) and an anaerobic zone at the bottom (12.4 L; 26% of the total volume). Geometrically, it consists of an internal square section of 0.20 x 0.20 m2 and a height of 1.30 m. A cross section of the detailed reactor geometry based on the square section side (D = 0.20 m) is shown in Fig. 1e. Mixing devices in the reactor consist in a Heidolph RZR2 000 impeller (100 rpm) (Fig. 1b) for the anoxic zone and a peristaltic pump Watson Marlow 313U for the continuous internal recycle of the anaerobic zone. Besides, a complex 3D geometry polyethylene (PE) baffle of 0.039 m height is placed (Fig. 1c) between the anoxic and clarification zones in order to improve suspended solids retention inside the reactor and reduce the up-flow velocity. Finally, a polyvinyl chloride (PVC) deflector of 4 cm width along the wall is introduced (Fig. 1d) to enhance the hydraulic separation between anoxic and anaerobic environments. Figure 1. (a) 3D scheme of the bench scale AnoxAn reactor, (b) Impeller, (c) Baffle between anoxic and clarification zones, (d) Deflector between anaerobic and anoxic zones and (e) Detailed cross section geometry based on the square section side (D = 0.20 m) 2.1.2. Experimental RTD conditions Two experimental and a simulated tracer tests in clean water were performed in AnoxAn to characterise the liquid phase flow pattern. All details for the experiments are shown in Table 1 and Fig. 2. The AnoxAn reactor was designed for a hydraulic residence time (HRT) up to 5 h (depending on the organic load applied) and for all the experiments performed the inlet stream flow Qin is 10.4 L/h, internal recycle rate (ratio between internal recycle stream flow and inlet stream flow, RIR = QIR/Qin) is 5.77 and nitrate recycle rate (ratio between nitrate recycle stream flow and inlet stream flow, RNR =QNR/Qin) is 2.98. Two pulse RTD tests were carried out for hydraulic characterization of AnoxAn: (i) RTD1 injecting the tracer through the inlet stream (see RTD1 in Table 1 and Fig. 2a) and (ii) RTD2 injecting the tracer through the nitrate stream (see RTD2 in Table 1 and Fig. 2b). For both pulse experiments, a concentrated solution of sodium chloride (NaCl, 350 g/L) was used as tracer. Moreover, in RTD1, the volume injected was 58 mL and RTD2 was 40 mL. Output concentration was estimated measuring the conductivity through a linear relationship between them (Tang et al. 2004, Martín-Dominguez et al. 2005) with a Hach CDC40103 probe connected to a HQ30d meter. In addition to pulse experiments, and in order to better assess the hydraulic separation between the anoxic and anaerobic zones, a step RTD test (RTD3) was simulated using the calibrated and validated hydraulic compartment-based model described in Diez-Montero et al. (2015) and presented as supplementary material. In RTD3, a constant concentrated solution of tracer (10 mg/L) is continuously injected in the nitrate recycle stream and tracer concentration is continuously measured in the outlet and in both anaerobic and anoxic zones (see RTD3 in Table 1 and Fig. 2c). Table 1. RTD tests conditions RTD experiment Type Tracer injection location (m) Tracer injection duration Tracer concentration measurement (m) RTD1 Pulse. Experimental Inlet stream Pi (-0.1, 0, 0.30) tpulse = 3 s Outlet Pout (-0.1, 0, 1.295) RTD2 Pulse. Experimental. Nitrate recycle stream PN (0.1, 0, 0.65) tpulse = 3 s Outlet Pout (-0.1, 0, 1.295) RTD3 Step. Simulated. Nitrate recycle stream PN (0.1, 0, 0.65) Continued injection Anaerobic zone Panae (0, 0, 0.30) Anoxic zone Panox (0, 0, 0.80) Outlet Pout (-0.1, 0, 1.295) Figure 2. Schematic diagram of pulse tracer tests (a) RTD1 (b) RTD2 (c) RTD3 2.2. Numerical model setup In this section, the computational setup of the CFD model is presented. The model is validated based on RTD experiments described in the previous section. Based on the results from the numerical simulations, a comprehensive hydrodynamic analysis for the different zones of the reactor is performed. For that purpose, a transient flow analysis is needed. At this respect, unlike RTD analysis, transient CFD models provide high spatial flow and tracer concentration resolution, along with time history of the latter. Besides hydrodynamics, reactor elements are also modelled. At this aim, the local mixing effect of the impeller is reproduced by means a flat disk approach (See section 2.2.1.3) and the baffle situated between anoxic and clarification environments is simulated as a porous media (See section 2.2.1.4). Both approaches led to a significant computational cost saving. 2.2.1. Computational Fluid Dynamics. Governing equations and models. 2.2.1.1. Hydrodynamic model Hydrodynamics of AnoxAn are simulated solving Navier Stokes (NS) equations for turbulent and incompressible flow. Numerical resolution of turbulent flows can be achieved through different approaches with several approximation degrees. Reynolds Averaged Navier Stokes (RANS) simulation is the most widely used approach in engineering due to its relative simplicity and lower computational cost. In RANS simulation, all turbulent scales are simulated by modelling, introducing a time averaging to the variables in order to separate their ensemble value and the fluctuant one. RANS equations include continuity (Eq. 1) and momentum conservation (Eq. 2) equations, linking the pressure and the velocity. 𝜕𝑣 𝜕𝑥=0 (1) 𝜕𝑣 𝜕t + 𝑣𝜕𝑣 𝜕𝑥=−1 𝜌𝜕P 𝜕𝑥+(µ+µ) 𝜕 𝜕𝑥󰇧 𝜕𝑣 𝜕𝑥+𝜕𝑣 𝜕𝑥 󰇨+ 𝑔+𝑓 (2) where 𝑣 is the ensemble velocity vector, ρ is the fluid density, μl is the fluid dynamic viscosity, μt is the eddy viscosity, p is the pressure, g is the gravitational acceleration and fbaffle the resistant force (Fbaffle) produced by the baffle per unit of volume normalised by density. In addition, model closure equations are needed for the turbulent stress tensor: In this work standard k-ε model (Launder and Spalding, 1972) is used (Eqs. 3-4): 𝜕(𝜌𝑘) 𝜕𝑡 +𝜕 𝜕𝑥(𝜌𝑘𝑣)=𝜕 𝜕𝑥󰇩µ+µ 𝜎𝜕𝑘 𝜕𝑥󰇪 + 𝐺− 𝜌𝜀 (3) 𝜕(𝜌𝜀) 𝜕𝑡 +𝜕 𝜕𝑥(𝜌𝜀𝑣)=𝜕 𝜕𝑥󰇩µ+µ 𝜎𝜕𝜀 𝜕𝑥󰇪+𝐶𝜀 𝑘𝐺−𝐶𝜌𝜀 𝑘 (4) in which k is the turbulent kinetic energy, ε is the turbulent kinetic energy dissipation rate, G is the production rate of turbulent kinetic energy, σk, σε, C1ε and C2ε are the k-ε model standard constants , and µt is the turbulent viscosity defined by Eq.5: µ= 𝜌𝐶𝑘 𝜀 (5) Numerical values for the constants are: Cµ = 0.09 C1ε=1.44 C2ε=1.92 σk=1.0 σε = 1.3 2.2.1.2. Tracer transport model For tracer test simulation, the resolution of a transport equation without chemical reaction is modelled in the software. The expression is obtained by means of averaging the general transport equation for turbulent flow (Eq. 6): 𝜕(𝜌𝑚) 𝜕𝑡 +𝜕 𝜕𝑥(𝜌𝑣𝑚)=𝜕 𝜕𝑥󰇩𝜌𝐷+𝜇 𝑆·𝜕𝑚 𝜕𝑥󰇪 (6) For the k-esim specie. 𝐷 is the self-diffusion coefficient of the tracer. In order to only reproduce the convective transport, this value must be very low, at least 10-10 m2/s (Fernandez 2012). In this work 𝐷 is set to 10-20 m2/s to ensure to convective behaviour.   term represents the turbulent diffusion, in which Schmidt number appears (Sct =  · ), being µt turbulent viscosity and Dt turbulent diffusivity. The value of the Schmidt number used for the tracer transport is 0.8 as it provided good agreement with experimental results, being this parameter typically between 0.7 and 1 (Brannock 2003, De Clercq 2003). The resolution method for tracer transport uses the concept of fluid mixture. First continuity, momentum and energy equations are solved in the mixture field. Next, tracer transport is performed. 2.2.1.3. Impeller model The impeller plays a key role in the hydrodynamics of the reactor because it is the main source of inner momentum supply. Because the purpose of this work is not to define in detail the hydrodynamic patterns produced around by impeller blades and because the rotation speed of the impeller is not very high creating a very significant turbulent flow patterns, a simplified modelling of the impeller has been performed. To do this, the flat disk hypothesis has been used (Jasak et al. 2019, Seb 2017), where the momentum induced by the impeller to the fluid is introduced into the model from a compound velocity field at an axial, radial and tangential velocity (i.e.: va, vr and vt), which are defined according to the rotation speed and the relative position of the different nodes that define the disk with respect to its centre (see Fig. 3). Those velocities are used a boundary conditions at the flat disk surface. This approach has been successfully used and validated in wastewater research (Brannock 2003, Climent et al. 2018, Climent et al. 2019, Rehman 2016). This approach has the advantage that it has a lower computational cost but allows an adequate flow characterization in the reactor. Figure 3. Schematic scheme of the impellers flat disk approach Velocity boundary condition is set for three impeller velocities (Eqs. 7-9) as: 𝑣 󰇍 󰇍 󰇍 󰇍  =𝑈, · 𝑎 (7) 𝑣 󰇍 󰇍 󰇍  =𝑈, · 𝑟 󰇍 󰇍 󰇍  (8) 𝑣 󰇍 󰇍 󰇍  =𝑈, ·𝑡 (9) being 𝑎 the vector normal to the impeller, 𝑟 the radial vector and 𝑡 tangential vector respectively. 𝑈, represent the module of three characteristic velocities. The impeller has a rotational velocity of 100 rpm and a radial flow behaviour. The following values are set in the model for the defined velocity vectors as they provided good agreement with experimental results: 𝑣 󰇍 󰇍 󰇍 󰇍  =𝑈,·𝑎=0· 𝑎=0 (10) 𝑣 󰇍 󰇍 󰇍  =𝑈, · 𝑟 󰇍 󰇍 󰇍  =0 · 𝑟 󰇍 󰇍 󰇍  =0 (11) 𝑣 󰇍 󰇍 󰇍  =𝑈,·𝑡=100 𝑟𝑝𝑚 ·𝑡= 10,5𝑟𝑎𝑑 𝑠·𝑅·𝑡 (12) Being R the radial distance to each node in the impeller flat disk to the impeller centre in meters. Rotational speed of the impeller is 100 rpm, the diameter is 0.08 m and the kinematic viscosity of the fluid is 10-6 m2/s. Taking these into account, the Reynolds number yields 10 670, and in consequence, fully developed turbulent flow can be considered in that zone (see Couper at al. 2010). 2.2.1.4. Baffle model The baffle between the anoxic and clarification zones consists of a plastic frame with a 3D complex geometry (see Fig. 1c). The realistic definition of its real geometry introduces a complexity in the computational mesh that results in a 30% increase in the number of elements of the computational mesh, with the consequent increase in the computational cost. To avoid this increase, the baffle is modelled as a porous media that simulates a pressure drop in the velocity field as a momentum sump. It is modelled by means of a Darcy type flow model defined in Equation 13: fbaffle = k·|U 󰇍 󰇍  |; ∀ Baffle region (13) Being k the friction coefficient, whose value is set to 500 000 kg/s (Mihovilovic 2010). To calculate Reynolds number inside the baffle, expression from Burcharth et al. (1995) and Losada et al. (2016) has been followed. In this regard, Reynolds number of the fluid inside the baffle is between 40 and 75, which is consistent with the range of applicability of Darcy’s flow modelling. Comparing the experimental RTD results (circles) to the CFD model simulations (green lines) it can be observed that the predictions of the model agree with the experimental measurements. This verifies that the NS solutions for the turbulent flow and diffusion-convection equations reliably represent tracer transport within AnoxAn. In addition, the coefficient of determination R2 has been determined in order to quantify the fit between the simulated and experimental results (Diez-Montero et al. 2015, Fang et al. 2011, López et al. 2010, Makinia et al. 2006), obtaining high values in all the cases (Table 4). Table 4. R2 coefficient of determination for different CFD models RTD Curve R 2 RTD1 0.98 RTD2 0.96 RTD3 – Anoxic curve 0.99 RTD3 – Clarification curve 0.99 RTD3 – Anaerobic curve 0.99 Moreover, experimental HRT in pulse tests and numerical HRT taken from CFD models are compared for further validation (Brannock et al. 2010a-b, Plascencia-Jatomea et al. 2015, Climent et al. 2018) and results are shown in Table 5. Both HRT have been calculated cutting the curve in the time corresponding to the last experimental measurement. The HRT numerical result is observed to have a difference of less than 7% in both cases. Therefore, modelling approaches followed for single elements, i.e. impeller and deflector, are considered satisfactory. Besides, the difference between real (experimental) and theoretical HRT is also discussed. It is observed that real HRT is higher than the theoretical value, as expected. This means that not the overall volume is being used in the mixing process, resulting in dead volumes (Eq. 20) or stagnant zones (Climent et al. 2018). It should be highlighted that cutting the curve in the time corresponding to the last experimental measurement, it could lead to slightly underestimated HRT, due to the tracer mass not taken into account in the tail of the curve. Therefore, the real dead volumes could be slightly lower than the calculated ones, as shown in Table 5. Vd = 󰇡 1−   󰇢 · 100% (20) Table 5. HRT comparison for theoretical, experimental and CFD models Experiment HRTtheo HRTexp HRTCFD Dead Volume RTD1 124 min 97 min 96 min ≤22% RTD2 50 min 44 min 47 min ≤12% 3.2. Hydrodynamic analysis based on RTD curves Analyzing the experimental RTD measurements (Fig. 8, black circles), the following can be stated: With respect to the pulse tracer tests (Fig. 8a and 8b), the non-ideal AnoxAn flow pattern is analysed. The time evolution of the RTD curve allows to confirm that they are between the Continuous Stirred Tank Reactor (CSTR) and a Plug Flow Reactor (PFR). This deviation in ideal flow patterns is a consequence of the presence of some preferential flow and channelling zones in the AnoxAn. In addition, in both RTD1 and RTD2 (Fig. 8a-b) experiments a remarkable tailing is observed. This is the result of the presence of stagnant or dead flow zones, where tracer transport is low, resulting in higher tracer concentration. The existence of stagnant zones means that the entire reactor volume is not being used efficiently, which may lead to a decrease in the actual HRT compared to the design value. However, the RTD curves do not provide information on the location and size of the channelling zones or dead volumes. This information can be crucial for the optimization of the AnoxAn design and its scalability. Thus, using a calibrated and validated CFD model, an additional hydrodynamic analysis can be performed. Detailed velocity fields and tracer concentration can be tracked in order to better understand the flow and mixing mechanisms within the reactor and identify dead flow zones. Regarding step tracer test (Fig. 8c), a remarkable hydraulic separation between anoxic and anaerobic zones is confirmed. In fact, the steady-state tracer concentration reached in anaerobic zone is 25% of the concentration observed in anoxic and clarification zones. In addition, a delay in the stabilization of the concentration in the clarification zone is observed compared to the anoxic one. According to compartmentbased hydraulic models built (Diez-Montero et al. 2015), this delay is due to the influence of the baffle, which in theory, reduces the up-flow velocity. However, a further hydrodynamic analysis based in CFD techniques are needed to study that hypothesis. 3.3. Hydrodynamic analysis based on CFD simulations Dead flow zones and short-circuiting can be studied based on velocity field analysis. In this sense, CFD model results regarding velocity vectors and/or streamlines have been previously used in water treatment studies (Arnaldos et al. 2018, Brannock et al. 2010, Climent et al. 2018, Plascencia-Jatomea et al. 2015, Rehman 2016). In this regard, the velocity magnitude (Fig. 9a, Fig. 9d and Fig. 9g), the vertical velocity component (Fig. 9b, Fig. 9e and Fig. 9h) and flow streamlines (Fig. 9c, Fig. 9f and Fig. 9i) for the most representative zones in AnoxAn are presented in Fig. 9. A cross section XZ is represented together with horizontal cross sections XY at different levels to better visualise the magnitude and vertical component of the velocity fields. For the streamlines, a full 3D representation of the main flow patterns is presented. Note that for velocity magnitudes, a logarithmic scale has been used for better visualization. For vertical velocities, negative downstream values are shown in blue and positive upstream values are shown in red. Different vertical velocity scales have been used for the flow lines in order to obtain more detailed information for each zone. Figure 9. Velocity fields in AnoxAn: (a) velocity magnitude in anoxic-clarification transition zone, (b) vertical velocity in anoxic-clarification transition zone, (c) streamlines in anoxic-clarification zone (d) velocity magnitude in the main anoxic zone, (e) vertical velocity in the main anoxic zone, (f) streamlines in the main anoxic zone, (g) velocity magnitude in anaerobic-anoxic transition zone, (h) vertical velocity in anaerobic-anoxic transition zone and (i) streamlines in anaerobic-anoxic transition zone. Experimental conditions: ωimp = 100 rpm, Qin = 10.4 L/h, Qnit = 31.0 L/h (COLOURED) 3.3.1. Anoxic-clarification transition zone Velocity fields for the upper anoxic and clarification zones which are separated by the baffle (in grey) are shown in Fig. 9a-c. The highest velocity profiles in that section are noted in the outlet. Similarly to other research (Arnaldos et al. 2018, Climent et al. 2018, Plascencia-Jatomea et al. 2015), those high velocities enhance a preferential flow channelling through it, creating a stagnant zone in the opposite corner of the outlet. In that zone, near zero velocities are observed. The aforementioned preferential path formation is represented in Fig. 9c, where it is observed that principal streamlines avoid outlets opposite corner. Besides, the baffle seems not to have influence in the hydrodynamic behaviour. This is attributed to the low velocity in this zone, caused by the limited influence of the impeller. 3.3.2. Main anoxic zone Fig. 9d, Fig. 9e and Fig. 9f show the values of the magnitude of the velocity, the vertical velocity and the flow lines in the main anoxic zone. It can be observed that the velocity magnitude profile (Fig. 9d) in the different XY sections reproduces a flow pattern with high rotationality as a consequence of the action of the impeller. It is also observed that the highest values of the velocity module are found on the sides and in the central part of the reactor. As for the vertical component of the velocity (Fig. 9e), preferential upward flow zones are observed near the reactor walls, with a flow channel on the outside of AnoxAn. The latter is also observed in Fig. 9f, where the flow lines describe an upward flow path formed through the walls. Consequently, most of the massive transport in the main anoxic zone occurs through the reactor walls. As a result, downward flow velocity profiles are found mainly in the central part of the reactor. In addition, the limited influence of the impeller is shown in Fig. 9d-f. Both the magnitude of the velocity and the vertical value of the velocity decrease with the height of the reactor. This coincides with what was observed and discussed for the anoxic-clarification transition zone in Fig. 9a-c. Consequently, the zone of influence of the impeller barely reaches the baffle in the clarification zone, which implies a delay in the homogenization of the concentration of the anoxic and clarification zones, as observed for RTD3 in section 3.2. 3.3.3. Anaerobic-anoxic transition zone Finally, the magnitude of the velocity, vertical velocity and flow lines in the volume around the deflector between the anaerobic and anoxic zones are shown in Fig. 9g, Fig. 9h and Fig. 9i respectively. The deflector, with a width of 4 cm from the reactor walls, is represented in grey. As in the main anoxic zone (Fig. 9d), the velocity magnitude profile (Fig. 9g) also reproduces the rotational flow patterns induced by the impeller. Due to the influence of the deflector, it is observed that the highest values of the velocity modulus are located in the inner part of the reactor while near the walls the velocity magnitude is lower. For the vertical component of the velocity (Fig. 9h), and due to the influence of the deflector, the highest positive value is clearly concentrated in the central part of the reactor, forming preferential upward flow patterns in the inner part of the section. This generation of preferential flow patterns can also be observed at Plascencia-Jatomea et al 2015. Once the influence of the impeller on the anoxic volume has been reached, an upward flow channel is formed through the walls as explained for the main anoxic zone in Fig. 9d-f. As shown in Fig. 9h, there are no vertical velocities around the deflector, confirming the presence of dead flow zones. Similar effect occurs for sectional changes in Climent et al. (2018). The latter is represented in Fig. 9i, where the flow lines avoid the outer corners of the deflector, which improves the behaviour of the stagnant flow. In addition, downward flow can be observed near the internal walls of the deflector (Fig. 9hi), in the anoxic zone. Although these observed velocity values are very small and not very significant in this case, depending on the geometry of the reactor they can contribute to creating a downward flow channel, causing unwanted mixtures or contamination from the anoxic zone to the anaerobic zone. Figure 10. Hydrodynamic scheme of AnoxAn. Experimental conditions: ωimp = 100 rpm, Qin = 10.4 L/h, Qnit = 31.0 L/h (COLOURED) All the main phenomena previously explained in the hydrodynamic analysis are resumed in Fig. 10, where green arrows represent main up flow zones and the red ones down flow zones. Dashed line zones mean potential stagnant or dead flow volumes. The thickness of the dates in the figure indicates the intensity of the flow and therefore the potential mass transport. 3.4. Tracer transport analysis based on CFD simulations Tracer concentration evolution is analysed next, analysing pulse RTD2 and step RTD3 tests for the anoxicclarification transition, main anoxic and anaerobic-anoxic transition zones. RTD analysis combined CFD and velocity field analysis to analyze dead volumes and chanelling have been already used in water treatment field (Brannock et al. 2010a, Climent et al. 2018, Wei et al. 2019). 3.4.1. Anoxic-clarification transition zone Figure 11. Tracer concentration field in anoxic-clarification transition zone for different time steps in RTD2 (a) 5 min, (b) 10 min, (c) 15 min, (d) 20 min, (e) 25 min and (f) 30 min. Experimental conditions: ωimp = 100 rpm, Qin = 10.4 L/h, Qnit = 31.0 L/h, Ctracer = 350 g/L (Pulse tracer test, 3s nitrate stream injection). (COLOURED) Fig. 11 shows the concentration range of the tracer for the different time steps in the RTD2 pulse test. It is observed that the tracer slowly reaches the effluent outlet. It takes between 5 and 10 minutes to reach the clarification zone due to the limited influence of the impeller in the upper part of the anoxic zone. This is consistent with the results obtained in the hydrodynamic analysis of the fluid in the reactor. In addition to this and as the simulation progresses, the presence of a preferential flow pattern towards the outlet is observed. A greater concentration of tracer is found around the outlet. This pattern is developed after 15 minutes of simulation, confirming the existence of a dead volume in the opposite corner of the outlet, as it was shown in Fig. 9a-c. Figure 12. (a) Scheme of tracer concentration measurement points in anoxic-clarification transition zone (dimensions in meters) and (b) Tracer concentration evolution in the outlet (P1) and its opposite corner (P2) for RTD2. Experimental conditions: ωimp = 100 rpm, Qin = 10.4 L/h, Qnit = 31.0 L/h, Ctracer = 350 g/L (Pulse tracer test, 3s nitrate stream injection). Fig. 12 shows the evolution of the concentration of the tracer for two different points at the exit and its opposite corner (P1 and P2). Firstly, it is observed that the maximum concentration value of the tracer at the exit (P1) is slightly higher than that of its opposite corner (P2). The delay between the two peaks is 2.5 minutes (5% of the total RTD2 HRT). Although P1 is located further from the tracer injection point (PN) than P2, the tracer first arrives at the outlet due to the channelling zone observed in Fig.9a-c and 11c-d and similar to Climent et al. (2018). The difference in tracer concentration at both points is less than 5% after 20.5 minutes of experiment (44% of total RTD2 HRT). After 26 minutes of experiment, the concentration of the tracer in the opposite corner (P2) is higher than in the output (P1) for the first time, reaching a maximum of 5% higher at t=66 minutes. This confirms that tracer dilution first happens through the preferential flow formed around the outlet (P1), generating a zone of stagnant behaviour in its opposite corner (P2) as shown in Fig. 9a-c and Fig. 11. 3.4.2. Main anoxic zone Figure 13. Tracer concentration field in main anoxic zone for different time steps in RTD2 (a) 5 min, (b) 10 min, (c) 15 min and (d) 20 min. Experimental conditions: ωimp = 100 rpm, Qin = 10.4 L/h, Qnit = 31.0 L/h, Ctracer = 350 g/L (Pulse tracer test, 3s nitrate stream injection). (COLOURED) Fig. 13 shows the tracer concentration fields for the different time steps in the RTD2 pulse test. It is observed that, first of all, the main tracer transport exists in the zones with higher upper velocity, i.e. the preferential flow channels near the walls. The latter matches with the observed in velocity field analysis for Fig 9d-f. As the simulation time progresses, it is observed that due to a higher vertical velocity, the tracer is firstly diluted in the mentioned outer part of AnoxAn, which is attributed to the action of the impeller. On the other hand, the tracer remains stagnant in the internal part of the system. This is attributed to a very low value of the previously observed flow velocities. Figure 14. (a) Scheme of tracer concentration measurement points in main anoxic zone (dimensions in meters) and (b) Tracer concentration evolution in the central part (P3) and near the walls of the reactor (P4) for RTD2. Experimental conditions: ωimp = 100 rpm, Qin = 10.4 L/h, Qnit = 31.0 L/h, Ctracer = 350 g/L (Pulse tracer test, 3s nitrate stream injection). Fig. 14 shows the evolution of the tracer concentration for the central points and near the wall (P3 and P4). Firstly, it is observed that the maximum tracer concentration value in P4 is 30% higher than in the internal part of the reactor (P3), with the delay in time between the two peaks being 5 minutes (10% of the total RTD2 HRT). These differences in the value and time of the tracer concentration peaks reveal that the mass transport is greater through the walls (P4) than in the central part of AnoxAn (P3), confirming the channelling phenomena observed in Fig. 9d-f. and Fig. 13 and also noticed for Climent et al. (2018). In addition, it is observed that a difference between the tracer concentration of both points is less than 5% after 7.5 minutes (16% of the total RTD2 HRT). After 9 minutes, the concentration of the tracer in the central part (P3) is higher than that of the walls (P4), reaching a maximum of 10% higher in t=58 minutes. This is due to the preferential flow that is formed through the walls. In consequence, dilution occurs first in the outer part of AnoxAn, while in the central part of the reactor a stagnant zone is formed, as also indicated in Fig. 9d-f and Fig. 13 of the analysis. 3.4.3. Anaerobic-anoxic transition zone Figure 15. Tracer concentration field in anaerobic-anoxic transition zone for different time steps in RTD2 (a) 5 min, (b) 10 min, (c) 15 min and (d) 20 min. Experimental conditions: ωimp = 100 rpm, Qin = 10.4 L/h, Qnit = 31.0 L/h, Ctracer = 350 g/L (Pulse tracer test, 3s nitrate stream injection). (COLOURED) The concentration range of the tracer for different time steps in the RTD2 simulated pulse test is shown in Fig. 15. Fig. 15c-d clearly shows the shape of the deflector in the concentration field. The lower concentration profiles are observed in the central part of the reactor, coinciding with the higher velocities of flow rise. In addition, a higher concentration of the tracer is observed in the zones where zero or low velocities are registered, especially in the zone above the deflector. This confirms that dead volumes and stagnant areas are found in these sections with the highest concentration of tracers. Furthermore, it shows that the tracer does not reach the anaerobic zone due to the presence of the deflector, suggesting that this element is crucial to achieve the desired anoxic-anaerobic hydraulic separation as also observed for a different multi-environmental reactor in Calder et al. (2013). Finally, as the time of the experiment progresses, a complete dilution in the reactor is observed. Figure 16. (a) Scheme of tracer concentration measurement points in anaerobic-anoxic transition zone (dimensions in meters) and (b) Tracer concentration evolution in the upper deflector zone (P5) and in the deflector (P6) for RTD2. Experimental conditions: ωimp = 100 rpm, Qin = 10.4 L/h, Qnit = 31.0 L/h, Ctracer = 350 g/L (Pulse tracer test, 3s nitrate stream injection). Fig. 16 shows the evolution of the tracer concentration for two different points near the deflector (P5 and P6). First, it is observed that the maximum tracer concentration value at P6 is 20% higher than near the deflector (P5). However, almost from the beginning of the experiment, the tracer concentration remains at least 5% higher at P5 than at P6 until complete dilution. The latter confirms that a stagnant zone is formed under the influence of the deflector as shown in Fig. 9g-i and Fig. 15. 3.4.4. Overall reactor Figure 17. Tracer concentration field in AnoxAn for different time steps in RTD3 (a) 20 min, (b) 40 min, (c) 80 min and (d) 150 min. Experimental conditions: ωimp = 100 rpm, Qin = 10.4 L/h, Qnit = 31.0 L/h, Ctracer = 10 mg/L (Step tracer test, constant nitrate stream injection). (COLOURED) Levenspiel, O. (1999) Chemical Reaction Engineering. 3rd Edition. John Wiley & Sons, New York, 54. Le Moullec, Y.; Potier, O.; Gentric, C.; Leclerc, J.P. (2008). Flow field and residence time distribution simulation of a cross-flow gas-liquid wastewater treatment reactor using CFD. Chem. Eng. Sci., 63(9), 2436-2449. https://doi.org/10.1016/j.ces.2008.01.029 Liu, M., Yang, M., Chen, M., Yu, D., Zheng, J., Chang, J., Wang, X., Ji, C., Wei, Y. (2018). Numerical optimization of membrane module design and operation for a full-scale MBR by computational fluid dynamics. Bioresour. Technol. 269, 300-308. https://doi.org/10.1016/j.biortech.2018.08.089 Lopez, I.; Borzacconi, L. (2010). UASB reactor hydrodynamics: residence time distribution and proposed modelling tools. Environ. Technol., 31:6, 591-600. https://doi.org/10.1080/09593331003646638 Losada, I.J.; Lara, J.L.; Del Jesus, M. (2016). Modeling the interaction of water waves with porous coastal structures. J.Waterway, Port, Coastal, Ocean Eng., Vol 142(6). https://doi.org/10.1061/(ASCE)WW.19435460.0000361 Martin-Dominguez, A.; Tzatchkov, V.G.; Martin-Dominguez, I.R.; Lawler, D.F. (2005). An enhanced tanks-in-series model for interpretation of tracer tests. J. Water. Supply. Res. T. 54, 435-448. https://doi.org/10.2166/aqua.2005.0041 Martin, K.J.; Nerenberg, R. (2012). The membrane biofilm reactor (MBfR) for water and istewater treatment: principles, applications, and recent developments. Bioresour. Technol. 122, 83–94. https://doi.org/10.1016/j.biortech.2012.02.110 Michalopoulos, I., Kamperidis, T., Seintis, G., Pashos, G., Lytras, C., Papadopulou, K., Boudouvis, A.G., Lyberatos, G. (2018). Experimental and numerical assessment of the hydraulic behaviour of a pilot-scale Periodic Anaerobic Baffled Reactor (PABR). Comp. Chem. Eng. 111, 278-287. https://doi.org/10.1016/j.compchemeng.2018.01.014 Makinia, J.; Rosenwinkel, K.H.; Spering, V. (2006). Comparison of Two Model Concepts for Simulation of Nitrogen Removal at a Full-Scale Biological Nutrient Removal Pilot Plant. J. Environ. Eng., Vol. 132, No. 4, 476-487. http://doi.org/10.1061/(ASCE)0733-9372(2006)132:4(476) Mihovilovic, B. (2010). Optimización del proceso BLAS® en alta carga: Modelización del atascamiento (Optimization of BLAS® process for high loads: Modelling of binding), Ph.D Diss, University of Cantabria, Spain. Nemade, P.D., Dutta, S.M., Shankar, H.S. (2010). Residence time distribution and oxygen transfer in a novel constructed soil filter. J Chem. Technol. Biot. 85, 77–84. https://doi.org/10.1002/jctb.2269 Olivet, D., Valls, J., Gordillo, M.A., Freixo, A., Sánchez, A. (2005). Application of residence time distribution technique to the study of the hydrodynamic behaviour of a full-scale wastewater treatment plant plug-flow bioreactor. J. Chem. Technol. Biotechnol. 80, 425–432. https://doi.org/10.1002/jctb.1201 Pereira, J., Karpinska, A., Gomes, P., Martins, A., Dias, M.&Santos, R. (2011). Activated Sludge Models coupled to CFD simulations In: Single and two-Phase Flows on Chemical and Biomedical Engineering. Bentham Science Publishers, Dubai, U.A.E. Plascencia-Jatomea, R., Almazán-Ruiz, F.J., Gómez J., Rivero, E.P., Monroy, O., Gonzalez, I. (2015). Hydrodynamic study of a novel membrane aerated bioflim reactor (MABR): tracer experiments and CFD simulation. Chem. Eng. Sci. 138, 324-332. https://doi.org/10.1016/j.ces.2015.08.004 Potier, O.; Leclerc, J.P.; Pons, M.N. (2005). Influence of geometrical and operational parameters on the axial dispersion in an aerated channel reactor. Water Res., Vol. 39, pp 4454-4462. https://doi.org/10.1016/j.watres.2005.08.024 Qi, W., Guo, Y., Xue, M., Li,Y. (2013). Hydraulic analysis of an upflow sand filter: tracer experiments, mathematical model and CFD computation. Chem. Eng. Sci. 104, 460–472. https://doi.org/10.1016/j.ces.2013.09.035 Qi, W.K.; Hojo, T.; Li, Y.Y. (2013). Hydraulic characteristics simulation of an innovative self-agitation anaerobic baffled reactor (SA-ABR). Bioresour. Technol., Vol. 136, pp 94-101. https://doi.org/10.1016/j.biortech.2013.02.033 Rehman, U. (2016). Next generation bioreactor models for wastewater treatment systems by means of detailed combined modelling of mixing and biokinetics; Ph.D Diss, BioMath, Ghent University, Belgium. Renade, V.V., (2002).Computational Flow Modeling for Chemical Reactor Engineering. (Academic Press, New York). Sarathai, Y., Koottatep, T., Morel, A. (2010). Hydraulic characteristics of an anaerobic baffled reactor as onsite istewater treatment system. J. Environ. Sci. 22(9), 1319–1326. https://doi.org/10.1016/S10010742(09)60257-6 Seb, B. (2017). Numerical characterisation of a ship propeller. Master Thesis, University of Zagreb. Tang, D.; Jess, A.; Ren, X.; Bluemich, B.; Stapf, S. (2004). Axial dispersion and wall effects in narrow fixed bed reactors: a comparative study based on RTD and NMR measurements. Chem. Eng. Technol. 27(8), 866-873. https://doi.org/10.1002/ceat.200402076 Tejero, I.; Diez, R.; Esteban, A.L.; Lobo, A.; Temprano, J.; Rodríguez, L. (2010). Reactor biológico anóxico-anaerobio para la eliminación de nutrientes de aguas residuales. Spanish Patent ES2338979B2 (in Spanish). Tejero, I.; Eguía, E. (1991). Proceso de tratamiento biológico de aguas residuales basado en biopelícula desarrollada sobre un soporte permeable a gases. Spanish Patent ES2038556 (in Spanish). Terashima, M.; Goel, R.; Komatsu, K.; Yasui, H.; Takahashi, H.; Li, Y.Y.; Noike, T. (2009). CFD simulation of mixing in anaerobic digesters. Bioresour. Technol., Vol 100, 2228-2233. https://doi.org/10.1016/j.biortech.2008.07.069 Trad, Z., Vial, C., Fontaine, J.P., Larroche, C. (2015). Modelling of hydrodynamics and mixing in a submerged membrane bioreactor. Chem. Eng. J. 282, 77-90. https://doi.org/10.1016/j.cej.2015.04.119 Wei. P., Mudde, R.F., Uittewaal, W., Spanjers, H., van Lier J.B., de Kreuk, M. (2019). Characterising the two-phase flow and mixing performance in a gas-mixed anaerobic digester: Importance for scale-up applications. Water Res. 149, 86-97. https://doi.org/10.1016/j.watres.2018.10.077 Weller, H.G.; Tabor, G.; Jasak, H.; Fureby, C. (1998). A tensorial approach to computational continuum mechanics using object oriented techniques. Comput. Phys. 12, 620-631. https://doi.org/10.1063/1.168744 Wicklein, E., Batstone, D., Ducoste, J., Laurent, J., Griborio, A., Wicks, J., Saunders, S., Samstag, R., Potier, O., Nopens, I. (2015). Good modelling practice in applying computational fluid dynamics for WWTP modelling. Water Sci. Technol. 73(5), 969-982. https://doi.org/10.2166/wst.2015.565 Wu, B. (2009). CFD Analysis of mechanical mixing in anaerobic digesters. Structures & Environment Division of ASABE 52 (4), 1371–1382. https://doi.org/10.13031/2013.27786 Yan, X., Xiao, K., Liang, S., Lei, T., Liang, P., Xue, T., Yu, K., Guan, J., Huang, X. (2015). Hydraulic optimization of membrane bioreactor via baffle modification using computational fluid dynamics. Biores. Technol. 175, 633-637. https://doi.org/10.1016/j.biortech.2014.10.133 Yerushalmi, L.; Alimahmoodi, M.; Mulligan, C.N. (2011). Performance evaluation of the BioCAST technology: a new multi-zone istewater treatment system. Water. Sci. Technol. 64(10), 1967-1972. https://doi.org/10.2166/wst.2011.776 Yerushalmi, L., Alimahmoodi, M., Behzadian, F., Mulligan, C.N. (2013). Mixing characteristics and liquid circulation in a new multi-environment bioreactor. Bioprocess Biosyst. Eng. 36, 1339–1352. https://doi.org/10.1007/s00449-012-0836-8 Zeng, Y., Mu, S.J., Lou, S.J., Tartakovsky, B., Guiot, S.R., Wu, P. (2005). Hydraulic modelling and axial dispersion analysis of UASB reactor. Biochem. Eng. J. 25(2), 113–123. https://doi.org/10.1016/j.bej.2005.04.024 Zhang, J., Tejada-Martinez, A.E., Zhang Q. (2016). Rapid analysis of disinfection efficiency through computational fluid dynamics. J. Am. Water Works Assoc. 2018, E50-E59. https://doi.org/10.5942/jawwa.2016.108.0005