Modeling Anaerobic Co-digestion of Food Wastes and Cattle Manure in an Industrial Plant: A System Dynamic Approach

Abstract


Introduction
Recognition of the negative impact of fossil fuels and biodegradable wastes on our environment has prompted research on the production of biogas through anaerobic digestion (AD).AD of biodegradable wastes is a biological process that produces methane, carbon dioxide, and other gases in trace amounts.Co-digestion takes place when 2 or more substrates are digested simultaneously.Co-digestion is known to be suitable for C/N ratio enhancement, methane production kinetics improvement, and inhibitory effect mitigation [1,2].Due to the potential of biogas production for waste management and source of alternative energy, the number of biogas plants is increasing [3].AD processes have been widely considered for biogas production from biodegradable wastes and the major biochemical stages recognized in AD processes are hydrolysis, acidogenesis, acetogenesis, and methanogenesis.These stages involve a different group of micro-organisms [4].
Different substrates can be co-digested to increase the quantity and quality of biogas production.However, a considerable reduction in biogas production due to process instability can occur as a result of inapt selection of co-substrates and co-substrate composition [5].Instability of the digester system as a result of an increase in volatile fatty acid (VFA) subsequently followed by a reduction in pH and biogas production has been reported by Donoso-Bravo et al. [6].To understand the system behaviour thereby applying process optimization and control, a large number of biodigester models [7][8][9] have been conducted.The Hill model with four state variables is visualized as a too simple model for representing the complex behaviour of a full-scale digester and therefore discarded in this study.The ADM1 was developed by International Water Association (IWA) Task Group to describe 3 physicochemical, 19 biochemical processes, and 7 different bacterial populations within the digester system.Four stages of anaerobic digestion along with the death and growth of separate biomass fractions are represented in the model.The model has 35 state variables and more than 85 model parameters.
The ADM1 is mostly accepted and extensively used among researchers [10] and it has manifested itself to be reliable in its demonstration in an agricultural biogas plant [11].The diagrammatic representation of these processes is shown in Fig. 1 [12].It was demonstrated by Biernacki et al. [13] that ADM1 had the capacity of describing an industrial biogas power plant in terms of biogas output and methane concentration.However, the structure of the ADM1 involves a large number of kinetic and stoichiometric parameters identification and full substrate characterization and therefore difficult to be utilized for optimization and control design [14].This has motivated many researchers to search for a simple but accurate model involving a few numbers of processes and model parameters.The anaerobic digestion model number 2 (AM2) which is a modification from Andrews and Graef [15] and developed for the design and control of a sludge digester system is being considered in this study.The model which contains 6 state variables and 13 model parameters was designed to be suitable for methane, biomass, VFA, alkalinity (Z), and organic matter prediction.Numerous studies have shown the applicability of the AM2 for digesters fed with different substrates.
Bernard et al. [16] developed an AM2 for an anaerobic wastewater treatment plant.Parameters identification during steady states and model validation were performed on the system.The model identification carried out under steady states proves to be efficient in dynamic conditions.Ficara et al. [17] carried out a comparative study of the anaerobic digestion model: ADM1 and AM2 to describe the anaerobic degradation of wasteactivated sludge.They concluded that AM2 could sufficiently replace ADM1 by including the hydrolysis stage in the model.The modeling of anaerobic digestion of maize silage using both AM2 and ADM1 was also conducted by Arzate et al. [4].A comparative study of the simulation results was analyzed.The analysis showed that consideration should be given to the implementation of the AM2 in an adaptive framework of a biogas plant.Delgadilo-Mirquez et al. [7] applied AM2 for modeling anaerobic digestion processes from organic wastes.Sensitivity analysis of parameters and model validation was conducted.Attar and Haugen [18] implemented AM2 on an anaerobic digestion reactor of water resources recovery facility in Norway.Model identification and validation were carried out.They all concluded that AM2 was sufficient to model accurately the behaviour of digester systems.The combination of food wastes and cattle manure has been recognized as a good substrate for anaerobic co-digestion because of the ease of biodegradability of food wastes and the high buffering capacity provided by cattle manure.However, the application of AM2 to simulate co-digestion processes of food wastes and cattle manure in an industrial plant has been rarely reported in the literature.Of the published AD models in the literature, process simulation was accomplished by SIMBA under MATLAB/SIMULINK [19], AQUASIM v 2.1b® [20,21], and Aspen plus [22,23].To overcome the challenge of interpreting and understanding the dynamic model of a complex biological system such as the AD process due to a lack of good background in mathematics, to make a good decision, SD modeling of biogas production system was considered.
SD is a methodology and mathematical modeling technique to manage, understand and discuss complex feedback systems.Several practical applications of the system dynamics model on different issues have been established in the literature such as solid waste management [24][25][26], wastewater treatment [27,28], assessment of lake eutrophication [29], wetland study [30,31], and environmental impact assessment of coalfields [32].
Modeling and simulation of the co-digestion process in an industrial plant using a SD approach are very scarce in the literature.Therefore, modeling and simulation of co-digestion processes in a SD modeling software using AM2 would not only improve the understanding of the behaviour and structure of the model but also serve as a valuable tool with less computation power and small model parameters for monitory and prediction of biogas production in an industrial plant.This research work is aimed at developing a SD model of the co-digestion process of food wastes and cattle manure in an industrial biogas plant's digester based on modified AM2 to display its capability of representing the dynamic behaviour of the system.The process stability of the co-digestion process and sensitivity analysis of biogas flow and concentrations to some important process variables were also investigated.The proposed model for the system under consideration would give insight into system dynamics thereby assisting in optimizing the rate of biogas output.

System Description
For the plausibility of the proposed model to be determined by direct comparison with the existing models in the literature, the operating system and the conditions must be the same.The systems considered for the simulation had been previously modeled by two different versions of ADM1.The system was the EWE biogas plant which was built in 1996 in Wittmud Lower Saxony, Germany.The system had been previously modeled by Biernacki et al. [13] using modified ADM1without including lactate and Satpathy et al. [33] using lactate including modified ADM1.The data from the plant is worthy of consideration because of its purpose, relevance, correctness, and source.The biogas plant was designed to consist of two parallel farmenters, each of 3500 m 3 , with an average retention time of 20 d, an underground tank of 1900 m 3 where substrates are collected, and a mixing tank (620 m 3 ), where two or more substrates are mixed to obtain a homogeneous consistency and an optimal buffer/substrate ratio [33].A detailed description of the plant can be found in the work of Biernacki et al. [13].During the data collection of the plant, it was fed with cattle manure of 180 m 3 d -1 and food wastes of 100 m 3 d -1 on average.The period of data collection was 28 d.The measurement of biogas produced was measured from the recorded cumulative gas flow and infrared sensors were used to monitor biogas composition.
The biogas produced is converted to heat (3.4 MW) and electricity (2.5 MW) in a combined heat and power (CHP) unit.The total biogas produced from the existing biogas power plant after 28 d was 127,711 m 3 with an average methane volume fraction of 66.85%.

Substrates Characterization
Detailed mathematical characterization is required because of the strong impact of substrate characterization on the composition of biogas.The substrates were fractionated into an easily biodegradable fraction (Ss), a slowly biodegradable fraction (Xs), and an inert fraction (Xi).Xs which comprises carbohydrates (Xch), proteins (Xpr), lipids (Xli) fraction and Xi were computed from ADM1 fractions provided by Biernacki et al. [13].These ADM1 fractions were calculated from the characteristics of the substrates obtained through Weender analysis and Van Soest estension.Ss such as monosaccharide (Ssu), amino acids (Saa), long-chain fatty acids (Sfa), total acetate (Sac), inorganic carbon (Sic), total valerate (Sva), total butyrate (Sbu), total propionate (Spr), and hydrogen carbonate salt ( "#$ % ) was determined by multiplying the total COD of the considered substrates with percentages calculated from the work of Fisgativa et al. [34].Table 1 gives an overview of the calculated basic data needed for estimating input data.Using the methodology proposed by Ficara et al. [17], the data in Table 1 was transferred to AM2 environment using Eqs.(1)(2)(3)(4).The influent compositions in terms of organic substrate, S1 (kgCOD m -3 ), VFA, S2 (mol m -3 ), alkalinity, Z (mol m -3 ) and inorganic carbon, C (mol m -3 ) were calculated for cattle manure and food wastes 2 = 1000 5

The ADM1
ADM1 is a model including 19 differential and 12 algebraic equations.It involves 4 stages of anaerobic biodegradation: hydrolysis, acidogenesis, acetogenesis, and methanogenesis.It has been proved by researchers that ADM1 is the most appropriate model for describing biogas production in a biogas power plant.[13,33].
The applicability of the modified ADM1 to the EWE biogas power plant was investigated by Biernacki et al. [13].The value of biogas and methane production from the simulation results deviated from the measured value by 1.3% and 1.7% respectively.Satpathy et al. [33] realized the limitation of ADM1 applicability to a digester Biogas overloaded with carbohydrates and included lactate into the ADM1.The lactate-included ADM1 proved better in reliability when applied to the EWE biogas power plant.

Digestate
Where D is the dilution rate expressed as the ratio of the flow rate of the substrate and effective volume of the digester,  is the fraction of bacteria in the liquid phase, K1 is the yield for substrate degradation, K2 is the yield for VFA production, K3 is the yield for VFA consumption, K4 is the yield for CO2 production,  > and  9 stand respectively for the specific growth rate of X1 and X2.
Monod kinetics given in Eq. ( 11) and Haldane kinetics given in Eq. ( 12) were respectively utilized for modeling X1 and X2 kinetics.
Where  >Z+[ and  9Z+[ are the maximum growth rate for  > and  9 respectively;  (> is the half-saturation constant associated with S1, ,  (9 and  c are the half-saturation constant and inhibition constant associated with S2 respectively. The performance of AD processes is affected by the pH and temperature of the digester content [35].In this model, the expression considered for modeling pH is given in Eq. ( 13) The alkalinity ratio (AR) which measures the reactor stability [36] is modeled using Eq.(14).
Due to the low solubility of CH4, it is assumed that CH4 produced is in the gas phase, the methane flow rate is estimated using Eq. ( 15) Where, K6 is the yield for CH4 production.
The flow of inorganic carbon from the liquid phase to the gas phase is calculated using Eq. ( 16) following Henry's law.
Where  p  is the liquid-gas transfer coefficient,  e is Henry's constant, and  W is the CO2 partial pressure.
The temperature dependency of  p  given in Eq. ( 16) [37] was considered for the estimation of  p  in this study.
p  = 0.56 + 27.9 Where T o C is the temperature of the reactor assumed to be constant.

The modified AM2
The methane gas produced is assumed to be strongly dependent on the conditions of the methanogenic reaction.
Therefore, the power of the inhibiting substrate of the third term in Eq. ( 12) was changed to 4 from 2 as shown in Eq. ( 20).This was done to further address the complexity involved in the process.To account for the death of acidogens and methanogens during the co-digestion process, the death rates with the value 5.3% and 5.3% in the AM2 [4] were incorporated into the biomass balance model and Eqs.(5 and 6) become Eqs.(18 and 19) respectively.

System Dynamics Model
A SD modeling has been identified as a methodology for studying and managing complex feedback systems.It involves the construction of a stock and flow diagram (SFD) or casual loop diagram to form the model [38].As one of the most advanced graphical system programming dynamic software packages, STELLA was used to develop the mathematical model for an industrial biogas production to address the mechanistic processes.The graphical model of the system was created using the four fundamental tools: stock, flow, converter, and connector.This was closely followed by automatic interpretation of the model in form of difference equations in the equation layer of the package.

Parameter Estimation and Sensitivity Analysis
The model developed needs to fit experimental data to estimate the kinetic parameters and yield coefficients.The initial values of the model parameters were based on the calibrated data of the VEAS water resources recovering facility model in salemmested [16].The parameters were estimated by minimizing the root mean square error (RMSE) given in Eq. ( 21) between the calculated model results and experimental results obtained from Biernack et al. [13].This was carried out in Premium Solver incorporated in Microsoft excel Where  ˆ‰Š 2 is the ith experimental value  (2Z 2 is the ith simulated value and n is the number of data points. The AM2 reconstructed in STELLA software with optimized parameters was then simulated for biogas production.The sensitivity analysis was performed on biogas flow and methane content by considering variations in biodegradable fraction (Xs and Ss), initial concentration of biomass, and substrate fed into the reactor.The stability of the co-digestion process in the digester was also investigated using the validated model results

Results and Discussion
The AM2 was modified by changing the kinetic parameters, yield coefficients, and the variables related to microbial communities.The model was then reconstructed in system dynamics modeling software.Following the model development, the model validation was executed.The simulations of the model in both excel 2007 and SD modeling software generated data that was compared with the data obtained from an industrial biogas plant's digester.
The illustration of the assumed biochemical processes explaining the behaviour of the digester in terms of X1, X2, S1, S2, C, and Z is depicted on the SFD (Fig. 2).The model considered six state variables and 13 parameters.The flows were all expressed on the unit mol m -3 d -1 except flows for estimation of (X1 and S1) and biogas flow which were measured in kgCOD m -3 and m 3 d -1 respectively.The unit of all the entities represented by the stock is mol m -3 .The rectangular box (stock) represents accumulations or state variables, the rates otherwise known as flows are represented by arrows while the circular symbol stands for factors affecting the inflow and outflow of the system.For simplicity, the interrelated reactions leading to biogas production in the SFD were divided into eight different sections.Fig. 2a presented the estimation of X1 while that of X2 is depicted in Fig. 2b.They are both affected by birth rate, death rate, and dilution rate represented in the SFD as flows.Variables affecting each of these flows were attached to the flows.The X1 and X2 generated in Fig. 2a and Fig. 2b were used in Fig. 2c and

Parameters Estimation
To improve the prediction capacity of the model, the measured biogas volume from the biogas plant's digester was used for model calibration.Of all the model parameters, twelve parameters were optimized while others were obtained from the literature.The value of model parameters before and after calibration is presented in Table 3.With the values obtained after calibration, a considerable fit between the simulation results and the experimental data is attained (Fig. 3) Table 3 The values of model parameters before and after calibration

The prediction of the modified AM2
The modified AM2 was run to predict the volumetric flow rate and methane concentrations of the biogas leaving the digester.Fig. 4a shows the comparison between the predicted values and experimental values.
Before the 15th day of the experiment, the predicted biogas output rate was approximately constant while the rest of the day witnessed a significant variation with a trend similar to the experimental data.As expected, the measured daily output rate varied tremendously throughout the period of the experiment since the data were taken from the full-scale industrial digester, probably experiencing variation in substrate concentration.This could subsequently cause instability in the system.The variation of measured data was more noticeable than the model results because the model was based on average substrate compositions and many assumptions.
Nonetheless, there was a quite good correlation between the simulated results and measured data.The model overestimated the total biogas production after 28 d by 775 m 3 (0.6%).As regards methane concentration, the modified AM2 was able to mimic satisfactorily the dynamic results of biogas methane content leaving the biogas plant's digester as indicated in Fig. 4b.Both the measured and the predicted values fluctuated within the narrow range though the predicted concentration (60.01%-73.84%) was slightly higher than the measured values (61.15% -69.00%).The average methane volume fraction was 67.80% which was 1.42% deviated from the experimental value.A relative error of 8% of the measuring systems in the plant was considered.

\ (a) (b)
Fig. 4 Comparison of the simulation of modified AM2 and the measured values from an industrial biogas plant of capacity 7000 m 3 with a varying input composition for a period of 28 days.A relative error of 8% in the measuring system was considered.(a) biogas production and (b) methane concentration.

Model Prediction of the original AM2
For comparison, the original AM2 without any modification was simulated using the calibrated parameters data obtained in the literature.As shown in Fig. 5a, the simulated results were not in good agreement with the experimental data considering biogas flow from the plant.The cumulative biogas volume after 28days was underestimated by 19533.8m3(15.30%).Also, an obvious disparity between the simulated and experimental results was observed in methane concentration as shown in Fig. 5b.Throughout the period considered for the simulation, there was an underestimation of methane content.At the average level, a difference of -15.82% (23.66%) was obtained.This deviation might be a result of the exclusion of representation

Model comparison
An inherent objective of this work was to investigate the applicability of AM2 which is reduced and tractable to a biogas plant's digester previously modeled using ADM1.Biernacki et al. [13] and Satpathy et al. [33] had already simulated the EWE biogas plant using optimized values for the model parameters.Biernacki et al. successively simulated the digester using ADM1 without including lactate.The total biogas flow was overestimated by 1.5% and the methane content was underestimated by 1.74%.Satpathy et al. [33] improved the prediction capacity of the ADM1 by including lactate in the model thereby reducing the error between the measured values and the simulated results to approximately 0.8% and 0.3% for the total biogas production and methane concentration, respectively.The erratic behavious of biogas production and methane concentration from the digester were adequately predicted by the two modified ADM1 as shown in Fig. 3a and 3b [33].The comparison of the different simulation results and the quantitative overall difference between the original and modified models are presented in Table 4.It was shown in Table 4 that the prediction capability of the modified AM2 conformed to that of modified ADM1.Despite the accurate mimicry of biogas quantity and quality of EWE biogas plant using ADMI, the complexity involved in substrate characterization and the high computational effort required especially during model calibration and optimization are of great concern.These setbacks were not only overtaken by AM2 but also predict the total biogas production and methane content with some level of accuracy.It is expected that the process behaviour of some variables in the digester and prediction capacity of this model for biogas quantity and quality can be further improved by considering the varying nature of input variables caused by the varying composition of food wastes and full substrate characterization.To further prove the plausibility of the model, the maximum absolute percentage error (MAPE) of the modified AM2 and the modified ADM1 developed and simulated by Satpathy et al. [33] was compared as shown in Table 5.It was observed that there was an insignificant difference between the MAPE values of both models which indicates that AM2 could be used in place of ADM1 for control and design purposes.The erratic nature of the daily biogas production (3280 m 3 d -1 -6150 m 3 d -1 ) from the biogas plant's digester as indicated in Fig. 4a could be a sign of process instability in the system [39].This might be the result of fluctuations in the substrate properties ( S1, S2, Z, and C) due to continuous changes in food waste composition.
The methane concentrations (61.15% -69.00%) were within the acceptable range and that is an indication that the inhibitory effect did not reflect in the methane content leaving the digester.Using the simulation results, process parameters for monitoring the stability of the processes in the digester were explored.Fig. 6 shows the process parameters (pH, Z, S2, and AR) analysis in the digester.Within the speculated time, the decrease in Z from 145.53 mol m -3 to 110 mol m -3 was accompanied by an increase in S2 from 27.65 mol m -3 to 58.36 mol m -3 .The nitrogen content of food wastes and cattle manure might contribute to the high values of Z.
There was a continuous accumulation of S2 which might cause an imbalance in the system.The critical concentration of S2 seems impracticable to specify since it depends on the substrate concentration and it varies from one AD plant to another [40].Considering the results of substrate characterization, the dominant VFA is acetic acid and the stability threshold of 2400 g m -3 suggested by Wang et al. [41] was considered in this study.
The VFA accumulation had an inhibitory effect after the 14th day since the S2 values were above the critical value.Following the increase in S2, the pH dropped from 7.84 to 6.55.Considering the finding of Weiland [42] and Kamar and Li [43] who opined that the microorganism in charge of methane production is highly sensitive to pH and flourishes well when the pH lies between 7.0 and 8.0, then, the AD plant was under stable operation up to 14th day of the experiment.The pH changes were small, even when there was a disturbance in the system as a result of VFA accumulation.This could be ascribed to the high buffering capacity provided by cattle manure.To authenticate the condition of the system, the AR values were calculated and presented in Fig. 6b.
The buffering capacity in the digester could not withstand VFA concentrations after 14th day thereby causing instability in the system (AR > 0.3) [39].The instability was reflected in Fig. 4a

3..6 Sensitivity Analysis
Sensitivity analysis is required to make a model-based prediction of the digester's effluent properties by changing some process parameters which in reality undergo variation.This will provide a better understanding of how a change in model output could be ascribed to the change in model input.With the calibrated model, the total biogas production and methane content with ±20% fluctuations of initial values of the state variables and input variables are presented in Table 6.Some of the process parameters considered for sensitivity analysis include inlet biodegradable fraction (S1in, S2in, Cin, Zin), initial biomass concentrations (X1 and X2), and initial substrate concentrations (S1init and S2init).As indicated in Table 6, the change in the initial concentration of VFA in the digester (S2init), the concentration of VFA in the substrate, (S2in), and the initial concentration of the biomasses in the digester (X1init and X2init) led to remarkable changes in both biogas production and methane content.These are variables that are directly related to the conversion of VFA to CH4 and CO2.It can be inferred that the variables linked to methanogenesis are the most sensitive.The variables such as Cin and Zin which determine the CO2 flow and pH of the digester display a lower degree of sensitivity.S1init and S1in were discovered to be insensitive to the model output.It was shown in Table 6 that biogas production decreased when S1in, S2in, S1init and S2init increased.An increase in S1in and S1init could ameliorate the acidogenesis step thereby increasing the S2 formation.Therefore, at constant buffering capacity, an increase in S1in, S2in, S1init, and S2init might lead to more accumulation of S2 in the digester beyond the tolerance level (AR > 0.3).This might attack methanogens leading to a reduction in biogas output.As expected, there was a decline in biogas production when X1init and X2init were reduced.
Lower biogas production occurred as a result of a lower degradation process in the presence of a small number of acidogens and methanogens.An increase in X1init leading to lower biogas output (-60.42%)might be a result of the acidogenesis step being accelerated resulting in more production of VFA and ammonia.The solubilized protein from food wastes and urea from cattle manure could also combine to produce ammonia.A high concentration of ammonia and VFA beyond the allowable limit could contribute to the inhibition in the digester by attacking the methanogens.The rate of biodegradation of VFA would increase to produce more biogas when X2init was increased.

1 .Fig. 1
Fig. 1 Schematic diagram of anaerobic digester showing AM2 variables (input and output variables). .(S1, S2, C, Z)in indicate the variable compositions in the inflow substrate while (S1, S2, C, Z. ) are the concentration of the process variables inside the reactor.The biochemical reactions within the reactor led to the production of biogas in the form of CH4 and CO2.

Fig. 2 Table 2
Fig. 2 The SFD showing different biochemical reactions assumed by the modified AM2 for biogas production in the digester.These diagrams are the results of the transformation of the model equations in the model layer of the STELLA package.(a) acidogenic bacterial balance.(b) methanogenic bacterial balance (c) organic substrate balance.(d) VFA balance.(e) inorganic carbon balance.(f) alkalinity balance.(g) biogas estimation.(h) pH and AR

Fig. 3
Fig. 3 Comparison of the cumulative biogas production of the simulation results to measured values from the

Fig. 5
Fig. 5 Comparison of the simulation of original AM2 and the measured values from an industrial biogas plant of capacity 7000 m 3 with a varying input composition for 28 days.. (a) biogas production and (b) methane concentration.

Fig. 6
Fig. 6 Prediction of the behaviour of stability indicators of an industrial biogas plant fed with cattle manure and food wastes for28 days (a) the variation of predicted S2, Z , and pH (b) variation of predicted AR.

Table 1
Fundamental data used for substrate characterization

Table 4
Biogas production and methane content predicted for the EWE biogas plant with different models

Table 5
Comparison of different models using MAPE

Table 6
Sensitivity analysis of the proposed model to biogas output rate and methane content