332
Biotechnol. hog. 1993, 9, 332-343
Modeling and Simulation of Anaerobic Wastewater Treatment and Its Application to Control Design: Case Whey Gerhard B. Ryhiner,?Elmar Heinzle,*and Irving J. Dunn Biochemical Engineering Group, Chemical Engineering Department, Swiss Federal Institute of Technology, Ziirich, Switzerland
Modeling and simulation of the anaerobic treatment of whey wastewater was studied. A model including five groups of organisms, whey substrate, butyric, propionic, and acetic acids, COZ, Hz, and CHI,was developed for single- and two-stage anaerobic fluidized bed reactors. It included thermodynamic limitations for the acetogenesis reaction, acid-base equilibria of dissociating compounds (organic acids and COZ) for variable and constant p H situations, gas-liquid phase mass transfer, and measurement dynamics. The model was adequate to describe both the steady-state and dynamic behavior of the reactors and was used to study various control schemes. In all cases the feed flow rate was the manipulated variable. Controls of pH, dissolved hydrogen, and organic acids were all successful in simulation and experiment. Simulations were used to determine optimal control parameters for conventional proportional-integral (PI)or proportionalintegral-differential (PID) control. The experimentally found optimal values were in close agreement with predetermined simulated values. Simple control using gas analysis results did not give stable control, either in the simulation or in the experiment.
Introduction Anaerobic processes recently have been used more frequently in wastewater treatment. This is especially true for highly polluted wastewater with a great excess of carbon in relation to nitrogen and phosphorus. The main reasons for the increased application of such processes are as follows: sludge production is reduced approximately 10-fold;energy input is drastically reduced; useful amounts of biogas (methane and carbon dioxide) may be produced with highly loaded wastewater; and air pollution by odorous materials may be reduced by burning the off-gases. There are, however, also some distinct disadvantages. The most significant among them are the following: methane producers grow much more slowly, and the stability of anaerobic methane-producing processes may be worse than that of aerobic processes. Instabilities may be caused either by toxic substrates or by overloading. Some types of disturbances may cause an irreversible collapse of the system. To overcome these problems, anaerobic processes may be controlled to allow prolonged stable operation. To do so, modeling and simulation can aid in understanding process dynamics and evaluating control strategies. As the process is very complex, it is not easy to design a useful controller. Several control strategies are possible, whose performances can only be evaluated either by simulation or by experiment. Experimentation is rather time-consuming and may, especially in the case sf anaerobic wastewater treatment, lead to irreversible loss of biological activity. Therefore, it seems useful to minimize experimentation and to replace it as far as possible with simulation studies. This combined approach using simulation and experimentation to develop control strategies and to determine control parameters is studied in this work. More detailed information is contained in the Ph.D. thesis of Ryhiner (1989) and also in a related
* Author to whom correspondence should be addressed.
t Present address: Sulzer-Chemtech, CH-8401Winterthur, Switzerland.
paper on adaptive optimal control of this process (Ryhiner et al., 1992). Various reviews concerning modeling and control of anaerobic wastewater treatment are available in the literature (e.g., Schiirbiischer and Wandrey (1991), McCarty and Mosey (19911, and Heinzle et a1 (1992)).
Reactions and Stoichiometry in Anaerobic Digestion Anaerobic degradation is a very complex multisubstrate, multiorganism process. It occurs in basically four steps as can be seen in Figure 1. In a first step, polymer materials like carbohydrates, proteins, or lipids are hydrolyzed to yield monomer compounds such as amino acids, sugars, and fatty acids. In a second acidogenic step, these compounds are fermented to organic acids (mainly acetic, propionic, and butyric acids) and hydrogen. In the third acetogenic step, organic acids with more than three atoms of carbon per molecule are converted to acetic acid and hydrogen. In a last methanogenic step, the intermediates, acetic acid, hydrogen, and carbon dioxide, are converted to methane. In Table I the most important reactions are summarized. Nitrogen and sulfur are not included in this stoichiometric model. Three typical cases may be distinguished. If particulate matter has to be degraded, then the hydrolysis step may be rate limiting. The kinetic analysis of such a process may be rather simple. Overload of the reactor may occur if easily degradable material having no toxic compounds has to be treated. It seems to be evident that hydrogen plays an important role in the control of such a process. Due to thermodynamic reasons, reactions R2 and R3 in Table I can, in the presence of acetic acid, only proceed to the right side, if hydrogen partial pressure is sufficiently low. Experimental measurements of dissolved hydrogen in the low concentration range of interest are difficult to obtain. According to Thauer et al. (1977) and Archer (1983), the oxidation of propionic to acetic acid thermo~ 10"' bar. dynamically should only be possible at P H < Experiments of Kaspar (1977) and Denac et al. (1986)
875&7938/93/3009-0332$04.00/0 0 1993 American Chemical Society and American Institute of Chemical Engineers
333
Blotechnol. Rog., 1993, Vol. 9, No. 3
as a measured variable for control. It does, however,mean that the response level of Hz to a disturbance can vary depending on the contribution of both electron-transfer reactions. From a series of observations it seems to be evident that disturbances in the methanogenic step, which generally is considered to be the most sensitive one, lead to accumulation of Hz (McCarty and Smith, 1986; Harper and Pohland, 1986; Whitmore et al., 1985; Archer, 1983). Additionally, the state of an anaerobic reactor may be characterized by its wlatile fatty acid levels, from the CH,/C02 ratio and from the total gas production rate. Elemental balances for carbon hydrogen and oxygen were set up to estimate the stoichiometric coefficients of reaction R1:
Figure 1. Reaction scheme of anaerobic degradation: Poly, polymer material (proteins, fats, hydrocarbons, etc.); X H Y , biomass hydrolyzing polymer materials; Mono, monomeric materials from hydrolysis of polymer materials; X A G , acid generating biomass; HPr, propionic acid; Pr-, propionate; XI+, biomass growing on propionate; HBu, butyric acid; Bu-, butyrate; XB", biomass growing on butyrate; HAC,acetic acid; Ac-, acetate; XAc,biomass growing on acetate; XH,biomass growing on hydrogen and carbon dioxide. Dashed arrows indicate gaseous compounds transfer to the liquid phase. The solid arrows symbolize gaseous compounds transferred to liquid-gas phase. Table I. Most Important Reactions in the Anaerobic Degradation Process of Whey" Hydrolysis and Acidification
-
+ 0.095HZO 0.5CH,00,, CH,,8500,853 (whey) + 0.12CH,00,,,, +0.15CHZO +O.23CO2 +0.24H2
(R1)
Acetification CH,,(CH,),COO-+ 2H,O s 2CH,COO- + 2H,
+ H+
AG', = +76.7 kJ (R2)
CH,,CH,COO-+ 2H,O s CH,COO-
+ CO, + 3H, AG', = +48.1 kJ (R3)
Methanogenesis CH,,COO-+ H,O CO,
i=t
CH,
+ HC0,-
+ 4H, s CH, + 2H,O
0:
0.853 + v H , ~= 1/2vB,
CH,(CH,),COO-
+ H+ + 2H,O 2CH,
(R5)
a Reaction R1 is given on a C-mol basis (one carbon atom per molecule) (AG'o for T = 25 "C, pH 7, p = 1atm, pure water, C = 1 mol kg').
applying pure H2 did not show any significant inhibition by high external concentrations of H2. Kaspar (1977) and, more recently, Whitmore et al. (1987), who actually measured dissolved Hz down to the 1pM range, pointed out that this may be caused by diffusion-reaction phenomena. A locally high consumption rate of H2 in a biofilm or floc may reduce local H2 concentration to allow endergonic reactions R2 and R3 to proceed. Labib et al. (1992) observed only slight inhibition of butyrate conversion to acetate, when Hz was added. This only moderate inhibition was attributed to mass-transfer effects of H2. Boone et al. (1989) claimed that interspecies electron transfer was mainly via formate and not via H2. These findings do not touch on the potential importance of H2
(3)
where v are the stoichiometric coefficients from eq R1 on a C-mol basis. One Csmol has only 1mole of carbon. The six unknown stoichiometric coefficientsof the acidification of whey were determined after the measurements of acid production and COD reduction in the first acidification reactor of the two-stage system (Ryhiner, 1990). In this reactor no methane was produced. Reactions R2 and R3 are assumed not to proceed under these conditions. Accumulation of hydrogen causes a shift of thermodynamic equilibria to the reactands propionic acid and butyric acid. The coefficientshave the dimensions Csmol (C-mol)-' when the reaction contains carbon and mol (C.mol)-l when it does not contain carbon. The experimentally determined values are given in Table I. Kinetic studies from McCarty and Smith (1986) showed that 70% of methane is produced via reaction R4 and 30% via reaction R5, given in Table I. A similar estimation was made for this system. Combination of eqs R2-R5 gives
AGto = -31.0 kJ (R4)
AG', = -130.7 kJ
+ 2/3vp, + vAC + 2vCo2
CH,CH,COO-
-+
+ 2C0,
2H, (R6 = R2+R4)
+ + + - +
+ H+ + 2H20 CH, c
CH,COOCO,
+ 4H,
2C0,
H+
-
3H, (R7 = R3+R4)
CH,
CH,
CO,
+ 2H,O
(R4)
(R5)
Use of the assumption that all carboxylic acids and all HP molecules react to methane allows one to estimate the fractional contribution~ofreactions R4 and R5 for methane formation. The numerator of eq 4 is derived from reactions R4, R6, and R7 by considering the amount of each acid necessary for the formation of 1 mol of methane. The denominator is derived from the stoichiometry of reactions
334
Bbtechnol. Prog., lQQ3,Vol. Q, No. 3
Table 11. Stoichiometric and Kinetic Parameterse substrate lmax (h-]) KS (C-molm-3) K I (C-molm-3) 0.36 0.7 whev 0.25 0.25 0.40 0.000085' 0.17 butyric acid 0.0113 12.5 0.0096 0.33 0.0154 0.032'3d 0.011 0.0008' 0.22 propionic acid 0.0108 21.4 0.0033 1.6 0.013 0.003 0.005d 0.0074C3d 0.003c 0.33 acetic acid 0.0145 0.0142 15.6 0.0142 5.16 0.1c.e 1.4C8e 0.01e 0.1' 1.4' O.OOSd 0.0375 hydrogen 0.0583 0.038 0.0583 0.001 0.0002 0.001 0.058 e
YX,s(Cemol (Cvmol-1)) 0.10b 0.03 0.03 0.013 0.051 0.083 0.05 0.023 0.023 0.06 0.06 0.04 0.043 0.041 0.0567 0.05 0.021 0.0283 0.03 0.026
kd (h-1)
ref Kissalita et al.. 1989 Mather, 1986 0.0004 used here Dalla Torre and Stephanopoulos, 1986 0.0004 Bryers, 1985 0.00011 Gujer and Zehnder, 1983 0.0004 used here Dalla Torre and Stephanopoulos, 1986 0.0004 Bryers, 1985 0.0004 Gujer and Zehnder, 1983 Mather, 1986 0.0004 used here 0.018 Dalla Torre and Stephanopoulos, 1986 0.00013 Bryers, 1985 0.0006 Gujer and Zehnder, 1983 Wiesmann, 1988 0.0004 used here 0.0004 Bryers, 1985 0.0004 Gujer and Zehnder, 1983 Mosey, 1983 Mather, 1986 0.0004 used here
0 Carbon molecular mass of biomass: 29.5 g (C.mol)-l. Assumption 50% content of C in biomass. Free acid. Experimentally determined. Mean value.
R5-R7. (4)
8 ' 4 ' 4
Dynamic Liquid Phase Balance Equations The reactors used in this study were fluidized bed biofilm reactors. The flow characteristics of these reactors are nearly plug flow, but because of the high recirculation rate of the reactor system and the corresponding small differences between inlet and outlet concentrations, they can be considered well-mixed. The basic balance equation for a well-mixed liquid phase is (5)
where V Lis the liquid volume (m3),FL is the liquid flow rate (m3 h-l), CL, is the concentration in the feed stream (C-molor mol m-3), CL,is the actual concentration in the liquid phase (Cmmol or mol m-3), ri is the reaction rate (C-mol or mol m-3 h-I), and Ni is the mass-transfer rate from the gas into the liquid phase (mol m-3 h-l). The amount of organic substances is always given in C.mol(1 mol of carbon in 1 Camel). Biomass. It was assumed that the biomass was homogeneously suspended in the reactor. The biomass retention in the biofilm was modeled using a biomass retention time (TX) which was always much larger than the hydraulic retention time. The resulting biomass balance is
biomass growing on substrate i (C-molm-9, pi is the specific growth rate on substrate i (h-l), k d , is the death rate for organism growingon substrate i (h-l), and Txis the biomass retention time (h). Individual biomasses were not measured. These were determined by taking kinetic parameter values from the literature as given in Table 11. Then individual biomass concentrations were estimated to obtain the best fit of simulated data with experimental data. T X was adjusted to fit effluent solid TOC attributed to biomass (Table 111). Substrates, Products, and Intermediates. Since substrates in the waste treatment system discussed here are not volatile, eq 5 can be simplified to
Each substrate may be consumed by several reactions, which is taken into account by the summation. The relation of substrate consumption to biomass formation is
,,I
-cL,xi
(7)
rs, = (9) yx,/s, where Y x , ~is, the yield of biomass (Xi) produced from the substrate (Si)(C-mol (C-mol)-l). To be able to closely solve the carbon balance, it was assumed that all of the carbon released during the death or maintenance processes remained in the system as C02. A substrate requirement for maintenance was not explicitly included in the model. But as all carbon released by death was modeled to create CO2, it is equivalent to using a maintenance term. Endogenous metabolism is also included in this model. A similar formulation was used by McCarty and Mosey (1991). The corresponding balances for intermediates (Hz, carboxylic acids) and products (COZ,CHd are also derived from eq 5 and include the mass-transfer term for the volatile compounds:
where i indicates the various substrates for biomass formation (whey, butyric, propionic, and acetic acids, hydrogen) (C-mol m-3), rx, is the rate of production of biomass growing on substrate i (C-mol m-3 h-l), Xi is the
Cri VLalso includes the rates for the formation of a product
where
Wtechnol. Prog, 1993, Vol. 9, No. 3
395
Table 111. Physical and Thermodynamic Parameters Darameter KD,,,,
KD, ,
z = kcations TX
descrbtion
value
dissociation constante
1.54 X mol m-3 1.34 X mol m-3 1.76 X mol m-3 3.02 X mol m-3 3.0 mol m-3 60-1440 h 360 h 24.5 mol m-3 bar1 1.12 mol m-3 bar-' 0.74 mol m-3 bar-' 2.0 h-l 0.8 m 100 m h-l
total concentration of cations biomass retention time henry Coefficients mass-transfer coefficient reactor diameter
UB
R PH,O
T
P TL
VGH
source
bubble rising velocity gas constant partial pressure of water at 310 K temperature total pressure time constant liquid analysis reactor headspace
Weast. 1987
measured Hall,1987 this work Mather, 1986 this work this work Aiba et al., 1973
8.314 X bar m3mol-' K-l 0.07 bar 310 K 1 bar 0.3 h 0.2 x 10-3 m3
Weast, 1987 measured measured determined determined
. .-
j from substrate i in the following form:
1 .o
I
(11) 0.8
with Yp,,, as the stoichiometric coefficient for formation of product (Pj)from substrate (Si)(C-mol(C.mol)-l) in the corresponding reaction. These values are identical to vi in R1. In R2, YAc/Bu = 1 and YH2/Bu = 0.5. In R3, YAc/Pr = 0.5, YHr/pr = 1, and YC02/pr = 0.333. In R4, YCH4/Ac = 0.5 and Y c o ~ / A ~0.5, and in R5 YCH4/COn = 1 and YCH4/H2 = 0.25.
Kinetics The kinetics of biomass growth are described by simple Monod relationships for whey, butyric acid, propionic acid, and hydrogen:
or with an additional substrate inhibition term (acetic acid)
where pi," is the maximal specific growth rate (h-l), Ksi is the Monod constant for the half-maximal growth rate (C-mol m-3), KI,is the inhibition constant (C-mol m-3), and EFi is the equilibrium factor. EFp, and EFB, are defined in Figure 2 and eqs 14-19 (EFH,= 1 and EFwy = 1). S and P are always related to one individual reaction, Le., the product in one reaction may be the substrate in another one. In many reactions the concentration of the actual substrate (usually the free acid form of carboxylic acid) is pH dependent. In some papers the acidic form is considered to be the kinetically determining form (Mtirkl et al., 1983; Dalla Torre and Stephanopoulos, 1986). This procedure was also followed in this work. Literature values for stoichiometric and kinetic coefficients and data used for simulation are summarized in Table 11. The stoichiometric coefficients are given in C-mol. In the acetogenic step, acetic acid, hydrogen, and carbon dioxide are produced from propionic and butyric acids. In additional to the kinetic data given in Table 111, the thermodynamic limits for these reactions are incorporated
LL
u
0.6 0.4 0.2
0.0
".- ,
;
1
0
3
K'IK
Figure 2. Equilibrium factors (EF) for reactions from eqs 14 and 15: K*, values from eq 18 or 19;K, equilibrium constant; EF, equilibrium factor to be used in eq 12.
by estimating the chemical equilibrium limits for butyric acid: CH,(CH,),COO- + 2H20 e 2CH3COO- + 2H, + H+ AG; = 48.3 kJ (14) From this one obtains an equilibrium constant of K B=~ 2.02 X 10-3 mol4 m-12.
For propionic acid we similarly obtain CH,CH,COO-
+ 2H,O
e 2CH3COO- + 3H,
+
AG; = 76.1 kJ (16) CO, and a equilibrium constant of Kp, = 1.35 X 10-l2 mol4 m-12.
The factor 4i3 is necessary because concentrations here are given in C-mol. Observed reactions in the backward direction were very slow (Smith and McCarty, 1990). They are not considered to be important and are therefore not included in the
Biotechnol. hog.., 1993, Vol. 9, No. 3
336
model. When the reaction is far from equilibrium, eq 1 2 describes the kinetics very well, whereas they are not well known closer to the equilibrium. An empirical approach was therefore chosen. As the reactions slowed down on approaching equilibrium, they were not all_owedto proceed to the right-hand side when the equilibrium condition was reached. From the actual concentrations, KB"* and Kpr* were estimated as follows.
KPr*=
4CAc2CH23CC02
3C~r
The fractions KB,*IKB, and Kp,*lKp,, whose values are >1 if the equilibrium has not yet been reached, were determined according to Figure 2. In many literature
models the equilibrium is not incorporated into the model because of inherent problems. An exact formulation of the equilibrium reaction requires the reaction rates of R2 and R3 to be equal in each direction under equilibrium conditions. This can easily be done at equilibrium conditions using eqs 15 and 17 combined with eqs 9 and 12, but there is no information available somewhat farther away from these equilibrium conditions. The experimental data describing the reverse reactions given by Smith and McCarty (1990) indicate that these reactions are not important for the process dynamics because of their slow rates. The most important facts are, therefore, described by the curve given in Figure 2, which is a somewhat arbitrary polygon starting from EF = 0 at K*/K < 1 and rising to EF = 1at K*/K > 2. The reaction to the right stops when the equilibrium is reached, but there is no reverse reaction when concentrations of acetate, hydrogen, COz, and H+ exceed the equilibrium values caused by other reactions. The reaction proceeds in an irreversible manner farther away from the equilibrium (K*/K > 2). The resulting factor EF is used in eq 12. A similar approach was taken by Mather (19861, who also used eq 7 to describe specific growth rate, but set pm,, = KAG. Here AG was the actual free reaction enthalpy at actual concentration values and K an empirical constant. From this, pmaxincreased monotonously with increasing substrate and decreasing product concentrations. No real saturation values for growth rates were obtained. Ion Charge Balance. The dissociation equilibria combined with an ion charge balance yields
where KA,is the acid dissociation constant (e.g., K A ~ KB, ), is the base dissociation constant (e.g., K N H JKW , is the dissociation constant of water, CB,,,, , is the total concenis the tration of base i (i indicates ammonia), and CA,,,,, total concentration of acid i (i indicates COz, acetic, propionic, and butyric acids, and hydrogen sulfide). From eq 20 the pH can be estimated for any situation by solving the nonlinear equation, provided the total concentrations of weak acids (CA,,,, ,), weak bases (CB,,,, ,), cations of strong bases (CK+),and anions of strong acids (CA) are known. From the above relations it is seen that COz is always important. If at 1bar of total pressure 30% of the gas is
COZ,the equilibrium liquid phase concentration will be around 4 mol m-3. At pH values below 7.25, less than 1% of the ammonia nitrogen exists in the free base form. As the total concentration of ammonia is rather low (X1.5 mol m-3, the corresponding terms can be neglected. The H2S-HS- system has its highest buffering capacity at pH 6.9. The experiments described here were typically around pH 6, which reduces the importance of this system for buffering. In whey wastewater the concentration of sulfur is quite low. The measured gas concentrations were usually below 2%. The corresponding HBS concentration was therefore always