- Email: [email protected]
Analysis of a model of hollow-fiber bioreactor wastewater treatment
Chemical En@werh# Science, Vol. 48. No. 6, pp. 1069-1075, Printed in Great Britain.
ANALYSIS
1993.
OF A MODEL OF HOLLOW-FIBER WASTEWATER TREATMENT R. N. PAUNOVI6,
Faculty
am-2509/93 s6.00 + om 0 1993 Pqamon Press Ltd
of Technology,
Z. Z. ZAVARGO
and M. N. TEKIC’
University of Novi Sad, Bul. Avnojs
(Firsr received 7 January
1992 accepted
in
BIOREACTOR
I, Novi Sad 21000, Yugoslavia
revisedform 28 May 1992)
Abstract-h this work, the effect of neglecting fiber lumen radial concentrations gradients in the onedimensional mathematical model of isothermal anaerobic acetate biodegradation in hollow-fiber bioreactor has been analyzed. The results show a significant increase of the error resulting from neglection of radial lumen concentration gradients with increase of hollow-fiber diameter value. The limiting fiber diameter value above which the inclusion of lumen concentration profile is recommendable is proposed. To provide a stable convergence in the whole domain of real model solutions, a self-starting procedure based on a shooting method has been developed.
INTRODUCTION
Modeling is an important step in developing biochemical reactors, not only for experimental data analysis but also for scaling-up purposes. Between the demand to include many hydrodynamic and kinetic parameters and the availability of the mathematical solutions of practical importance, the compromise is always to be made when constructing the model. Recently, Dall-Bauman et al. (1990) presented a one-dimensional model describing steady-state isothermal anaerobic biodegradation of acetate in a hollow-fiber bioreactor (HFBR). The authors have neglected the radial concentration gradients in the region of fiber lumen. The errors introduced by such an assumption could become significant at larger hollow-fiber diameters and high biofilm densities. In this paper, an analysis based on comparisons of the numerical solutions of the model proposed by Dall-Bauman et al. (1990) and the revised model which includes the radial concentration gradients in fiber lumen is presented.
In the model proposed by Dall-Bauman et ai. (1990) (model l), the radial concentration gradients in phase 1 are also neglected, i.e.
6”
= c,,
O
(la)
where r. = aR is the inner fiber radius, and C, is the concentration value at r0 = 0. In the revised model (model 2) a parabolic concentration profile in phase I is assumed: C”’ = A,?
+ C o,
OGrbr,
WI
and the boundary condition is C(l) = Co
at r = 0.
The parameter AI can be related to the substrate flux, F,, at the interface r = ro:
The concentration profile in phases 2 and 3 is described by the following differential equations: dCc2’
V-X dr
(3)
MODELS
Assuming that a module consists of N identical and independent hollow fibers, a single fiber shown in Fig. 1 could be considered as being representative. The fiber feed, qf, is, thus, equal to the one-Nth of the module feed, Q,. Each fiber can be considered to be composed of three phases (Dali-Bauman et nl., 1990); the fiber lumen (l), the fiber wall (2) and the biofilm (3). The parameter p of the biofilm thickness is determined by the shell-side volume, the effective reactor length and the number of fibers. At high recycle ratio, p, values, the axial concentration gradients in the substrate could be neglected (Dall-Bauman et al., 1990).
-
kX,Cr3’ K, + c’“”
(4)
with the following boundary conditions: At r=ro, C”’ = C’Z’
@a)
(5b) At r = R,
0’3’
‘Author to whom correspondence
R
should be addressed.
-.
dc’3’ dr
1069
NW
R.
1070
N. PAUNOVIC
et al.
Assuming the parabolic axial velocity profile u=“o[l
-(JJ]
one obtains from eq. (14):
co
c* =
riA,/3
-1; Co
for model 1 for model 2.
(15a) (15b)
ALGORITHM
The problem of determination of the concentration profile in a hollow-fiber reactor is a boundary value problem. The methods of orthogonal collocation (OC) and orthogonal collocation on finite elements (OCEF) proved to be efficient to solve diffusion-reaction problems (Finlayson, 1980). Dall-Bauman er al; (1990) applied the OCFE method for the solution of the hollow-fiber reactor model, but, due to the nonconvergence of the Newton-Raphson technique applied for the solution of the resulting nonlinear system of equations, they had to modify the procedure by decoding the system. The material balance equation is solved for the etIIuent concentration, CO, in the outer iteration loop, and the remaining equations in the inner loop. By decoupling, the potential rate of convergence of the original OCFE method is signifieantly decreased, approximately by a factor equal to the necessary number of iterations in the outer iteration loop. In addition, to make the convergence easier, the influent concentration, Ci, is used as a specified variable, instead of “natural” specification of the feed concentration, C,, the latter being calculated from
Fig. 1. Schematic diagram of single-hollow-fiberbioreactor. At r = /3R, _
D'3'dCo
= ,-,,
(7)
dr
Seemingly, there is a surplus of boundary conditions, but one should recall that the first two conditions [eqs (6a) and (6b)] are shared with the profile in phase 1. If the radial velocity, D, is approximated as (Dall-Bauman et al., 1990) v=v,-,
KR
ro-crcj7R
r
eq. (3) has the following analytical solution: (-y
=
-
&,
f’@’ + B2
and the values of the parameters A2 and B2 follow from the boundary conditions (5a) and (5b):
c, A, = +,
,.6’-14~‘)
Co + Fob,
(10)
for model 1
(1 la)
for model 2
(1 lb)
The recycle ratio, p, introduced by the proposed formulation (Dall-Bauman et al., 1990) has, however, no meaning in a model that ignores the axial concentration gradients. In order to avoid problems accompanying the application of the OCFE method, in this paper a robust (in terms of convergence characteristics) and simpleto-use method for the solution of the two comparative models of the hollow-fiber reactor is proposed. The method is based on the shooting technique (Finlayson, 198Q Carnahan et al., 1969)--an iterative solution of the corresponding initial value problem. A shooting procedure could be formulated as follows. Starting from an estimate of the concentration Co, the flux F, can be calculated from the material balance equation (13):
B2 = Co + :
(1 - 0.5Pe”))
where Per’) and Pe(*) are the Peclet numbers for phases 1 and 2, respectively: Pet’) = v,r,/D(‘)
> i = 1, 2.
(12)
One more constraint should be included in the model and that is the material balance equation: 2nroLFo + u,C(r)(r,)
= q,(C,
The average concentration
- C*).
(13)
in the lumen, C*, is
F
=
0
defined as 9JC * = 2x
OJro(Cf
CO)/~
-
(0.5ro(Cf - C,)/r
-
V&O - o,C,)~
for model 1 for model 2
(16a) (16b)
where
0J s0
= (p + 1)Ci - pco.
u(r)C(l)(r)r
dr
(14)
where u(r) is the axial velocity profile in the lumen.
@= and tg) and
7
1 1 - @‘/(12r)
- 0.5Pe’l
are the diffusion and contact
(17) time,
1071
Analysis of a model of hollow-fiber bioreactor wastewater treatment
volume, V, are
respectively: t’,” = t$D”’ r = L/VI.
(18)
s = 2nRL
(19)
v = R(/32- 1)RZL.
Subsequently, the parameters Al, Az and B, in the analytical solutions (lb) and (9) are calculated from eqs (2), (lo), (11 a) and (11b) and the differential equation (4) is numerically integrated, starting from the initial conditions (6a) and (6b), which take the following form: C”’ = GoCo + G,(C, - Co) dCt3’
~
at r = R
= GzCo - G,(C,- Co) at r = R
dr
(20a) (2’3
After introduction of eqs (20a) and (2Ob), one obtains the following quadratic equation: oC~+bCo+c=O
(26a)
G1)(G2 + G3)Dc3’Cf
a = (Go -
b = (Go - GI )[O.S(p’ - l)RkX, - G3D’3’C,] +(K I/C,+ G,)(GZ+ G&Y3'C / c = 0.5(/Y*-
l)RkX,G,
- (K./C1
WC)
+ GI)GJD%,-.
where
(26d)
Go =
c1- pe’2’
=
- rg’/lZr]@
P&2
[a-
1 I
-030
G. -1
(2W
(1
-
for model 1 (21a) for model 2 (21b)
In the domain of realistic parameter values, eq. (26) has two real roots of different signs, the positive one having a physical meaning. Thus, the solution of the
for mode1 1 (22a)
tl-Pe”‘)
VW5
02
(1
. _
_
O-Pe”’
- O.We”))0
for mode1 2
(22b)
original problem obeys the relation
u,
for model 1
(_p -Pe”’
0’3’
GZ =
v,
El --Pe”’0
0’3’
for model 2
C,, > (-
(23a) (23b)
b + ,/=)/2a.
dC” ---CO dr
0.5r0 &_l _
-PIP’
D’%
for model 1
(244
for model 2.
(24b)
(27)
A closer upper bound of the initial interval may be obtained from the condition atr=R
resulting in [after introduction of eq. (2Ob)]
G3 = 0.5ro fi=
i -ppm
Q,
If the boundary condition (7) is satisfied, the problem is solved. Otherwise, it is necessary to repeat the procedure with a new guess for Co. Starting from a correct interval estimate for C,,, one might apply the half-interval method or the Regula falsi method to converge to the solution. A possible starting interval for the described iteration procedure is (0, C,), but it is a very rough estimate and, from the view point of computational time consumption, a closer one would be preferable. For the lower limit of CO, the solution of the flat concentration profile case in phase 3 (no resistance to diffusion) could be adopted, as follows. For the case of uniform concentration in phase 3, the substrate flux from phase 2 to phase 3 satisfies the relation
where
L)‘*’
dc’*’ dr
-
S =
kXfC”’ V K, c 03’
atr=R
(25)
the interfacial surface area, S, and the biofilm
The computer algorithm can be summarized as follows:
(1) Input all the specifications:
(a) The design parameters: R, a, B, L, X/, 2 and
Cf (b) The transport and kinetic parameters: u,, D”’ (i = 1, 3), k, K,. (2) Define the initial interval for Co [eqs (27) and
WI.
(3) Convert the second-order differential equation (4) into two first-order equations and intcgrate them by fourth-order Runge-Kutta method. An automatic adjustment of the integration step by successive bisections to get the prescribed accuracy is recommended. (4) Check the boundary condition (7). If necessary, correct Co and return to step 3. Otherwise continue. (5) Determine the average effluent concentration and the extent of biodegradation of acetate.
R. N.
1072
PAUNOVIC et al.
For correction of Co values, the half-interval method is used until the interval [C,.,j,,. C,,.,,,J, which satisfies the following conditions is located: For
CO
=
Grnin
GO
atr=flR
(29)
at r = DR.
(30)
For Co = CO.mS.:
Afterwards, the secant method is applied to obtain Co with the prescribed convergence criterion. A real solution to the problem (positive’ concentration values) need not exist throughout a wide region of specification values. A necessary and sufficient condition to obtain a real solution is the existence of an interval CCO.min,CO, mi ] that satisfies conditions (29) and (30). The absence of a real solution for the formulated models indicates that the substrate is completely consumed in a part of biofilm, not reaching its outer border. To describe such a case, one additional variable is to be introduced: the distance r* (R -c r* c BR) up to which the substrate penetrates through the biofilm. At the “moving” moundary r*, the following two conditions are valid: C’S’ _ -
_dcY’ dr
at r = r*_
= 0
This extension of model is not included in the existing computer program.
RESULTS AND DISCUSSION
is obvious that the introduction of parabolic concentration profile in phase 1 (model 2) affects the boundary conditions at the interface between phases 2 and 3. Such a “perturbation” of boundary conditions is a function of the following dimensionless groups: the ratio of diffusion to contact time [eqs (18) and (19)], the Peclet numbers for phases 1 and 2 [eq. (12)], as well as of the geometric parameters Q. The effect of perturbing the boundary conditions on the predicted average effluent concentration depends on the parameters of the dimensionless biofilm model: It
(31) where R ‘=+ y = @/c, Pd3)
= u,ro/D~5)
(314 @lb)
(31c)
(314
To conclude, the effect of neglecting the radial coucentration gradients in the reactor lumen is dependent on the model geometric parameters o! and 8, simplexes tb”/r and KS/C, and four dimensionless groups, Pet) (i = l-3), Da”. To perform a quantitative analysis for design purposes, the influences of some selected most important single design parameters (instead of dimensionless groups) is studied. The included design paratIEters are as follows: inner/outer fiber diameters, biofilm density, X,-, feed concentration, C,, contact time, 7, and biofilm thickness parameter, 8, we qualitatively predict the influences of these single dimensionless parameters in the light of their physical meanings. In view of the significance of the unremoved quantity of polluting substance (acetate), the relative deviation between the average reactor effluent concentrations predicted by models 2 and 1 6 = cp/c;
(32)
is adopted as a measure of error resulting from neglecting the lumen radial concentration gradients. The results of the numerical analysis are shown in Figs 2-5. All the calculations were performed with the following transport and kinetic parameter values (Dall-Bauman et al., 1990): V, = 2.2 cm/d, D(l) = Df2) = 1.09 crn2/d, Dc3) = 0.87 cm3/d, k = 2.5 mg acctatejmg cells/d, K, = 0.01 mg/cm’. Figure 2 shows the relative deviation of the predicted effluent concentrations, 6, as a function of the inner hollow-fiber diameter, &, and biomass density, X,. It is clear that the error resulting from the neglect of the radial lumen concentration gradients becomes significant as the inner fiber diameter increases and this trend is more pronounced for higher biomass density. Tolerating the differences between effluent acetate concentrations predicted by models 1 and 2 up to 20% (6 4 1.2), it appears that one should take into account the radial lumen gradients for inner fiber diameters above 0.03 cm, at biofilm density of 1 mg/cm3. The values of the remaining design parameters are indicated in the figure. At a biofihn density of 100 mg/cm’, the limiting value of the fiber diameter is, however, reduced to about 0.01 cm. Figure 3 illustrates a decrease of d with increase in feed concentration, CJ, and shows that the limiting hollow-fiber diameter is almost insensitive to the C, variation. The influence of biofilm thickness is presented in Fig. 4. The curves /I = 4 and fi = 2 tend to intersect, i.e. the endpoint (boundary of the real model solution) of the curve /I = 4 lies approximately on the /? = 2 curve. Quite a similar trend is exhibited by the curves B = 2 and B = 1.3. As far as the limiting fiber diameter vaues are concerned, the influence of /I is similar to that of XJ. Further, one can see that the curve 0 = 1.3 has a maximum at the fiber diameter of approximately 0.22 cm. This could be explained by the influence of the decreasing extent of acetate biodegradation with the increase of do. Namely, as one can see from
Analysis
1.00
0.00
of a model of hollow-fiber
0.65
I
0.10
bioreactor
wastewater
I
0.15
treatment
1073
I
0.20
0
5
&/cm/
Fig. 2. Relativedeviationbetweenthe averaage effluentconcentrationpredictionsas a functionof the inner hollow-fiber diameter and the biofilm density.
l.00
0.00
.‘“,.‘..,“.‘,‘I’l,‘II’ 0.05
0.10
0.15
0.20
0
‘5
%/cm/
Fig. 3. Relativedeviationbetweentheaverageeffiuentconcentrationpredictionsas a functionof the inner hollow-fiberdiameterand the feed concentration.
Fig. 5, the relative deviation, S, decreases with decrease of the extent of biodegradation. In Figure 5, the extent of acetate biodegradation predicted by the two models is presented as a function of residence time, t. The difference between the extents of conversion predicted by models 1 and 2 appears to decrease with increase in t, while the relative deviation between the predicted effluent concentrations increases. The computing algorithm proved to be reliable and provided convergence in the whole range of the tested model parameter values, in which the models have real solutions. However, one should expect a con-
siderably higher computing time consumption in comparison to a method based on OC or OCFE techniques. CONCLUSIONS
The error resulting from neglection of the radial lumen concentration profile increases with increase of inner hollow-fiber diameter. This trend is more noticeable for higher values of biofilm density and thickness. As a limiting value, an inner fiber diameter of 0.01 cm may be suggested, above which the inclusion of the radial concentration profile in the fiber lumen model is recommendable.
R. N. PAUNOVIC et
1074
I 0.05
o.io
o.i5
al.
o.lo
I
0.25
0.30 dJcm/
i5
Fig. 4. Relative deviation between the average effluent concentration predictions as a function of the inner hollow-fiber diameter and the biofilm thickness: (0) endpoint of the curve @= 4; (A) endpoint of the curve fl = 2.
0.90
0.70 0.50 0.30
Fig. 5. Deviations
between the predictions of models 1 and 2. 6: relative deviation of the average effluent concentrations; X, : extent of acetate biodegradation predicted by model 1;x2: extentof acetatebiodegradation predictedby model 2.
The formulated computing algorithm for bioreactor simulation proved to be reliable, assuring stable convergence in the whole range of the tested parameter values, in which real solutions of the formulated models exist. The algorithm is self-starting, i.e. automatic generation of an initial guess for the reactor effluent concentration is included. NOTATION Al
parameter in parabolic substrate concentration profile in fiber lumen defined by eq. (1b)
A2 a B* b c; C”’
integration parameter in the analytical solution for concentration profile in the fiber wall parameter in quadratic equation (26) integration parameter in the analytical solution for concentration profile in the fiber wall parameter in quadratic equation (26) average substrate effluent concentration predicted by model i, mg/cm3 substrate concentration in phase i, mg/cm’
Analysis
c* Cf
co C D(i)
Dan
do Fo
Gi k
KS L N PC?@)
Qf q/ R r
SO
S
us
46:6-E
of a model of hollow-fiber bioreactor
average substrate concentration in the fiber lumen, mg/cm3 feed substrate concentration, mg@n3 substrate concentration at fiber axis mg/cm3 parameter in quadratic equation (26) diffusion coefficient of the substrate in phase i, crn2/d Damkohler second number, defined by eq. (31d) inner fiber diameter, cm substrate diffusion flux between the fiber lumen and the fiber wall, m&m2 d parameters defined by eqs (21a)-(24b) maximum rate of substrate utilization, mg substrate/mg cells/d Monod half-velocity coefficient, mg/cm” effective length of the fiber, cm number of fibers in the bioreactor Peclet number for phase i defined by eqs (12) and (31~) reactor feed rate, cm”/d hollow-fiber feed rate, cm3/d outer fiber radius, cm radial coordinate inner fiber radius, cm surface area of the interface between the fiber wall and the biofilm, cm’
V V
ur u, Xf
wastewater
treatment
1075
biofilm volume, cm3 radial velocity, cm/d average feed velocity, cm/d wall Permeation velocity, c/d biomass density, mg cells/cm3
Greek letters ratio of inner to outer fiber radius ratio of outer film radius to outer fiber ; radius relative deviation between the reactor aver6 age effluent concentrations predicted by models 2 and 1, defined by eq. (32) T contact time, d @ dimensionless factor defined by eq. (17)
REFERENCES
Carnahan,B., Luther,H. A. and Wilkes, J. O., 1969,AppIied NumericalMethods, PP. 405-415. Wiley, New York. Dali-Bauman,L., Ilias,S. and Govind, R.: 1990,Analysis of hollow fiber bioreactor wastewater treatment. Biotechnol. Bioengng 35, 837-842. Finlayson, B. A., 1988,Nonlinear Analysis in Chemical Engineering, pp. 73-164.McGraw-Hill, New York. Ralston,A. and Wilf, H. S., 1960,Mathematical Methods for Digital Computers, pp. 110-120. Wiley, New York.
ANALYSIS
1993.
OF A MODEL OF HOLLOW-FIBER WASTEWATER TREATMENT R. N. PAUNOVI6,
Faculty
am-2509/93 s6.00 + om 0 1993 Pqamon Press Ltd
of Technology,
Z. Z. ZAVARGO
and M. N. TEKIC’
University of Novi Sad, Bul. Avnojs
(Firsr received 7 January
1992 accepted
in
BIOREACTOR
I, Novi Sad 21000, Yugoslavia
revisedform 28 May 1992)
Abstract-h this work, the effect of neglecting fiber lumen radial concentrations gradients in the onedimensional mathematical model of isothermal anaerobic acetate biodegradation in hollow-fiber bioreactor has been analyzed. The results show a significant increase of the error resulting from neglection of radial lumen concentration gradients with increase of hollow-fiber diameter value. The limiting fiber diameter value above which the inclusion of lumen concentration profile is recommendable is proposed. To provide a stable convergence in the whole domain of real model solutions, a self-starting procedure based on a shooting method has been developed.
INTRODUCTION
Modeling is an important step in developing biochemical reactors, not only for experimental data analysis but also for scaling-up purposes. Between the demand to include many hydrodynamic and kinetic parameters and the availability of the mathematical solutions of practical importance, the compromise is always to be made when constructing the model. Recently, Dall-Bauman et al. (1990) presented a one-dimensional model describing steady-state isothermal anaerobic biodegradation of acetate in a hollow-fiber bioreactor (HFBR). The authors have neglected the radial concentration gradients in the region of fiber lumen. The errors introduced by such an assumption could become significant at larger hollow-fiber diameters and high biofilm densities. In this paper, an analysis based on comparisons of the numerical solutions of the model proposed by Dall-Bauman et al. (1990) and the revised model which includes the radial concentration gradients in fiber lumen is presented.
In the model proposed by Dall-Bauman et ai. (1990) (model l), the radial concentration gradients in phase 1 are also neglected, i.e.
6”
= c,,
O
(la)
where r. = aR is the inner fiber radius, and C, is the concentration value at r0 = 0. In the revised model (model 2) a parabolic concentration profile in phase I is assumed: C”’ = A,?
+ C o,
OGrbr,
WI
and the boundary condition is C(l) = Co
at r = 0.
The parameter AI can be related to the substrate flux, F,, at the interface r = ro:
The concentration profile in phases 2 and 3 is described by the following differential equations: dCc2’
V-X dr
(3)
MODELS
Assuming that a module consists of N identical and independent hollow fibers, a single fiber shown in Fig. 1 could be considered as being representative. The fiber feed, qf, is, thus, equal to the one-Nth of the module feed, Q,. Each fiber can be considered to be composed of three phases (Dali-Bauman et nl., 1990); the fiber lumen (l), the fiber wall (2) and the biofilm (3). The parameter p of the biofilm thickness is determined by the shell-side volume, the effective reactor length and the number of fibers. At high recycle ratio, p, values, the axial concentration gradients in the substrate could be neglected (Dall-Bauman et al., 1990).
-
kX,Cr3’ K, + c’“”
(4)
with the following boundary conditions: At r=ro, C”’ = C’Z’
@a)
(5b) At r = R,
0’3’
‘Author to whom correspondence
R
should be addressed.
-.
dc’3’ dr
1069
NW
R.
1070
N. PAUNOVIC
et al.
Assuming the parabolic axial velocity profile u=“o[l
-(JJ]
one obtains from eq. (14):
co
c* =
riA,/3
-1; Co
for model 1 for model 2.
(15a) (15b)
ALGORITHM
The problem of determination of the concentration profile in a hollow-fiber reactor is a boundary value problem. The methods of orthogonal collocation (OC) and orthogonal collocation on finite elements (OCEF) proved to be efficient to solve diffusion-reaction problems (Finlayson, 1980). Dall-Bauman er al; (1990) applied the OCFE method for the solution of the hollow-fiber reactor model, but, due to the nonconvergence of the Newton-Raphson technique applied for the solution of the resulting nonlinear system of equations, they had to modify the procedure by decoding the system. The material balance equation is solved for the etIIuent concentration, CO, in the outer iteration loop, and the remaining equations in the inner loop. By decoupling, the potential rate of convergence of the original OCFE method is signifieantly decreased, approximately by a factor equal to the necessary number of iterations in the outer iteration loop. In addition, to make the convergence easier, the influent concentration, Ci, is used as a specified variable, instead of “natural” specification of the feed concentration, C,, the latter being calculated from
Fig. 1. Schematic diagram of single-hollow-fiberbioreactor. At r = /3R, _
D'3'dCo
= ,-,,
(7)
dr
Seemingly, there is a surplus of boundary conditions, but one should recall that the first two conditions [eqs (6a) and (6b)] are shared with the profile in phase 1. If the radial velocity, D, is approximated as (Dall-Bauman et al., 1990) v=v,-,
KR
ro-crcj7R
r
eq. (3) has the following analytical solution: (-y
=
-
&,
f’@’ + B2
and the values of the parameters A2 and B2 follow from the boundary conditions (5a) and (5b):
c, A, = +,
,.6’-14~‘)
Co + Fob,
(10)
for model 1
(1 la)
for model 2
(1 lb)
The recycle ratio, p, introduced by the proposed formulation (Dall-Bauman et al., 1990) has, however, no meaning in a model that ignores the axial concentration gradients. In order to avoid problems accompanying the application of the OCFE method, in this paper a robust (in terms of convergence characteristics) and simpleto-use method for the solution of the two comparative models of the hollow-fiber reactor is proposed. The method is based on the shooting technique (Finlayson, 198Q Carnahan et al., 1969)--an iterative solution of the corresponding initial value problem. A shooting procedure could be formulated as follows. Starting from an estimate of the concentration Co, the flux F, can be calculated from the material balance equation (13):
B2 = Co + :
(1 - 0.5Pe”))
where Per’) and Pe(*) are the Peclet numbers for phases 1 and 2, respectively: Pet’) = v,r,/D(‘)
> i = 1, 2.
(12)
One more constraint should be included in the model and that is the material balance equation: 2nroLFo + u,C(r)(r,)
= q,(C,
The average concentration
- C*).
(13)
in the lumen, C*, is
F
=
0
defined as 9JC * = 2x
OJro(Cf
CO)/~
-
(0.5ro(Cf - C,)/r
-
V&O - o,C,)~
for model 1 for model 2
(16a) (16b)
where
0J s0
= (p + 1)Ci - pco.
u(r)C(l)(r)r
dr
(14)
where u(r) is the axial velocity profile in the lumen.
@= and tg) and
7
1 1 - @‘/(12r)
- 0.5Pe’l
are the diffusion and contact
(17) time,
1071
Analysis of a model of hollow-fiber bioreactor wastewater treatment
volume, V, are
respectively: t’,” = t$D”’ r = L/VI.
(18)
s = 2nRL
(19)
v = R(/32- 1)RZL.
Subsequently, the parameters Al, Az and B, in the analytical solutions (lb) and (9) are calculated from eqs (2), (lo), (11 a) and (11b) and the differential equation (4) is numerically integrated, starting from the initial conditions (6a) and (6b), which take the following form: C”’ = GoCo + G,(C, - Co) dCt3’
~
at r = R
= GzCo - G,(C,- Co) at r = R
dr
(20a) (2’3
After introduction of eqs (20a) and (2Ob), one obtains the following quadratic equation: oC~+bCo+c=O
(26a)
G1)(G2 + G3)Dc3’Cf
a = (Go -
b = (Go - GI )[O.S(p’ - l)RkX, - G3D’3’C,] +(K I/C,+ G,)(GZ+ G&Y3'C / c = 0.5(/Y*-
l)RkX,G,
- (K./C1
WC)
+ GI)GJD%,-.
where
(26d)
Go =
c1- pe’2’
=
- rg’/lZr]@
P&2
[a-
1 I
-030
G. -1
(2W
(1
-
for model 1 (21a) for model 2 (21b)
In the domain of realistic parameter values, eq. (26) has two real roots of different signs, the positive one having a physical meaning. Thus, the solution of the
for mode1 1 (22a)
tl-Pe”‘)
VW5
02
(1
. _
_
O-Pe”’
- O.We”))0
for mode1 2
(22b)
original problem obeys the relation
u,
for model 1
(_p -Pe”’
0’3’
GZ =
v,
El --Pe”’0
0’3’
for model 2
C,, > (-
(23a) (23b)
b + ,/=)/2a.
dC” ---CO dr
0.5r0 &_l _
-PIP’
D’%
for model 1
(244
for model 2.
(24b)
(27)
A closer upper bound of the initial interval may be obtained from the condition atr=R
resulting in [after introduction of eq. (2Ob)]
G3 = 0.5ro fi=
i -ppm
Q,
If the boundary condition (7) is satisfied, the problem is solved. Otherwise, it is necessary to repeat the procedure with a new guess for Co. Starting from a correct interval estimate for C,,, one might apply the half-interval method or the Regula falsi method to converge to the solution. A possible starting interval for the described iteration procedure is (0, C,), but it is a very rough estimate and, from the view point of computational time consumption, a closer one would be preferable. For the lower limit of CO, the solution of the flat concentration profile case in phase 3 (no resistance to diffusion) could be adopted, as follows. For the case of uniform concentration in phase 3, the substrate flux from phase 2 to phase 3 satisfies the relation
where
L)‘*’
dc’*’ dr
-
S =
kXfC”’ V K, c 03’
atr=R
(25)
the interfacial surface area, S, and the biofilm
The computer algorithm can be summarized as follows:
(1) Input all the specifications:
(a) The design parameters: R, a, B, L, X/, 2 and
Cf (b) The transport and kinetic parameters: u,, D”’ (i = 1, 3), k, K,. (2) Define the initial interval for Co [eqs (27) and
WI.
(3) Convert the second-order differential equation (4) into two first-order equations and intcgrate them by fourth-order Runge-Kutta method. An automatic adjustment of the integration step by successive bisections to get the prescribed accuracy is recommended. (4) Check the boundary condition (7). If necessary, correct Co and return to step 3. Otherwise continue. (5) Determine the average effluent concentration and the extent of biodegradation of acetate.
R. N.
1072
PAUNOVIC et al.
For correction of Co values, the half-interval method is used until the interval [C,.,j,,. C,,.,,,J, which satisfies the following conditions is located: For
CO
=
Grnin
GO
atr=flR
(29)
at r = DR.
(30)
For Co = CO.mS.:
Afterwards, the secant method is applied to obtain Co with the prescribed convergence criterion. A real solution to the problem (positive’ concentration values) need not exist throughout a wide region of specification values. A necessary and sufficient condition to obtain a real solution is the existence of an interval CCO.min,CO, mi ] that satisfies conditions (29) and (30). The absence of a real solution for the formulated models indicates that the substrate is completely consumed in a part of biofilm, not reaching its outer border. To describe such a case, one additional variable is to be introduced: the distance r* (R -c r* c BR) up to which the substrate penetrates through the biofilm. At the “moving” moundary r*, the following two conditions are valid: C’S’ _ -
_dcY’ dr
at r = r*_
= 0
This extension of model is not included in the existing computer program.
RESULTS AND DISCUSSION
is obvious that the introduction of parabolic concentration profile in phase 1 (model 2) affects the boundary conditions at the interface between phases 2 and 3. Such a “perturbation” of boundary conditions is a function of the following dimensionless groups: the ratio of diffusion to contact time [eqs (18) and (19)], the Peclet numbers for phases 1 and 2 [eq. (12)], as well as of the geometric parameters Q. The effect of perturbing the boundary conditions on the predicted average effluent concentration depends on the parameters of the dimensionless biofilm model: It
(31) where R ‘=+ y = @/c, Pd3)
= u,ro/D~5)
(314 @lb)
(31c)
(314
To conclude, the effect of neglecting the radial coucentration gradients in the reactor lumen is dependent on the model geometric parameters o! and 8, simplexes tb”/r and KS/C, and four dimensionless groups, Pet) (i = l-3), Da”. To perform a quantitative analysis for design purposes, the influences of some selected most important single design parameters (instead of dimensionless groups) is studied. The included design paratIEters are as follows: inner/outer fiber diameters, biofilm density, X,-, feed concentration, C,, contact time, 7, and biofilm thickness parameter, 8, we qualitatively predict the influences of these single dimensionless parameters in the light of their physical meanings. In view of the significance of the unremoved quantity of polluting substance (acetate), the relative deviation between the average reactor effluent concentrations predicted by models 2 and 1 6 = cp/c;
(32)
is adopted as a measure of error resulting from neglecting the lumen radial concentration gradients. The results of the numerical analysis are shown in Figs 2-5. All the calculations were performed with the following transport and kinetic parameter values (Dall-Bauman et al., 1990): V, = 2.2 cm/d, D(l) = Df2) = 1.09 crn2/d, Dc3) = 0.87 cm3/d, k = 2.5 mg acctatejmg cells/d, K, = 0.01 mg/cm’. Figure 2 shows the relative deviation of the predicted effluent concentrations, 6, as a function of the inner hollow-fiber diameter, &, and biomass density, X,. It is clear that the error resulting from the neglect of the radial lumen concentration gradients becomes significant as the inner fiber diameter increases and this trend is more pronounced for higher biomass density. Tolerating the differences between effluent acetate concentrations predicted by models 1 and 2 up to 20% (6 4 1.2), it appears that one should take into account the radial lumen gradients for inner fiber diameters above 0.03 cm, at biofilm density of 1 mg/cm3. The values of the remaining design parameters are indicated in the figure. At a biofihn density of 100 mg/cm’, the limiting value of the fiber diameter is, however, reduced to about 0.01 cm. Figure 3 illustrates a decrease of d with increase in feed concentration, CJ, and shows that the limiting hollow-fiber diameter is almost insensitive to the C, variation. The influence of biofilm thickness is presented in Fig. 4. The curves /I = 4 and fi = 2 tend to intersect, i.e. the endpoint (boundary of the real model solution) of the curve /I = 4 lies approximately on the /? = 2 curve. Quite a similar trend is exhibited by the curves B = 2 and B = 1.3. As far as the limiting fiber diameter vaues are concerned, the influence of /I is similar to that of XJ. Further, one can see that the curve 0 = 1.3 has a maximum at the fiber diameter of approximately 0.22 cm. This could be explained by the influence of the decreasing extent of acetate biodegradation with the increase of do. Namely, as one can see from
Analysis
1.00
0.00
of a model of hollow-fiber
0.65
I
0.10
bioreactor
wastewater
I
0.15
treatment
1073
I
0.20
0
5
&/cm/
Fig. 2. Relativedeviationbetweenthe averaage effluentconcentrationpredictionsas a functionof the inner hollow-fiber diameter and the biofilm density.
l.00
0.00
.‘“,.‘..,“.‘,‘I’l,‘II’ 0.05
0.10
0.15
0.20
0
‘5
%/cm/
Fig. 3. Relativedeviationbetweentheaverageeffiuentconcentrationpredictionsas a functionof the inner hollow-fiberdiameterand the feed concentration.
Fig. 5, the relative deviation, S, decreases with decrease of the extent of biodegradation. In Figure 5, the extent of acetate biodegradation predicted by the two models is presented as a function of residence time, t. The difference between the extents of conversion predicted by models 1 and 2 appears to decrease with increase in t, while the relative deviation between the predicted effluent concentrations increases. The computing algorithm proved to be reliable and provided convergence in the whole range of the tested model parameter values, in which the models have real solutions. However, one should expect a con-
siderably higher computing time consumption in comparison to a method based on OC or OCFE techniques. CONCLUSIONS
The error resulting from neglection of the radial lumen concentration profile increases with increase of inner hollow-fiber diameter. This trend is more noticeable for higher values of biofilm density and thickness. As a limiting value, an inner fiber diameter of 0.01 cm may be suggested, above which the inclusion of the radial concentration profile in the fiber lumen model is recommendable.
R. N. PAUNOVIC et
1074
I 0.05
o.io
o.i5
al.
o.lo
I
0.25
0.30 dJcm/
i5
Fig. 4. Relative deviation between the average effluent concentration predictions as a function of the inner hollow-fiber diameter and the biofilm thickness: (0) endpoint of the curve @= 4; (A) endpoint of the curve fl = 2.
0.90
0.70 0.50 0.30
Fig. 5. Deviations
between the predictions of models 1 and 2. 6: relative deviation of the average effluent concentrations; X, : extent of acetate biodegradation predicted by model 1;x2: extentof acetatebiodegradation predictedby model 2.
The formulated computing algorithm for bioreactor simulation proved to be reliable, assuring stable convergence in the whole range of the tested parameter values, in which real solutions of the formulated models exist. The algorithm is self-starting, i.e. automatic generation of an initial guess for the reactor effluent concentration is included. NOTATION Al
parameter in parabolic substrate concentration profile in fiber lumen defined by eq. (1b)
A2 a B* b c; C”’
integration parameter in the analytical solution for concentration profile in the fiber wall parameter in quadratic equation (26) integration parameter in the analytical solution for concentration profile in the fiber wall parameter in quadratic equation (26) average substrate effluent concentration predicted by model i, mg/cm3 substrate concentration in phase i, mg/cm’
Analysis
c* Cf
co C D(i)
Dan
do Fo
Gi k
KS L N PC?@)
Qf q/ R r
SO
S
us
46:6-E
of a model of hollow-fiber bioreactor
average substrate concentration in the fiber lumen, mg/cm3 feed substrate concentration, mg@n3 substrate concentration at fiber axis mg/cm3 parameter in quadratic equation (26) diffusion coefficient of the substrate in phase i, crn2/d Damkohler second number, defined by eq. (31d) inner fiber diameter, cm substrate diffusion flux between the fiber lumen and the fiber wall, m&m2 d parameters defined by eqs (21a)-(24b) maximum rate of substrate utilization, mg substrate/mg cells/d Monod half-velocity coefficient, mg/cm” effective length of the fiber, cm number of fibers in the bioreactor Peclet number for phase i defined by eqs (12) and (31~) reactor feed rate, cm”/d hollow-fiber feed rate, cm3/d outer fiber radius, cm radial coordinate inner fiber radius, cm surface area of the interface between the fiber wall and the biofilm, cm’
V V
ur u, Xf
wastewater
treatment
1075
biofilm volume, cm3 radial velocity, cm/d average feed velocity, cm/d wall Permeation velocity, c/d biomass density, mg cells/cm3
Greek letters ratio of inner to outer fiber radius ratio of outer film radius to outer fiber ; radius relative deviation between the reactor aver6 age effluent concentrations predicted by models 2 and 1, defined by eq. (32) T contact time, d @ dimensionless factor defined by eq. (17)
REFERENCES
Carnahan,B., Luther,H. A. and Wilkes, J. O., 1969,AppIied NumericalMethods, PP. 405-415. Wiley, New York. Dali-Bauman,L., Ilias,S. and Govind, R.: 1990,Analysis of hollow fiber bioreactor wastewater treatment. Biotechnol. Bioengng 35, 837-842. Finlayson, B. A., 1988,Nonlinear Analysis in Chemical Engineering, pp. 73-164.McGraw-Hill, New York. Ralston,A. and Wilf, H. S., 1960,Mathematical Methods for Digital Computers, pp. 110-120. Wiley, New York.