In this work, the integration of dynamic bioenergetic calculations in the IWA Anaerobic Digestion Model No. 1 (ADM1) is presented. The impact of bioenergetics on kinetics was addressed via two different approaches: a thermodynamic-based inhibition function and variable microbial growth yields based on dynamic Gibbs free energy calculations. The dynamic bioenergetic calculations indicate that the standard ADM1 predicts positive reaction rates under thermodynamically unfeasible conditions. The dissolved hydrogen inhibition approach used in ADM1 is, however, deemed as adequate, offering the trade-off of not requiring dynamic bioenergetics computation despite the need of hydrogen inhibition parameters. Simulations of the model with bioenergetics showed the low amount of energy available in butyrate and propionate oxidation, suggesting that microbial growth on these substrates must be very limited or occur via alternative mechanisms rather than dissolved hydrogen.
Anaerobic treatment of wastewater is a well-established technology (McCarty 2001). Traditional anaerobic digestion (AD) is used in wastewater treatment plants for waste stabilisation (decrease in organic content) and for methane production (Donoso-Bravo et al. 2011; Kleerebezem et al. 2015). Methane can be used as an energy source once the biogas is upgraded to biomethane. Biomethane, however, is uncompetitive compared to natural gas, unless it is subsidised (Kleerebezem et al. 2015). Therefore, products from the AD process with higher economic value (compared to methane) have been targeted. Examples of those alternative products are hydrogen, bioplastics (Kleerebezem & van Loosdrecht 2007) or carboxylic acids (Agler et al. 2011).
AD modelling has been a very active area of previous work and especially after the publication of the IWA Anaerobic Digestion Model No. 1 (ADM1) (Batstone et al. 2002) – by far the most widely applied and studied AD model to date. The modelling structures used in most of the existing AD models and particularly in the ADM1 typically involve a large set of parameters. Several of these parameters have limited mechanistic interpretation or actual importance (Rodríguez et al. 2006). Models such as ADM1 and alternative AD models are also highly centred on the kinetics, a legacy of previous wastewater treatment models such as the ASM (activated sludge model) series (Henze et al. 2000).
In AD, however, reactions do run close to thermodynamic equilibrium as opposed to aerobic systems. Therefore, bioenergetic factors need to be taken into consideration in anaerobic processes (Kleerebezem & Stams 2000; Rodríguez et al. 2008). Previous studies on the bioenergetics of microbial reactions in anaerobic environments recognise the very low metabolic energy available for propionate and butyrate oxidation, forcing these organisms to live under syntrophic associations (Schink 1997; Junicke et al. 2016). Product formation in fermentative environments is affected by environmental conditions such as pH and temperature (Zoetemeyer et al. 1982a, 1982b). Models based on bioenergetics were developed to predict the variable product formation in these fermentative environments (Rodríguez et al. 2006; González-Cabaleiro et al. 2015; Regueira et al. 2018). Kinetic rates also appear to be affected by thermodynamics when reactions run close to equilibrium (Rottenberg 1973; Westerhoff et al. 1982; Hoh & Cord-Ruwisch 1996).
In the ADM1 model (Batstone et al. 2002), bioenergetics are not explicitly incorporated. The thermodynamic limitations on long chain fatty acids (LCFA) and volatile fatty acid (VFA) oxidation are modelled by proxy via an empirical inhibition function based on dissolved hydrogen concentration. This approach requires an extra inhibition parameter (KI,h2) for each reaction that contains hydrogen as a product, as it is the main contributor for the inhibition of the reaction. This approach currently used in ADM1 may lead to scenarios in which reactions have a positive rate under thermodynamically unfavourable conditions (Kleerebezem & van Loosdrecht 2006), therefore breaking the second law of thermodynamics.
Additionally, the ADM1 model relies on the use of fixed yields obtained from a large compilation of data (Pavlostathis & Giraldo-Gomez 1991; Batstone et al. 2002). Microbial growth yields, however, can be correlated via bioenergetic calculations (von Stockar et al. 2008) and do vary as a function of the catabolic energy available (Pirt 1987). The most important correlation methods to estimate microbial growth yields are the TEEM (thermodynamic electron equivalents model) (McCarty 1965, 2007; Rittmann & McCarty 2001) and the Gibbs energy dissipation method (Heijnen & Van Dijken 1992). The estimation of yields in anaerobic environments, however, appears as not sufficiently accurate to date (Liu et al. 2007).
In this work, the effect of bioenergetics on microbial reactions is incorporated into the ADM1 model via two approaches, namely (i) ADM1 with thermodynamic-based inhibition function and (ii) ADM1 with microbial growth yields as a function of the Gibbs free energy available dynamically calculated. The results obtained from each implementation are compared against the simulation results obtained with the standard ADM1. A more mechanistic model with more interpretable results regarding reaction inhibitions is sought through the integration of bioenergetics in the ADM1 model. A secondary objective is to decrease the large number of parameters required to be calibrated by correlating parameters such as inhibition factors and/or microbial growth yield to the metabolic energy available for each microbial functional group.
Aspects of ADM1 implementation
The ADM1 model (Batstone et al. 2002) was implemented in an Excel/Simulink framework described by Rodríguez et al. (2009). Molar concentration units were used to ensure that mass balances are closed (Kleerebezem & van Loosdrecht 2006) and to facilitate the dynamic computation of both pH and thermodynamics at each time step. The reference parameter set for ADM1 obtained from the BSM2 (Benchmark Simulation Model no. 2) report (Rosen & Jeppsson 2006) was used in all simulations.
ADM1 with thermodynamic inhibition
In the standard ADM1, the effect of hydrogen on VFA oxidation is modelled via an apparent inhibition by the dissolved hydrogen concentration. In the ADM1-Ith implementation, the apparent inhibition term is substituted by the aforementioned thermodynamic inhibition term (Equation (3)). A comparison of the two kinetic models for propionate oxidation used is shown in Table 1.
|Dissolved hydrogen inhibition (ADM1) .||Thermodynamic-based inhibition (ADM1-Ith) .|
|Dissolved hydrogen inhibition (ADM1) .||Thermodynamic-based inhibition (ADM1-Ith) .|
The difference between the two inhibition functions for a microorganism oxidising propionate to acetate is shown in Figure 1.
A clear advantage of the thermodynamic-based inhibition approach is that no parameter is required, as opposed to the apparent inhibition function (for which different KI,H2 are needed for each microbial functional group). The thermodynamic-based inhibition also provides a mechanistic consistency as the model only allows for reactions that are thermodynamically favourable to run. In Figure 1, the apparent inhibition function taking positive values even when the ΔG of the reactions is positive is shown (as previously reported by Kleerebezem & van Loosdrecht (2006)). By applying the thermodynamic-based inhibition function, this error can be avoided.
ADM1 with variable yields
In this alternative implementation of the ADM1 model (ADM1-VY), microbial growth yields (YXS) are dynamically calculated as a function of the available free energy in the catabolic and anabolic reaction. The Gibbs energy dissipation method (Heijnen & Van Dijken 1992) was used to calculate the required dissipated energy (ΔGDiss) per mol of biomass carbon for each carbon source. This energy needs to be met, providing a value for the microbial growth yield function only of the bioenergetics available at each time step. The assumed biomass composition for the calculation of the dynamic microbial growth yields was that reported in ADM1 (CH1.4O0.4N0.2).
Experimental scenario for model evaluations
A highly dynamic AD experiment was selected to conduct the comparative simulations between the standard ADM1 model and the two alternative implementations (ADM1-Ith and ADM1-VY) considering bioenergetic calculations (as described in the two previous sections). The experimental scenario selected was for a laboratory-scale 2 L reactor fed with distillery wastewater (Zaher et al. 2004).
The organic loading rate (OLR) profile of the reactor is shown in Figure 2.
The initial biomass composition for each one of three simulation sets (standard ADM1, ADM1-Ith and AMD1-VY) was estimated by running a simulation in steady state. An influent composition was selected at a time when the reactor was considered to be operationally stable and representative of the biomass of the reactor (t = 2,250 h). The water inflow composition used to simulate the steady state is shown in Table 2.
|Influent composition at t = 2,250 h .|
|Influent composition at t = 2,250 h .|
After the estimation of the initial biomass composition, dynamic simulations were run for all three ADM1 variants: standard ADM1, ADM1 with thermodynamic-based inhibition (ADM1-Ith) and ADM1 with variable yields calculated (ADM1-VY). The model results were plotted against experimental data and analysed. It is important to note that no parameter fit to the data was conducted, since the objective of this work was to compare between the different ADM1 simulations. The complete simulation results and experimental data are provided in the Supplementary Material (available with the online version of this paper).
RESULTS AND DISCUSSION
ADM1 with thermodynamic inhibition
The bioenergetics for all ADM1 reactions are shown dynamically in Figure 3(a) for the standard ADM1 and Figure 3(b) for the ADM1 with thermodynamic inhibition (ADM1-Ith). Sugars and amino acids degradation and acetoclastic methanogenesis are reactions that run far from equilibrium with values lower than −240, −80 and −20 kJ/mol respectively. Therefore, they are unlikely to be ever limited thermodynamically. The rest of the reactions, however, show Gibbs free energy values closer to equilibrium and even become unfavourable at some stages during the simulation (e.g. at time circa 3,000 h). Those reactions that run close to equilibrium have hydrogen either as a substrate (hydrogenotrophic methanogenesis) or as a product (LCFA and VFA oxidation). The bioenergetics of the model in dynamic conditions highlight the key role of the hydrogen concentration on the thermodynamic limitation of the reactions.
The catabolic rates of VFA oxidation are shown in Figure 3(c) and 3(d). Impact of the bioenergetics on the ADM1-Ith is observed in Figure 3(d), where reactions are stopped at positive ΔGR. The standard ADM1, however, predicted at some stages (shaded in Figure 3(a) and 3(c) at t = 1,700–2,500 and 3,100–4,000 h) the oxidation of VFAs at thermodynamically unfeasible conditions, as previously reported by Kleerebezem & van Loosdrecht (2006). Both looking at the simulated bioenergetics (Figure 3(a) and 3(b)) and at the experimental data in Figure 4, it appears that butyrate oxidation was indeed not feasible and propionate oxidation is somehow very difficult to proceed under the current pathways used. These observations suggest that a reconsideration and evaluation of the oxidation pathways for both butyrate and propionate may be necessary.
The key simulation results for the proposed case study are presented in Figure 4. The ADM1 with a thermodynamic-based inhibition (ADM1-Ith) was compared against the standard ADM1 (with the apparent hydrogen inhibition function).
No major differences were observed in gas flow (data not shown) or pH simulation between the two model implementations. However, differences appeared to be more significant for VFA concentrations. The thermodynamic-based model was indeed able to capture in a similar way to ADM1 the concentration of propionate, although still overestimating to a lower degree butyrate concentrations. If temperature corrections are not incorporated, the model is not able to predict accurately VFA concentrations (data not shown). The need for those accurate calculations with temperature effect is reasonable thermodynamically considering that the entropy of VFA oxidation reactions increases. Therefore, by increasing the temperature, the energy available in the reaction will increase. It therefore becomes clear that if thermodynamic calculations are to be used to inhibit reactions, the effect of temperature must be included.
The results shown in Figure 4 validate the dissolved hydrogen inhibition factor approach used in the standard ADM1. It reproduces the actual thermodynamic inhibition quite accurately at typical AD process concentrations in a simple manner. Although it requires H2 inhibition parameters it has a big trade-off for not requiring the dynamic computation of Gibbs energies of each reaction.
ADM1 with variable microbial growth yields
The results of a dynamic simulation with ADM1-VY for key variables measured (pH, butyrate, propionate and acetate) are also shown in Figure 4. Moderate agreement is observed in Figure 4 between all models for both pH and acetate concentrations. However, simulations strongly differ for propionate and especially for butyrate concentrations. Propionate concentrations show an increase between t = 750–1,500 h that is not observed in the experimental data nor in any of the previous two ADM1 implementations. This is caused by the much lower initial population of VFA oxidisers due to the much lower yields predicted (compared to the default values in ADM1) using bioenergetic correlations for those microbial groups (Patón & Rodríguez 2019). In the case of butyrate, ADM1-VY widely overestimated its concentration. Calculated energetics in the simulated data show that the reaction should not proceed, suggesting again (as in the comparison with the standard ADM1 with the ADM1-Ith) that butyrate oxidation might occur through an alternative mechanism, such as a coupled syntrophic reaction where hydrogen is not an endproduct of the reaction or via conductive pili as observed with Geobacter microorganisms (Reguera et al. 2005; Shrestha & Rotaru 2014).
The dynamic simulation of microbial growth yields (Figure 5(b)) shows that the already low initial microbial growth yield values for VFA oxidisers become even lower, to the extent of not allowing microbial growth for some periods at reactor conditions. This occurs due to the observed low catabolic energies (Figure 5(a)) (−20 kJ/mol or more positive) and even becoming positive at some stages (e.g. time circa 2,850 h).
These findings suggest that for the experimentally observed VFA oxidation (see Figure 4) to proceed in the given conditions, alternative microbial functional groupings for VFA oxidisers involving wider metabolic diversities may be required to yield sufficient energy to sustain their population at reactor conditions. An alternative approach might consist of the inclusion of formate as an alternate electron transfer species that can enhance reactions proceeding forward and/or reduce hydrogen concentration in the liquid. Formate, however, has been reported to be thermodynamically equivalent to hydrogen when accounting for reaction bioenergetics (Batstone et al. 2006). Other possibilities include reported syntrophisms between butyrate or propionate oxidisers and methanogenic hydrogenotrophs (Schink 1997; de Bok et al. 2004; Stams & Plugge 2009). Syntrophic electron transfer could even occur via conductive pili (Reguera et al. 2005; Shrestha & Rotaru 2014).
The main conclusions of this work are as follows.
The hydrogen inhibition factor approach used in ADM1 to mimic the thermodynamic limitation of VFA oxidations appears as adequate. Equivalent dynamic results are obtained to those using a thermodynamic-based inhibition. A trade-off between the two approaches exists in the sense that a hydrogen inhibition parameter is required to be defined in exchange for not requiring the dynamic computation of Gibbs free energies.
The low amount of catabolic energy available for VFA oxidisers theoretically predicts that their sustained growth is very unlikely. The inconsistency between commonly used growth yields and energy available suggests that microbial growth on butyrate and propionate might involve an alternative mechanism for interspecies electron transfer (such as conductive pili or other species for electron transfer). A thorough evaluation of the catabolic pathway of the oxidation of syntrophic organisms is recommended.
This publication is based upon work supported by Khalifa University's Award No. CIRA-2018-84 and the Government of Abu Dhabi.