Novel dissolved oxygen distribution model of a full-scale aeration tank for a municipal wastewater treatment plant in a cold climate region

a

gas−liquid interfacial area (m²/m³)

AOR

actual oxygen requirement (mol/(m³·s))

AS

activated sludge

d_B

air bubble diameter (m)

BOD

biological oxygen demand (g/m³)

c

dissolved gas concentration in liquid (mol/m³)

c*

equilibrium concentration of dissolved gas in liquid (mol/m³)

CFD

computational fluid dynamics

COD

chemical oxygen demand (g/m³)

D

diffusion coefficient (m²/s)

DO

dissolved oxygen concentration (mg/L or mol/m³)

D_h

hydraulic diameter of a rectangular tank (m)

f_D

drag force (N/m³)

F

additional volume force, e.g., virtual mass force (N/m³)

ϕ

phase volume fraction (m³/m³)

ϕ_g

gas phase volume fraction (m³/m³)

ϕ_l

liquid phase volume fraction (m³/m³)

GHU

gas hold-up (%)

g

gravity (m/s²)

H

molar concentration Henry coefficient (Pa·m³/mol)

I

identity matrix

k_L

liquid side mass transfer coefficient (m/s)

k_La

volumetric mass transfer coefficient (L/s)

L

length of aeration tank (m)

mass transfer rate from the gas to the liquid (kg/(m³·s))

M

molecular weight (kg/mol)

MLSS

mixed liquor suspended solids (kg/m³)

MVLSS

mixed liquor volatile suspended solids (kg/m³)

n

bubble number density (L/m³)

NH₄−N

ammonium nitrogen concentration (mg/L)

NO₃−N

nitrate nitrogen concentration (mg/L)

O₂

oxygen gas

OF

hydraulic overflow over the aeration tank baffle walls (m)

OTR

oxygen transfer rate (mol/(m³·s))

p

pressure (Pa)

p_ref

reference pressure (Pa)

Q

volumetric flow rate (m³/s)

Q_g

volumetric air flow rate (m³/s)

Q_l

volumetric water flow rate (m³/s)

RANS

Reynold's average Navier–Stokes equations

Re

Reynold's number (-)

SCADA

supervisory control and data acquisition

ρ

density (kg/m³)

ρ_g

air density (kg/m³)

ρ_l

water density (kg/m³)

T

temperature (°C)

t

time (s)

u

flow velocity rate (m/s)

u_slip

slip velocity rate (m/s)

water dynamic viscosity (Pa·s)

turbulent viscosity (Pa·s)

W

width of aeration tank (m)

WWTP

wastewater treatment plant

Note(s) The subscripts ‘l’ and ‘g’ denote quantities related to the liquid phase and the gas phase, respectively. Bolded abbreviation refers to a vector.

Aeration can make up 45–75% of the activated sludge (AS) plant's total energy demand (Rosso et al. 2008). In the AS process, aeration of wastewater is needed to provide enough dissolved oxygen (DO) for microorganisms that metabolize organic matter and convert ammonia to nitrate (nitrification). If there is not enough DO available, the microorganisms cannot metabolize organic matter or do nitrification (Åmand et al. 2013).

At the moment, the automated aeration control is usually done based on one DO measurement point (Fan and Boshnakov 2010) or no measurement in the aerated zone but based on an assumed oxygen consumption pattern (Åmand et al. 2013). Automated aeration control keeps the DO levels in the aerated reactors at the desired setpoints. The oxygen requirement of the process is dynamic due to variations in the plant's loading, meaning that to sustain the desired DO setpoint, air flow to the system needs to be regulated (Gray et al. 2011). By using a dynamic, time-dependent, DO distribution model, it is possible to understand the relationships between the plant's loading, aeration equipment, control algorithms, process performance, and energy consumption, thus leading to a significantly more realistic prediction of aeration performance.

The use of computational fluid dynamics (CFD) modelling in wastewater treatment has been steadily increasing as a result of recent advancements in multiphase flow research. By coupling turbulence and multiphase modelling, it becomes possible to simulate the aeration system of the wastewater treatment plant (WWTP) with reasonable accuracy. CFD can predict the DO distribution and multiphase flow pattern in an aeration tank (Karpinska & Bridgeman 2016). When turbulent flows are considered, the modelling of the bubbly flow is based on the Eulerian two-fluid model derived from Reynold's average Navier–Stokes equations (RANS). This model is generally used to simulate bubbly flow reactors with a large number of bubbles (Fayolle et al. 2007,). The resulting equations describe the motion of the two fluids as if they were interpenetrating and interacting continuously, assuming that each element of the spatial domain contains a certain fraction ϕ_l of the continuous (liquid) phase and a fraction ϕ_g = 1 − ϕ_l of the dispersed (gas) phase. Since the resulting equations are formulated in the Eulerian frame of reference, such models are called two-fluid or Euler–Euler models (Sokolichin et al. 2004). The multi-fluid model is the most used multiphase flow model for the large bubble population in an aeration tank (Höhne & Mamedov 2020).

Smaller bubbles increase the oxygen mass transfer by increasing the interfacial gas–liquid bubble surface area. Fine pore membrane diffusion can produce smaller bubbles (Rosso et al. 2008). The flow structure in the bubble reactor is greatly influenced by gas–liquid properties, especially gas flow rate and bubble size (Sommerfeld & Broder 2009). Important multiphase parameters are the gas–liquid interfacial area (a) and the gas hold-up (GHU), which affect the inter-phase mass and momentum transfer. They depend on the physical properties of the liquid and the fluid flow regime. GHU is the volume fraction of gas (ϕ_g) in the total volume of gas–liquid phase in the bubble reactor. It is closely related to the d_B, to the volumetric mass transfer coefficient (k_La) and to the superficial gas velocity or the relative motion between the gas bubbles and the liquid (u_slip) (Saleh et al. 2018). The k_La is a parameter that determines the rate at which a gaseous compound can transfer between the gas phase and the liquid phase. These parameters are affected by the hydrodynamic characteristics, such as the aeration tank dimensions, the liquid height (Sasaki et al. 2017) and the liquid properties (Mouza et al. 2005). An increase in u_slip increases GHU and decreases the d_B, thus increasing the a and the k_La between the gas and liquid. The important forces that have been found to replicate the bubble-fluid interactions are drag, virtual mass, lift, turbulent dispersion, and wall lubrication forces. The viscous drag force (f_D) expresses the resistance experienced by a bubble as it moves relative to the surrounding liquid. The virtual or added mass force accounts for the work done by the bubbles to accelerate the surrounding liquid. A bubble travelling through a shearing flow will experience a lift force perpendicular to the direction of motion (Amol & Joshi 2005). The lift force models the contribution shearing motion on the momentum (Hidman et al. 2022). The turbulent dispersion force models the influence of random turbulent eddies on the flow field from the bubbles. Of the different interaction forces between gas and liquid, pressure and drag force are the most important (Sokolichin et al. 2004).

The aim of this study was to develop a new method for designing and quantifying hydraulic profiles and DO distribution patterns in an aerated tank for minimizing DO gradients. We addressed three different hydraulic flows and DO distribution patterns in a full-scale aeration tank to estimate the effects of variating hydraulic flow, common in cold climate regions. A single, full-scale aeration tank was studied using the CFD model created, supported by experimental data. This study aims to present an experimentally validated simulation predicting the hydraulic flow and the DO distribution in an as-build aeration tank, equipped with a bottom fixed aeration grid consisting of fine bubble disc diffusers.

The aeration tank flow profile, the geometry and the sampling point positions are presented in Figure 2. There is an aeration system on a tank bottom with one aeration air inlet, consisting of a total of 560 membrane fine bubble disc diffusers. See Table 2 for disc diffuser specification.

We made a CFD model for estimating DO distribution of three different hydraulic profile cases in Viikinmäki WWTP's aeration line 5 aerated zone 5. These are represented in this study as hydraulic overflow (OF) over the aeration tank baffle walls: 0.1 m (water depth 12.1 m), 0.35 m (water depth 12.35 m) and 0.6 m (water depth 12.6 m); see Table 1. Calculational domain boundary conditions were (see Figure 2(a)) were:

• Inlet air

• – Liquid boundary condition

• • ○ Air/normal inflow velocity (derived from measured air flow taken from plant SCADA), total air flow was divided evenly to 560 diffusers

• – Gas boundary conditions

• • ○ Air mass flux (derived value from measured air flow), total mass flux was divided evenly to 560 diffusers

  • ○ Air bubble density flux (derived value from total air flow and air bubble diameter), total bubble density flux was divided evenly to 560 diffusers

• Inlet water

• – Liquid boundary condition

• • ○ Normal inflow velocity (from plant SCADA), flow from upstream zone 4 entering zone 5 over cross-sectional baffle wall top and bottom openings area

• – Gas boundary condition

• • ○ No gas flux

• Inlet DO

• – DO (measured value, sampling point 1), flow from upstream zone 4, entering zone 5 over baffle wall top and bottom opening

• Outlet air

• – Liquid boundary condition

• • ○ Atmospheric pressure at the water surface (measured value)

• – Gas boundary conditions (gas outlet)

• • ○ Effective gas (O₂) density (derived value from measured air flow and from a separate mass balance calculation)

  • ○ Air bubble diameter and air density at the water surface

• Outlet water

• – Liquid boundary condition

• • ○ Normal outflow velocity (from plant SCADA), flow to downstream zone 6 leaving zone 5 over cross-sectional baffle wall top and bottom openings area

• – Gas boundary condition

• • ○ No gas flux

• Outlet DO

• – Outflow by convection

• Tank walls and bottom

• – No slip condition

The parameter values from Table 1 were used as inputs for simulations. It took a total of 128 calibration simulation runs for calibrating k_L and d_B. Calibration simulations were run until the simulated DO level in sampling point 2 matched with the measured level within simulation tolerance. The d_B was calibrated by changing its value by 0.5 mm intervals in a range of 0.5–3.5 mm. This range was selected based on disc diffuser specifications (Table 2). The k_L was calibrated by changing its value by ±50% from the initial value depending on whether the simulated DO level in sampling point 2 was over or underestimated. The initial value for k_L was taken from the literature (Mohan et al. 2021). The liquid side oxygen mass transfer coefficient was calibrated to k_L = 0.7 × 10⁻⁵ m/s. This k_L gives k_La = 2.1 × 10⁻⁴ L/s (0.75 L/h, 0.013 L/min) when a = 30 m²/m³ (Equation (12)) and is in the range of typical values for bubble columns (Wilhelms & Martin 1992; k_La = 0.72 L/h, Mohan et al. 2021; k_La = 0.013 L/min). The initial bubble size was calibrated to d_B = 0.0025 m. This diameter is typical for a fine bubble diffuser type used in this study, representing a typical mean diameter found in literature (Painmanakul et al. 2004; d_B = 0.001– 0.0035 m, Hasanen et al. 2006; d_B = 0.002–0.0035 m).

In addition to the three OF cases, the simulation for the OF 0.6 m case was separately validated with 4 March–24 May 2024 measured DO data collected from sampling points 1, 2, 3, 4, 5 and 6.

The aerated zone 5 average Reynold's number was Re = 23,042 (Equation (1)), highlighting the condition of the turbulent flow regime. Average measured and simulated DO levels and validation errors are presented in Table 4.

There is a good agreement between the simulated and the measured average volume DO. Validation errors were all less than 15% and in the range of 0.1– 0.3 mg/L. Measured DO and simulated DO sampling point correlations were 99% in lateral and 78% in vertical directions, respectively.

For sampling point 2, the measured and simulated DOs are very close, indicating good model performance. At the tank inlet, the simulated DO was overestimated and at the tank outlet, the simulated DO was underestimated. This could be caused by the high turbulence at the tank inlet and outlet. It was difficult to get reliable and repeatable measurements from sampling points 1 and 3. Going deeper in the middle of the tank (middle sampling points 4, 5 and 6) increases validation error. The simulated DOs were all overestimated. One possible reason for this is the assumption that the same DO level from the previous zone 4 enters zone 5 over cross-sectional baffle wall top and bottom openings area. This assumption might not be correct because the trend of the measured vertical DO profile from the tank bottom towards the surface is increasing, whereas the trend of the simulated vertical DO profile is flat. The vertical DO profile will be further studied in the coming studies. We assumed uniform oxygen uptake and uniform biomass distribution throughout the aeration tank. This uniformity assumption was supported by an MLSS study done in Viikinmäki WWTP aeration line 5 aerated zone 6 in 2020. Zone 6 has the same water depth, dimensions, baffle walls and similar disc diffuser layout as zone 5.

One reason for the validation errors is assuming that there is no bubble coalescence or break-up. This is a very significant assumption both for the modelling of the turbulent flow and for the O₂ mass transfer. This can be considered a major limitation of the model. In this model, bubbles can expand or shrink but not completely vanish, merge or split. It should be noted that in real cases, there is always a heterogeneous distribution of bubble sizes in bioreactors due to bubble coalescence and break-up. Another limitation is that only two phases (water and air bubbles) are included but the solid phase (biomass) is excluded. In real cases, there are three interacting phases, but including solids in this model would have increased the model complexity and calculational cost too much considering our resources.

The average u_slip and GHU in Table 5 are both in the range of typical values for aeration tanks for all three OF cases (Clift et al. 1978, Gresch et al. 2011).

Figure 6 shows how water circulates counterclockwise in a big vortex in the middle of the tank. Another, smaller vortex can be seen circulating clockwise on the lower right-hand side of the tank. These vortices are caused by the upward gas flow from the bottom diffusers and the hydraulic flow entering (left) and leaving (right) the tank. Matko et al. (2021) noticed that the upward flow from the bottom diffuser creates a local re-circulating radial flow pattern and that the maximum velocities are near the water surface. Non-homogeneous water velocity distribution can be observed over the aeration between the two vortices due to the chaotic motion of the air bubbles. A similar phenomenon was observed by Höhne & Mamedo (2020). The bigger vortex pushes oxygen-rich water to the surface and towards the beginning of the tank. There is a U-shaped short-circuiting going from the baffle wall top inlet down to the tank bottom, over the aeration grid, back to the surface and further to the baffle wall top outlet. The big vortex and short-circuiting probably decrease the oxygen transfer efficiency in zone 5.

The overall predicted water flow patterns are quite similar for all three OF cases, but there are noticeable variations in DO distribution patterns. Since the air flow rate from diffusers to zone 5 and the DO from the previous zone 4 are the same for all three OF cases, one possible reason for variations in DO distribution patterns could be the local differences in turbulent flow velocities. This seems to be caused by the three different water inlet and outlet velocities, 0.44, 0.32 and 0.22 m/s, to and from zone 5 (Table 1). The possible reason for the similar water flow patterns for all three OF cases could be a combination of a deep basin, dense aeration grid and baffle wall structure. Since similar combinations can be found in other Viikinmäki WWTP aerated zones, the results can be adapted to these too. In general, utilizing the results of this study in other aerated tanks would require a similar deep tank, dense aeration grid and baffle wall structure.

When comparing the average flow velocity rates in Table 5 for the OF 0.1 m and OF 0.35 m cases, it can be seen that for the OF 0.1 m case u_l is 6.0% higher, u_g is 3.0% higher, u_l range size is 3.7% higher and u_g range size is 7.8% lower, indicating that there is more turbulent flow in the tank for the OF 0.1 m case. In addition, the average DO is 5.9% lower and the DO range is 27.6% lower for the OF 0.1 m case, indicating that the DO is more evenly distributed in the volume.

When comparing the average flow velocity rates in Table 5 for the OF 0.35 and 0.6 m cases, it can be seen that for the OF 0.35 m case, u_l is 5.8% higher, u_g is 1.8% higher, u_l range size is 5.2% higher and u_g range size is 14.3% higher, indicating that there is more turbulent flow in the tank for the OF 0.35 m case. In addition, the average DO is 6.4% lower and the DO range is 15.2% lower for the OF 0.35 m case, indicating that the DO is more evenly distributed in the volume.

To summarize, lower OF creating more turbulence in the volume could be more effective levelling the DO gradients in zone 5. This was also observed by Shah et al. (2024) in a full-scale bottom-aerated tank equipped with disc diffusers. In addition, increased air flow improved the mixing and DO distribution became more uniform (Shah et al. 2024).

The biggest difference in simulated average volume DO was 0.3 mg/L between the OF 0.6 m to OF 0.1 m cases (Table 5). The case OF 0.1 m average DO 2.1 mg/L was a sufficient DO source for the estimated DO sink. Since the DO setpoint value was 3.0 mg/L during the measurement campaign, we estimated that a 10% lower DO setpoint value, 2.7 mg/L, could be sufficient for the nitrification process.

Based on a separate oxygen mass balance calculation, this would mean using 3% less air for the same treatment result, thus lowering aeration energy cost. For the aeration process dimensioning details, we refer the reader to a book on biological wastewater treatment (e.g., Grady et al. 2011). With more even DO distribution, the DO setpoint value could be even lower.

Because it is very difficult to control Viikinmäki WWTP's incoming hydraulic flow, one way of levelling DO gradients more effectively in a single aerated zone could be to modify the water inflow and outflow rates over and through the baffle walls. This could be done by changing the current baffle wall height, opening size and opening positioning (Table 1, Figures 2(a) and 3(a)). Another way could be to modify the existing aeration grid layout so that it is divided into two separate aeration groups, both with their own aeration air feeds and control valves. This would enable feeding different air flow to the 50% of tank bottom area. Then the big turbulent vortex seen in Figure 6 could be modified to increase the hydraulic retention time. Both of the modifications here presented should be CFD modelled first to see their effect on levelling the DO gradients.

The aim of this study was to develop a new method for designing and quantifying hydraulic profiles and DO distribution patterns in an aerated tank for minimizing DO gradients. We addressed three different hydraulic flows and DO distribution patterns in a full-scale aeration tank to estimate the effects of variating hydraulic flow, common in cold climate regions.

In addition to the three OF cases, the simulation for the OF 0.6 m case was separately validated with 4 March–24 May 2024 measured DO data collected from sampling points 1, 2, 3, 4, 5 and 6 (Figure 2(b)). There is a good agreement between the simulated and the measured average volume DO. Validation errors were all less than 15% and in the range of 0.1–0.3 mg/L. Measured DO and simulated DO sampling point correlations were 99% in lateral and 78% in vertical directions, respectively.

Based on the observed differences in DO, we estimated that a 10% lower DO setpoint could be applied. This would mean using 3% less air for the same treatment result, thus lowering aeration energy cost. With more even DO distribution, the DO setpoint value could be even lower. The results show that variating hydraulic flow causes noticeable DO level variation in zone 5. Lower OF, creating more turbulence in the volume, could be more effective levelling the DO gradients in zone 5.

This study presents a novel CFD DO distribution three-dimensional model for a full-scale, deep, bottom-aerated tank with an as-build layout consisting of disc diffusers. The CFD model predicts the hydraulic flow and the DO distribution pattern in a deep tank with dense aeration grid. The aerated tank is divided from the aeration line with two baffle walls. The results provide a numerical method to design and quantify aeration in a more adaptive and energy-efficient way.

The authors would like to thank Maa-ja Vesitekniikan Tuki ry (MVTT) for their financial support and Helsinki Region Environmental Services HSY for providing the data for this study.

All relevant data are included in the paper or its Supplementary Information.

The authors declare there is no conflict.