Abstract

The impact on the prediction of key process variables in anaerobic digestion (AD) when activity corrections are neglected (e.g. when ideal solution is assumed) is evaluated in this paper. The magnitude of deviations incurred in key variables was quantified using a generalised physicochemistry modelling framework that incorporates activity corrections. Deviations incurred on the intermediate and partial alkalinity ratio (a key control variable in AD) already reach values over 20% in typical AD scenarios at low ionic strengths. Deviations of moderate importance (∼5%) in free ammonia, hydrogen sulfide inhibition, as well as in the biogas composition, were observed. Those errors become very large for components involving multiple deprotonations, such as inorganic phosphorus, and their magnitude (∼40%) would impede proper precipitation modelling. A dynamic AD case simulation involving a series of overloads showed model underpredictions of the process acidification when activity corrections are neglected. This compromises control actions based on such models. Based on these results, a systematic incorporation of activity corrections in AD models is strongly recommended. This will prevent model overfitting to observations related to inaccurate physicochemistry modelling, at a marginal computational cost. Alternatives for these implementations are also discussed.

INTRODUCTION

Bioprocesses are characterised by (bio)chemical reactions that occur at different rates and with components in different phases (solid, liquid and gas) (Batstone et al. 2002). Bioprocess models typically describe microbial activity aspects meticulously but they often assume a low impact of the physicochemistry compared to the biological activity (González-Cabaleiro 2015).

Changes in physicochemistry can lead to major impacts in biological reactions in anaerobic digestion (AD) (Batstone 2009). For example, a pH drop leads to an accumulation of fully protonated fatty acids and results in microbial inhibition. Other common impact is a pH increase that leads to higher concentration of inhibitory free ammonia (Yenigün & Demirel 2013). pH is also strongly tied up with ionic precipitation (Li & Zhao 2001), which is not included in the widely used Anaerobic Digestion Model No. 1 (ADM1) because of its complexity (Batstone et al. 2002).

The accurate description of any physical or chemical reaction (and the ability to quantify its overall impact on the system), is a requirement of any mathematical model aiming at a reliable prediction of the behaviour of a system. An inaccurate physicochemistry modelling can lead to wasted efforts on calibration of microbial-related parameters such that erroneous conclusions from wrongly calibrated parameters are reached.

Physicochemistry reactions are well understood and can be accurately described. Several models are also available for activity corrections: Debye-Hückel (Holtzer 1954), the Davies equation (Davies 1962), the B-Dot model (Helgeson 1969) or Pitzer equations (Pitzer 1987).

Those models are included in various geochemistry modelling packages such as PHREEQC (Parkhurst & Appelo 2013) or MINTEQA2 (Allison et al. 1991), among others. However, a limitation is that those packages need to be coupled with other software(s) in order to run dynamic simulations. Another limitation when coupling external packages to compute physicochemistry is the large size of their database, which decreases simulation speed significantly (Vaneeckhaute et al. 2018).

Activity corrections have been neglected in the past by most AD models, such as those by Angelidaki et al. (1999), Siegrist et al. (2002) or ADM1 (Batstone et al. 2002). Ionic strength effect in AD models was initially incorporated by adjusting parameters such as Henry constants, gas–liquid transfer rate or biomass yields (Sötemann et al. 2005) or by correcting ionisation constants (ka) (Smith & Chen 2006; Nielsen et al. 2008).

Over the last few years, coordinated works started to address the need for a consensual, accurate and generalised framework to describe physicochemistry in biochemical models (Batstone et al. 2012). That led to several recent modelling efforts that describe physicochemistry more accurately in biochemical models. Those works were described in literature such as: (i) a framework for modelling ion speciation and ion-pairing in plant-wide modelling (Flores-Alsina et al. 2015), (ii) the evaluation of the effects of ionic strength and ion-pairing in AD in plant-wide modelling (Solon et al. 2015) and (iii) a framework with a systematic evaluation for precipitation modelling (Kazadi Mbamba et al. 2015; Kazadi Mbamba et al. 2016). Recent modelling efforts resulted in coupling ADM1 with PHREEQC (Huber et al. 2017).

Process control variables for AD include the so-called α-value (also known as ratio between intermediate alkalinity and partial alkalinity or IA/PA ratio). An anaerobic digester is believed to be operating in a stable manner as long as α <0.3 (Ripley et al. 1986). Measurements of total and partial alkalinity thresholds were defined also at 4.3 and 5.75 in a rather arbitrarily manner (Jenkins et al. 1983). Examples of the use of alkalinity in AD modelling studies include: (i) a model-based controller based on alkalinity using the ratio between volatile fatty acids (VFA) and total alkalinity (Batstone & Steyer 2007) and (ii) a simulation of the ratio VFA/Alkalinity (Schoen et al. 2009). Other variables that affect the AD process are the concentration of free ammonia (which is inhibitory for acetoclastic methanogens) and H2S (inhibitory for VFA degraders and methanogens). Phosphorus precipitation modelling has been of importance due to the value of struvite as a fertiliser and the need to model phosphorus in a plant-wide model context accurately (Flores-Alsina et al. 2016).

In this work, the impact of activity corrections on the accuracy of physicochemistry modelling in AD systems is assessed. That impact is evaluated in terms of key process variables. Also, a number of approaches to implement activity correction are presented and discussed.

MATERIALS AND METHODS

In this section, the generalised algorithm used for calculating pH and activity correction is described, along with the model selected for activity correction. The variables of interest and the computation of its errors are also presented.

Generalised solution algorithm for equilibrium ion speciation and pH calculation

Acid-base and hydration reactions can be assumed to occur instantaneously with respect to the biological reactions taking place in any bioprocesses. Also, some acid base equilibrium reactions, such as for inorganic carbon, involve a hydration and multiple deprotonations (Figure 1).

Figure 1

Equilibrium reactions for inorganic carbon in liquid solution.

Figure 1

Equilibrium reactions for inorganic carbon in liquid solution.

Analogous to ADM1 (Batstone et al. 2002), each dynamic variable of the model is considered as the total sum of concentrations of the chemical species in the same phase. As an example, the chemical species of inorganic carbon are shown in Table 1.

Table 1

Different species of total inorganic carbon state

StateSpecies
Total inorganic carbon (liquid phase) Not hydrated Fully protonated 1st deprotonation 2nd deprotonation 
 CO2 H2CO3 HCO3 CO32− 
StateSpecies
Total inorganic carbon (liquid phase) Not hydrated Fully protonated 1st deprotonation 2nd deprotonation 
 CO2 H2CO3 HCO3 CO32− 
Algebraic solution algorithms for pH and ionic equilibria have been widely applied in recent models (Batstone et al. 2002; Rosen et al. 2006), in which charge neutrality is always assumed (Equation (1)). This approach allows to reduce the number of dynamic variables of the model and thereby, the number of differential equations.  
formula
(1)
where Ci is the concentration and zi the charge for each chemical species.

Activity correction coefficient and generalised procedure to calculate pH, speciation and ionic strength

As described in the introduction, several activity correction models can be used, finding among them the Debye-Hückel equation, the Davies equation and the B-Dot model. The Debye-Hückel equation was not considered due to the limitations of application to ionic strengths solutions above 0.1 molal (Stumm & Morgan 1996; Bethke 2007). B-Dot model (or Pitzer equations) may provide more accurate options but at the expense of complexity with its associated higher computational cost and additional parameters. Therefore, the activity coefficients (γ) were calculated using Davies equation (Equation (2)), which is suitable in AD environments (Batstone et al. 2012) and applicable up to ionic strengths of 0.5 molal (Stumm & Morgan 1996).  
formula
(2)
where , ɛ is the dielectric constant of water corrected by temperature (Malmberg & Maryott 1956), zi is the charge component and IC is the ionic strength of the solution.
The generalised algorithm used to solve the overall physicochemistry is shown in Figure 2. Generalised equations for any number of deprotonations were used to calculate the activity of each species in solution. Equilibrium constants were corrected by temperature using the van 't Hoff equation (Equation (3)):  
formula
(3)
where ΔH0 is the enthalpy of formation of the species at 25 °C.
Figure 2

Generalised ionic speciation and pH solution procedure including ionic strength and its effect on chemical activity (adapted from González-Cabaleiro (2015)).

Figure 2

Generalised ionic speciation and pH solution procedure including ionic strength and its effect on chemical activity (adapted from González-Cabaleiro (2015)).

The procedure implies a double loop to converge both pH and ionic strength (Figure 2). We find this approach sufficiently general to allow it for flexibility, as it can be used for multiple states subject to speciation. This solution procedure was run using an Excel/MATLAB framework described by Rodríguez et al. (2009). To validate this approach, the results obtained with this solution procedure were compared with those obtained with VMINTEQ (Gustafsson 2014) (See Supplementary material, available with the online version of this paper).

Selected physicochemistry dependent variables of interest in AD modelling

Total, partial and intermediate alkalinities (TA, PA and IA), including their ratio (IA/PA), carbon dioxide in biogas (PCO2), the free ammonia concentration (SNH3), hydrogen sulphide (SH2S) and phosphate species concentrations (SH2PO4, SHPO42, SPO43) are key AD variables that depend highly on pH and ionic speciation. Alkalinity and biogas composition values are routinely used to monitor AD process operation while the free ammonia is an important inhibitor of acetoclastic methanogenesis. H2S is a corrosive gas that can also cause inhibition of VFA degradation in the AD process. Inorganic phosphorus (IP) is a variable of interest due to its role in precipitation and the interest of recovering struvite as a fertiliser. Some of the potential precipitated salts with phosphorus are summarised in Table 2.

Table 2

Possible precipitated salts containing inorganic phosphorus species

 Inorganic phosphorus species
IonsH2PO4HPO42−PO43−
Ca CaH2PO4 CaHPO4 CaPO4, Ca5(PO4)3OH (hydroxyapatite) 
Mg MgH2PO4 MgHPO4 MgPO4, NH4MgPO4·6H2O (struvite) 
Na  NaHPO4  
K  KHPO4  
Fe FeH2PO4 FeHPO4 FePO4, Fe3(PO4)2:8H2O (vivianite) 
 Inorganic phosphorus species
IonsH2PO4HPO42−PO43−
Ca CaH2PO4 CaHPO4 CaPO4, Ca5(PO4)3OH (hydroxyapatite) 
Mg MgH2PO4 MgHPO4 MgPO4, NH4MgPO4·6H2O (struvite) 
Na  NaHPO4  
K  KHPO4  
Fe FeH2PO4 FeHPO4 FePO4, Fe3(PO4)2:8H2O (vivianite) 

The free ions (e.g. Ca2+ or Na+ in Table 2) were not evaluated under this study, as ion-pairing is less important compared to ion activity corrections (Solon et al. 2015). Therefore, free ions were grouped into states such as cations (SCat) and anions (SAn), as considered in ADM1 (Batstone et al. 2002).

The errors incurred by any AD models using ideal solution assumption (such as the ADM1 (Batstone et al. 2002)) on the aforementioned variables (IA/PA, SNH3, PCO2, SH2S and phosphates) were assessed at different values of ionic strength. A range of ionic strength from 0 to 0.5 M was used to evaluate those errors, as it covers the range of typical ionic strength of AD influents (Batstone et al. 2012) and Davies equation is reasonably accurate within that range (Bethke 2007).

To compute the pH range at a constant IC, net strong acid or base monovalent ions were added to the solution while cation/anion pairs were subtracted to compensate for the IC change caused.

Alkalinities were computed by simulated titrations from the initial calculated pH of the wastewater to the reference pH values for TA (4.3) and PA (5.75) (Ripley et al. 1986) by adding a concentration of a strong acid (assumed to be HCl 1 M). The value of the ratio intermediate/partial alkalinity was calculated as:  
formula
(4)

Error calculations of key variables in physicochemistry simplifications in AD modelling

Errors were calculated as the difference between the values under ideal solution assumption (e.g. with activity corrections neglected) with respect to the values accounting for the effect of ionic strength on activities. The assessment was conducted on a typical AD scenario as defined by the steady state concentrations described for the Benchmark Simulation Model No. 2 (BSM2) (Rosen et al. 2006). BSM2 consists of a simulation of a digester treating typical municipal wastewater and operating at stable conditions. The wastewater described in BSM2 has an ionic strength of ∼0.17 M. To consider the errors on H2S speciation, a digester where sulfate-reduction had a relevant effect was evaluated (Barrera et al. 2015). For inorganic phosphorus modelling, a scenario was selected in which phosphorus was modelled in an anaerobic digester (Flores-Alsina et al. 2016).

Free ammonia and hydrogen sulphide are relevant inhibitory species in AD. The impact of their errors is reflected on the inhibitory terms commonly used in modelling acetoclastic methanogenesis and other AD reaction kinetics. Therefore, the errors for those two variables were also studied in terms of its inhibition terms. The inhibition terms for ammonia (INH3) and hydrogen sulphide (IH2S) were calculated as:  
formula
(5)
where KI,NH3 and KI,H2S are the inhibition constants for NH3 and H2S. A value of KI,NH3 of 1.8·10−3 M (Batstone et al. 2002) and KI,H2S of 7.5·10−3 M – an average among the numerous values of KI,H2S reported in literature (Fedorovich et al. 2003; Barrera et al. 2015) – were used to assess the magnitude of the committed errors. The magnitude of the errors in the calculated inhibition depends on the adopted values. Higher values of KI imply smaller error and vice versa.
Those terms (INH3 and IH2S) are used to capture the microorganisms inhibition due to a concentration of the selected components (NH3 or H2S). Other inhibition terms, such as inhibition at low pH, are used in AD models. For this work, the combinatorial errors on pH and an inhibition term are not evaluated as the maximum observed errors occur at neutral pH in most cases. Absolute and relative errors for each variable studied (x) were calculated as:  
formula
(6)
 
formula
(7)

Dynamic AD simulation

To evaluate the impact of neglecting activity correction by ionic strength in dynamic AD simulations, data from an AD reactor subject to a series of increasing organic loading rates were simulated. The data correspond to a 1 m3 UASB-UAF reactor fed with synthetic wastewater containing diluted wine as described by Rodríguez et al. (2006) and by Ruiz et al. (2004). Two comparative identical simulations using the same parameters were run with and without activity correction. In this case, the ionic strength of the solution was on the low range between 0.04 and 0.07 M. The system was simulated using the ADM1 (Batstone et al. 2002) under an Excel-MATLAB framework described by Rodríguez et al. (2009).

RESULTS AND DISCUSSION

The error assessment results for alkalinity, PCO2, ammonia, H2S and inorganic phosphorus are presented below.

Alkalinity

Figure 3 presents simulations of the errors incurred on alkalinities if activity correction is neglected. For this simulation, the BSM2 AD conditions (Rosen et al. 2006) were used.

Figure 3

(a) Effect of activity correction by ionic strength on the model predictions of alkalinity (PA, TA and IA/PA). Total and partial alkalinities are calculated through a simulated titration by adding cations until pH of 5.75 and 4.3, respectively, were reached. (b) Errors are computed with the non-ideal solution as reference (i.e. error indicates the overestimation when ideal solution is assumed).

Figure 3

(a) Effect of activity correction by ionic strength on the model predictions of alkalinity (PA, TA and IA/PA). Total and partial alkalinities are calculated through a simulated titration by adding cations until pH of 5.75 and 4.3, respectively, were reached. (b) Errors are computed with the non-ideal solution as reference (i.e. error indicates the overestimation when ideal solution is assumed).

The difference between ideal and non-ideal solution shown in Figure 3(a) are due to the incorporation of the activity correction (γ), which diverges from the ideal solution (γ = 1) with increasing ionic strengths. The results presented in Figure 3(b) show the significance of the errors on partial alkalinity predictions such that they cannot be generally neglected. This error amplifies in the prediction of IA/PA up to very significant values over 20% and at already typical ionic strength values for AD of municipal wastewater. Considering the importance of the IA/PA ratio for control purposes, this result cannot be ignored. It calls for better physicochemistry modelling in AD systems, with respect to, e.g., current ADM1 approach, or for a revision of the reference values for IA/PA for stable AD operation. The magnitude of the errors for this variable implies that neglecting activity corrections at ionic strengths lower than 0.2 M (Solon et al. 2015) should be further adjusted. The incurred errors for IA/PA ratio become relevant at IC > 0.01 M, highlighting the need to include activity corrections in biochemical model-based controllers.

Carbon dioxide in biogas

The accurate prediction of biogas composition is also a key objective for AD models. The impact of the ideal solution assumption on the predicted PCO2 was also evaluated for a simple example with glucose as sole substrate under typical AD operational conditions and assuming that biogas composition contained only CH4 and CO2 (see Supplementary material, available with the online version of this paper).

The moderately significant errors on the prediction of CO2 composition in the biogas (presented in Figure 4(a)) also relate to the errors in the alkalinity, but depend on other factors such as type of substrate and process HRT and mixing conditions. Figure 4(b) presents a simulation of the errors incurred on the predicted composition of carbon dioxide in the biogas if ideal solution is assumed. Moderately significant errors on the biogas composition ranging from 2–10% are incurred in the prediction of PCO2 already at typical IC and moderately alkaline pH values (7.5–8.5). The errors on CO2 become more relevant due to the multiple deprotonations of inorganic carbon (Figure 1). When modelling species that have multiple deprotonations, the value of the activity coefficient of each ion decreases exponentially with its charge (see Equation (2)). For instance, the activity coefficient of CO32− will be much lower (due to having a −2 charge) than the more protonated species HCO3 (which has a charge of −1). The more charged an ion is, the more the activity coefficient will diverge from the ideal solution (γ = 1). When the pH is more alkaline, more inorganic carbon is in the form of CO32−. Therefore, the error estimating the concentration of CO32− becomes larger, and that is reflected in the committed error in the CO2 biogas composition at moderate alkaline values.

Figure 4

(a) Effect of activity correction by ionic strength on a model prediction of PCO2 for an example case scenario with glucose as sole substrate. (b) Errors are computed with the non-ideal solution as reference (i.e. error indicates the overestimation when ideal solution is assumed).

Figure 4

(a) Effect of activity correction by ionic strength on a model prediction of PCO2 for an example case scenario with glucose as sole substrate. (b) Errors are computed with the non-ideal solution as reference (i.e. error indicates the overestimation when ideal solution is assumed).

Free ammonia and hydrogen sulfide

The error assessments for free ammonia (for BSM2 AD conditions) and hydrogen sulfide (for an AD reactor with sulfate modelled (Barrera et al. 2015)) are presented respectively in Figure 5(a) and 5(b).

Figure 5

Concentrations of (a) NH3 and (b) H2S and the value of its inhibitions (INH3, IH2S; top figures, right axis) for ideal and non-ideal solutions (with activity corrected by ionic strength). Absolute errors for species concentrations and absolute error for their inhibition values are also shown, taking the non-ideal solution as a reference (a negative value means that the ideal solution calculates a larger inhibition than the activity-corrected one).

Figure 5

Concentrations of (a) NH3 and (b) H2S and the value of its inhibitions (INH3, IH2S; top figures, right axis) for ideal and non-ideal solutions (with activity corrected by ionic strength). Absolute errors for species concentrations and absolute error for their inhibition values are also shown, taking the non-ideal solution as a reference (a negative value means that the ideal solution calculates a larger inhibition than the activity-corrected one).

The potential errors reported in the prediction of free ammonia if the ideal solution is assumed (Figure 5(a)) can become of moderate significance in pH ranges of interest in AD operation, especially in regards to the impact on inhibition terms. Ideal solution simplifications overpredict the inhibition by free ammonia to relevant extents already at typical AD ionic strength values (∼5% for IC = 0.17 M, which is the ionic strength at BSM2 conditions). These model errors can lead to inadequate model overfitting of the ammonia inhibition parameters to fit experimental data (García-Gen et al. 2013). Less important (∼2% under BSM2 conditions) are the errors committed in the prediction of H2S inhibition in models dealing with sulfate reduction in AD (Figure 5(b)). However, these errors could become more relevant if the adopted KI,H2S is smaller than the one used in this study (KI,H2S = 7.5·10−3 M). Despite H2S being a species subject to multiple deprotonations (similar to the inorganic carbon), the impact of activity correction is not that important because the second deprotonation (with a double negative charge) occurs only at very high pH (pKa ∼ 17). For practical purposes, H2S speciation can be considered as a single deprotonation. Such high pH values are not observed in AD environments.

Inorganic phosphorus

An evaluation of errors for the three main deprotonations is shown in Figure 6. The magnitude of the errors clearly indicate that the assumption of ideal solution is not feasible if inorganic phosphorus and its possible precipitations are to be described in the model.

Figure 6

Inorganic phosphorus species concentrations (SH2PO4, SHPO4 and SPO4) for ideal and non-ideal solutions (with activity corrected by ionic strength) (top), along with its absolute and relative errors. Positive relative errors indicate an overestimation of the species concentration by the ideal solution while negative values indicate an underestimation of the species concentration when using ideal solution.

Figure 6

Inorganic phosphorus species concentrations (SH2PO4, SHPO4 and SPO4) for ideal and non-ideal solutions (with activity corrected by ionic strength) (top), along with its absolute and relative errors. Positive relative errors indicate an overestimation of the species concentration by the ideal solution while negative values indicate an underestimation of the species concentration when using ideal solution.

Inorganic phosphorus, a component of interest in AD subject to multiple deprotonations function of pH, is particularly sensitive to the assumption of ideal solution due to its multivalency (see Figure 6). In the case of inorganic phosphorus, the errors are bigger compared to inorganic carbon or H2S because the predominant species at neutral pH is HPO42− (with a pKa ∼ 7.2). The activity coefficient of HPO42− is much lower (due to having a double negative charge) than the more protonated species H2PO4 (which has a single negative charge).

The most predominant species of inorganic phosphorus in typical AD conditions is H2PO4 (Figure 6(a)). The absolute errors under ideal solution assumption are high within typical AD pH ranges, and the magnitude of the relative errors is of high relevance (∼40% of overestimation at IC = 0.17 M, which corresponds to that of the BSM2 under steady state conditions), increasing at higher ionic strengths. For biochemical models in which phosphorus speciation is of importance, an accurate physicochemistry prediction cannot be expected without activity correction. The other major species of phosphorus show similar magnitude of errors. HPO42− (Figure 6(b)) is underestimated by ∼15% at BSM2 conditions, meaning that decreased precipitation (e.g. calcium phosphate, iron phosphate) would be predicted if the ideal solution is assumed. Although PO43− concentration (Figure 6(c)) in AD conditions is very low compared to the other phosphorus species, the relative errors of its underestimation reach ∼60% for phosphate under BSM2 conditions if the ideal solution is assumed. This error is even bigger compared to HPO43− because its charge (triple negative) impacts the activity coefficient to a larger extent. That error is magnified at alkaline pH, as it becomes closer to its pKa (∼12.3). These errors would translate into much less predicted precipitation (e.g. struvite or hydroxyapatite) when these species are involved.

Overall, the results obtained indicate that inorganic phosphorus modelling requires activity corrections. This aligns to previous recommendations from literature (Batstone et al. 2012). Furthermore, phosphorus species pair with other protonated species (such as Na+ or Ca2+). Therefore, an approach with ion-pairing implementation such as that described by Flores-Alsina et al. (2015) or Huber et al. (2017) should be used if phosphorus speciation is to be modelled.

Dynamic AD process simulation and implications for control

Comparative simulations using the ADM1 model for a dynamic AD case scenario (Ruiz et al. 2004) involving a series of increasing overloads were conducted.

The dynamic simulation of a series of overloads of an AD ethanol-rich water treatment illustrates some of the control implications of the ideal solution assumptions. Figure 7(a) shows how the model not considering activity corrections predicts a higher pH than the one that corrects for ion activity by ionic strength. This implies an underprediction of the reactor destabilisation risks by AD models not considering the ionic strength impact on ion activity.

Figure 7

Comparative dynamic ADM1 simulations adjusted to the experimental data from an AD system treating ethanol containing wastewater subject to a series of increasing overloads until system destabilisation (Ruiz et al. 2004). Simulated data considering ideal solution (dashed lines) and non-ideal solution (i.e. including activity correction) (solid lines). The presented variables are (a) pH, IA/PA ratio, and (b) concentrations of butyrate (Sbu), propionate (Spro) and acetate (Sac).

Figure 7

Comparative dynamic ADM1 simulations adjusted to the experimental data from an AD system treating ethanol containing wastewater subject to a series of increasing overloads until system destabilisation (Ruiz et al. 2004). Simulated data considering ideal solution (dashed lines) and non-ideal solution (i.e. including activity correction) (solid lines). The presented variables are (a) pH, IA/PA ratio, and (b) concentrations of butyrate (Sbu), propionate (Spro) and acetate (Sac).

Particularly important in control operation is the error committed on the α (IA/PA ratio) as it is a parameter largely used for AD stable control. Neglecting activity correction causes an underestimation of the acidification risk in the reactor, which carries problems if the models is used to tune AD process controllers for optimum process performance and can lead to less conservative controllers and increased destabilisation risks.

Alternatives for physicochemistry calculation in AD modelling

Existing and proposed recommended alternatives for the calculation of the ionic speciation in AD models are summarised in Figure 8. Current ideal solution approaches require a single loop iteration only if pH is calculated as an output (idS-(St) in Figure 8). However, pH is a commonly measured variable in AD that can also (and indeed should) be provided as an input when there is uncertainty in the water composition. In this case, no iterative loop calculation is required (idS(St+pH) case).

Figure 8

Alternatives for physicochemistry calculation in AD models: current ideal solution approaches (dashed frame top) and recommended non-ideal solution approaches (solid frame bottom); with pH calculated as an output (left) and alternatives with pH as an input (right). Ic(prev.) refers to an approximation to the ionic strength of the solution considering the value of an earlier time step.

Figure 8

Alternatives for physicochemistry calculation in AD models: current ideal solution approaches (dashed frame top) and recommended non-ideal solution approaches (solid frame bottom); with pH calculated as an output (left) and alternatives with pH as an input (right). Ic(prev.) refers to an approximation to the ionic strength of the solution considering the value of an earlier time step.

The adoption of non-ideal solution adds the ionic strength (IC) as a variable in the procedure. Four possible alternatives (presented in Figure 8) arise depending on whether pH and IC are provided as inputs or calculated (see Supplementary material for the alternative implementations, available with the online version of this paper).

Alternatives are possible if the implementation of loops for the IC wants to be avoided, namely the NidS-(St+Ic) and the NidS-(St+pH+Ic) with pH as input (see Figure 8). In these cases, the IC could be provided as input based on the value at the previous time step and activity corrections applied accordingly. This would ensure that most of the errors are minimised by assuming that IC does not change fast over time and would require no additional iterations for the ionic strength.

The most accurate approach with calculation of both pH and IC providing only the state variables (as total concentrations) as input is the NidS-St case in Figure 8 that has been used for the error analyses presented above. Although the NidS-St is the only one involving a double loop iteration, the additional computational costs observed were minimal due to the very fast convergence of the IC loop. It is therefore the recommended approach together with its variation with measured pH as input (the NidS-(St+pH) case with single-loop iteration for the IC alone).

The implementations described for activity corrections should be compatible with most widely used softwares by AD model users, such as Aquasim, SIMBA or WEST, by following the flow diagram provided in Figure 2 (and in the Supplementary material).

CONCLUSIONS

The comprehensive computation of the errors incurred by those physicochemistry models in AD that assume ideal solution demonstrates the need for the systematic incorporation of activity corrections in AD models, particularly the following:

  • The errors on the alkalinity ratio are very significant even at low ionic strengths (IC < 0.01 M) and can compromise model-based control actions. If the IA/PA ratio is used as a model parameter for control purposes, activity corrections must be included.

  • The inhibition factors of key components such as NH3 and H2S are also moderately impacted (∼ 4% at IC < 0.2 M) if activity corrections are not considered. It is not necessary to include activity corrections if those are the only variables of interest of the model or if the ionic strength of the wastewater is lower than 0.2 M.

  • The magnitude of the errors committed in modelling inorganic phosphorus speciation are unacceptably high and will impede precipitation modelling. Therefore, if phosphorus is modelled, activity corrections have to be included.

  • The Davies method appears as an adequate accurate option as it does not require additional parameters and involves very low computational cost.

ACKNOWLEDGEMENTS

This work has been supported by Masdar Institute of Science and Technology (grant number SSG2015-0057) and the Government of Abu Dhabi.

REFERENCES

REFERENCES
Allison
,
J. D.
,
Brown
,
D. S.
&
NovoGradac
,
K. J.
1991
MINTEQA2/PRODEFA2, A Geochemical Assessment Model for Environmental Systems: Version 3.0 User's Manual
.
Environmental Research Laboratory, Office of Research and Development, US Environmental Protection Agency
,
Athens, GA, USA
.
Angelidaki
,
I.
,
Ellegaard
,
L.
&
Ahring
,
B. K.
1999
A comprehensive model of anaerobic bioconversion of complex substrates to biogas
.
Biotechnology and Bioengineering
63
(
3
),
364
372
.
Batstone
,
D. J.
2009
Towards a generalised physicochemical modelling framework
.
Reviews in Environmental Science and Bio/Technology
8
(
2
),
113
114
.
Batstone
,
D. J.
&
Steyer
,
J.-P.
2007
Use of modelling to evaluate best control practice for winery-type wastewaters
.
Water Science & Technology
56
(
2
),
147
152
.
Batstone
,
D. J.
,
Keller
,
J.
,
Angelidaki
,
I.
,
Kalyuzhnyi
,
S. V.
,
Pavlostathis
,
S. G.
,
Rozzi
,
A.
,
Sanders
,
W. T. M.
,
Siegrist
,
H.
&
Vavilin
,
V. A.
2002
The IWA Anaerobic Digestion Model No. 1 (ADM1)
.
Water Science & Technology
45
(
10
),
65
73
.
Batstone
,
D. J.
,
Amerlinck
,
Y.
,
Ekama
,
G.
,
Goel
,
R.
,
Grau
,
P.
,
Johnson
,
B.
,
Kaya
,
I.
,
Steyer
,
J.-P.
,
Tait
,
S.
,
Takács
,
I.
,
Vanrolleghem
,
P. A.
,
Brouckaert
,
C. J.
&
Volcke
,
E.
2012
Towards a generalized physicochemical framework
.
Water Science & Technology
66
(
6
),
1147
1161
.
Bethke
,
C. M.
2007
Geochemical and Biogeochemical Reaction Modeling
.
Cambridge University Press
,
New York, NY, USA
.
Davies
,
C. W.
1962
Ion Association
.
Butterworths, Washington, DC, USA
.
Fedorovich
,
V.
,
Lens
,
P.
&
Kalyuzhnyi
,
S.
2003
Extension of Anaerobic Digestion Model No. 1 with processes of sulfate reduction
.
Applied Biochemistry and Biotechnology
109
(
033
),
33
45
.
FloresAlsina
,
X.
,
Kazadi Mbamba
,
C.
,
Solon
,
K.
,
Vrecko
,
D.
,
Tait
,
S.
,
Batstone
,
D. J.
,
Jeppsson
,
U.
&
Gernaey
,
K. V.
2015
A plant-wide aqueous phase chemistry module describing pH variations and ion speciation/pairing in wastewater treatment process models
.
Water Research
85
,
255
265
.
FloresAlsina
,
X.
,
Solon
,
K.
,
Kazadi Mbamba
,
C.
,
Tait
,
S.
,
Gernaey
,
K. V.
,
Jeppsson
,
U.
&
Batstone
,
D. J.
2016
Modelling phosphorus (P), sulfur (S) and iron (Fe) interactions for dynamic simulations of anaerobic digestion processes
.
Water Research
95
,
370
382
.
García-Gen
,
S.
,
Lema
,
J. M.
&
Rodríguez
,
J.
2013
Generalised modelling approach for anaerobic co-digestion of fermentable substrates
.
Bioresource Technology
147
,
525
533
.
González-Cabaleiro
,
R.
2015
Bioenergetics-based Modelling of Microbial Ecosystems for Biotechnological Applications
.
PhD, Grupo de Ingeniería Ambiental y Bioprocesos
,
Universidade de Santiago de Compostela, Santiago de Compostela
,
Spain
.
Gustafsson
,
J.
2014
Visual MINTEQ ver. 3.1. https://vminteq.lwr.kth.se/
.
Helgeson
,
H. C.
1969
Thermodynamics of hydrothermal systems at elevated temperatures and pressures
.
American Journal of Science
267
(
7
),
729
804
.
Huber
,
P.
,
Neyret
,
C.
&
Fourest
,
E.
2017
Implementation of the anaerobic digestion model (ADM1) in the PHREEQC chemistry engine
.
Water Science & Technology
76
(
5–6
),
1090
1103
.
Jenkins
,
S. R.
,
Morgan
,
J. M.
&
Sawyer
,
C. L.
1983
Measuring anaerobic sludge digestion and growth by a simple alkalimetric titration
.
Journal (Water Pollution Control Federation)
55
(
5
),
448
453
.
Kazadi Mbamba
,
C.
,
Tait
,
S.
,
Flores-Alsina
,
X.
&
Batstone
,
D. J.
2015
A systematic study of multiple minerals precipitation modelling in wastewater treatment
.
Water Research
85
,
359
370
.
Kazadi Mbamba
,
C.
,
Flores-Alsina
,
X.
,
Batstone
,
D. J.
&
Tait
,
S.
2016
Validation of a plant-wide phosphorus modelling approach with minerals precipitation in a full-scale WWTP
.
Water Research
100
,
169
183
.
Malmberg
,
C.
&
Maryott
,
A.
1956
Dielectric constant of water from 0 to 100 C
.
Journal of Research of the National Bureau of Standards
56
(
1
),
1
8
.
Nielsen
,
A. M.
,
Spanjers
,
H.
&
Volcke
,
E. I.
2008
Calculating pH in pig manure taking into account ionic strength
.
Water Science & Technology
57
(
11
),
1785
1790
.
Parkhurst
,
D. L.
&
Appelo
,
C. A. J.
2013
Description of Input and Examples for PHREEQC Version 3 – A Computer Program for Speciation, Batch- Reaction, One-Dimensional Transport, and Inverse Geochemical Calculations
.
US Geological Survey
, p.
497
.
Pitzer
,
K. S.
1987
A thermodynamic model for aqueous solutions of liquid-like density
.
Reviews in Mineralogy and Geochemistry
17
(
1
),
97
.
Ripley
,
L.
,
Boyle
,
W.
&
Converse
,
J.
1986
Improved alkalimetric monitoring for anaerobic digestion of high-strength wastes
.
Journal (Water Pollution Control Federation)
58
,
406
411
.
Rodríguez
,
J.
,
Kleerebezem
,
R.
,
Lema
,
J. M.
&
van Loosdrecht
,
M. C.
2006
Modeling product formation in anaerobic mixed culture fermentations
.
Biotechnology and Bioengineering
93
(
3
),
592
606
.
Rodríguez
,
J.
,
Premier
,
G. C.
,
Dinsdale
,
R.
&
Guwy
,
A. J.
2009
An implementation framework for wastewater treatment models requiring a minimum programming expertise
.
Water Science & Technology
59
(
2
),
367
380
.
Rosen
,
C.
,
Vrecko
,
D.
,
Gernaey
,
K. V.
,
Pons
,
M. N.
&
Jeppsson
,
U.
2006
Implementing ADM1 for plant-wide benchmark simulations in Matlab/Simulink
.
Water Science & Technology
54
(
4
),
11
19
.
Ruiz
,
G.
,
Rodrıguez
,
J.
,
Roca
,
E.
&
Lema
,
J.
2004
Modification of the IWA-ADM1 for application to anaerobic treatment of ethanolic wastewater from wine factories
. In:
Proceedings of 10th IWA World Congress on Anaerobic Digestion (AD10)
,
Montreal, Canada
, pp.
1341
1344
.
Schoen
,
M. A.
,
Sperl
,
D.
,
Gadermaier
,
M.
,
Goberna
,
M.
,
Franke-Whittle
,
I.
,
Insam
,
H.
,
Ablinger
,
J.
&
Wett
,
B.
2009
Population dynamics at digester overload conditions
.
Bioresource Technology
100
(
23
),
5648
5655
.
Siegrist
,
H.
,
Vogt
,
D.
,
Garcia-Heras
,
J. L.
&
Gujer
,
W.
2002
Mathematical model for meso- and thermophilic anaerobic sewage sludge digestion
.
Environmental Science & Technology
36
(
5
),
1113
1123
.
Smith
,
S. A.
&
Chen
,
S.
2006
Activity corrections for ionization constants in defined media
.
Water Science & Technology
54
(
4
),
21
29
.
Solon
,
K.
,
Flores-Alsina
,
X.
,
Mbamba
,
C. K.
,
Volcke
,
E. I. P.
,
Tait
,
S.
,
Batstone
,
D.
,
Gernaey
,
K. V.
&
Jeppsson
,
U.
2015
Effects of ionic strength and ion pairing on (plant-wide) modelling of anaerobic digestion
.
Water Research
70
,
235
245
.
Sötemann
,
S.
,
Van Rensburg
,
P.
,
Ristow
,
N.
,
Wentzel
,
M.
,
Loewenthal
,
R.
&
Ekama
,
G.
2005
Integrated chemical/physical and biological processes modeling Part 2-Anaerobic digestion of sewage sludges
.
Water SA
31
(
4
),
545
568
.
Stumm
,
W.
&
Morgan
,
J. J.
1996
Aquatic Chemistry: Chemical Equilibria and Rates in Natural Waters
.
Wiley
,
New York, NY, USA
.
Vaneeckhaute
,
C.
,
Claeys
,
F. H. A.
,
Tack
,
F. M. G.
,
Meers
,
E.
,
Belia
,
E.
&
Vanrolleghem
,
P. A.
2018
Development, implementation, and validation of a generic nutrient recovery model (NRM) library
.
Environmental Modelling & Software
99
,
170
209
.
Yenigün
,
O.
&
Demirel
,
B.
2013
Ammonia inhibition in anaerobic digestion: a review
.
Process Biochemistry
48
(
5–6
),
901
911
.

Supplementary data