- Email: [email protected]
Estimation of empirical coefficients of a granular filtration model
~
Pergamon
Pll:S0043-1354(96)00329-6
War. Res. Vol. 31, No. 5, pp. 1083-1091, 1997 © 1997 Elsevier Science Ltd. All rights reserved Printed in Great Britain 0043-1354/97 $17.00 + 0.00
ESTIMATION OF EMPIRICAL COEFFICIENTS OF A G R A N U L A R FILTRATION MODEL H. B. D H A R M A P P A
I**, H. P R A S A N T H F *, M. U. K R I S H N A 3 a n d Y. X I A O )
~Department of Civil and Mining Engineering, University of Wollongong, Northfields Avenue, Wollongong, NSW 2522, Australia, 2School of Civil Engineering, University of Technology, Sydney, PO 123, Broadway NSW 2007, Australia and 3Environmental Engineering Program, Asian Institute of Technology, GPO Box 2754, Bangkok 10501, Thailand (Received February 1996; accepted in revised form November 1996)
Abstract--Inthis paper a methodology comprising a non-linear least-squares algorithm for determining the empirical parameters of a filtration model is presented and is tested against literature data. Through the case study it is found that incorporation of a multiplication factor for the sum of the square of errors associated with effluent quality could lead to a more balanced estimation of parameters, while achieving a better fit for both effluent quality and headloss. In the case study it is found that the shape of the fitted effluent quality curve has significant impact on the values of the model parameters. The proposed algorithm appears to be well suited for obtaining the best fit for both effluent quality and headloss, simultaneously. The value of the multiplication factor should be fixed based on the accuracy required under both effluent quality and headloss predictions. The model parameters predicted under each of the case studies, correlated well with the experimental observations. © 1997 Elsevier Science Ltd Key words--filtration, modelling, Marquardt algorithm, particle separation
NOMENCLATURE A= co = co = d~ = dp = f= = g = hf =
I= J = K= k~ =
L = Np = Arc = no = n, = S~ = $2 = t=
coefficient in equation (6) influent concentration (mg 1-~) effluent concentration (mg 1-~) average diameter of filter medium grain (m) diameter of particles (m) porosity in the filter medium at any given time porosity of clean bed acceleration due to gravity (m s -2) total headloss across filter bed at any given time (rn) identity matrix hydraulic gradient (m m -~) Koze~y's constant multiplication factor for the sum of the square of errors associated with effluent quality depth of filter medium (m) total number of particles retained in a unit area and depth AL of the filter bed total n umber of filter grains in a unit area and depth AL of the filter bed influent number concentration entering each layer (number ml -~) number of particles in the effluent (number ml -I) shape factor of particles shape factor of filter grain time (:0
*Author to whom all correspondence should be addressed [Fax: (+61) 42 2;1 3238]. The Fortran 77 listing of the algorithm is also available from this author.
v. = filtration rate (m s -)) Y = dependent variable in equation (6) Greek symbols ~t = attachment factor for filter grain to particle ~p = attachment factor for particle to particle fl = fraction of retained particles which act as particle collectors /32 = coefficient of detachment fl' = headloss coefficient AL, Ahf, At = small increments in the respective parameters Ed = porosity of the deposit on the filter grain rl = initial contact efficiency of single filter grain r/p = contact efficiency of a retained particle r/r = single collector efficiency of filter grain 2 = factor for modifying normal equations /~ = dynamic viscosity (kg m ~ s i) p = density o f water (kg m - ' ) Subscripts/symbols L t = time steps Symbols _, t, A, A and * in equations (6)-(8) indicate that the variable is a vector, transposed, estimated, a deviation and obtained in the previous iteration, respectively.
INTRODUCTION F i l t r a t i o n is o n e o f the m a i n s o l i d - l i q u i d s e p a r a t i o n processes widely used in a typical w a t e r t r e a t m e n t system. V a r i o u s theoretical f o r m u l a t i o n s have b e e n p r o p o s e d b y different r e s e a r c h e r s in o r d e r to p r e d i c t its p e r f o r m a n c e . A semi-empirical m o d e l p r o p o s e d by
1083
1084
H. B. Table
Dharmappa et al.
I. Experimental conditions Experimental conditions
Case No.# I) Fulvic acid
2) KCI = 40 mM
3) KCI = 60 mM
4) KCI = 80 mM
5) U‘r= 0.21 x 10-2 m s-1 6) U. = 0.28 x IO-’ m s-’ 7) “a = 0.42 x 10m2m s-’
Filtration of haematite particles with fulvic acid = 2 mg I-‘, KCI = 1 mM, glass beads as medium, d = 0.175 x lo-1 m, L = 3 x IO-* m, f= 0.4. dp = 0.096 x IO-~ m, co = 37-38 mg I-’ as Fez01, pH = 3, cr = 0.043 x lo-* m s-1 Filtration of haematite particles with ionic solution, KCI = 40mM. glass beads as medium, d=0.175 x lOm)m, L=3 x IO-‘m,f=0.4, d,=O.l8 x 10-6m, co=34.>36mg -I as FelOl, pH = 3, Ca= 0.043 x IO-2tn s-1 Filtration of haematite particles with ionic solution, KCI = 60 mM, glass beads as medium, -I as FeJOl, pH = 3, dc = 0.175 x IO-‘m, L = 3 x 10-zm, /= 0.4, dp = 0.46 x 10-6m, co = 34.5-36mg I)*= 0.043 x IO-’ m s-1 Filtration of haematite particles with ionic solution, KCI = 80 mM, glass beads as medium, -I as Fe201, pH = 3, d=O.l75 x IOm’m, L=3 x IO-*m,f=0.4, dp=0.99 x 10-6m, co=34.5-36mg Ua= 0.043 x IO-2 m s-1 Filtration with Kaolin clay suspensions, glass beads as medium, d = 350 x 10-6m (300 x 10m6400 x lO-6 m). f= 0.39, dp = I2 x lO-L m, L = 2 x IO-* m, co = 80 mg I-’ us = 0.21 x IO-‘m s-’ Filtration with Kaolin clay suspensions, glass beads as medium, d, = 350 x 10-Lm (300 x IO-” 400 x 1O-6m), f= 0.39, dp = 12 x lO-6 m, L = 2 x IO-* tn. co = 80 mg I-‘, ~1= 0.28 x IO-> m s-’ Filtration with Kaolin clay suspensions, glass beads as medium, d = 350 x IO-“tn (300 x lO-4400 x 1O-6m), f= 0.39, do = 12 x lO-6 m, L = 2 x IO-‘m, ca = 80 mg I-‘, u, = 0.42 x IO-> m s-’
#First four cases taken from Prasanthi et al. (1994) and the last three from Chang (1989).
O’Melia and Ali (1978) is one of the commonly used models in granular filtration. This model was initially developed for predicting the ripening stage of granular filtration. Later, Vigneswaran and Chang (1986) modified O’Melia and Ali’s model (1978) to extend it to the entire cycle of filtration. The effluent quality across an elemental depth AL can be calculated by:
where, the single collector removal efficiency, v~,, is given by:
Similarly, the headloss can be (Vigneswaran and Chang, 1986):
computed
by
Out of these the value of qa is calculated experimental observations using: 2 rla=
-3
(3)
dcln5 = L 0co ' cl.
(5)
The particle contact efficiency Q is normally assumed to be unity. This leaves three parameters, namely a$, p2 and /I’, to be determined. These parameters are usually determined through a trial and error method (e.g. Darby et al., 1992; Vigneswaran and Chang, 1986). However, Tare and Venkobachar (1985) reported an application of the Marquardt algorithm for determining the above model parameters. But the results of parameter estimation are not discussed in detail by the authors. The objectives of this study include: (i) to develop an algorithm comprising standard numerical techniques for determination of the above (a,/?, j12and fi’) model parameters; (ii) to demonstrate the use of proposed algorithm with experimental data from the literature; and (iii) to interpret the results obtained from the algorithm. ALGORITHM
1 + P’(~,l~~)(~,/~~)‘(s,/sl)* 2 1 + (N,/lv,)(~,/~,)‘(S,/sl)s 1
1 (l-5)
from
FOR PARAMETER
ESTIMATION
Equations (lH4) are principally non-linear. The parameters in these equations can be determined through non-linear regression techniques. One can use three types of technique: (i) Gauss-Newton, (ii) steepest descent and (iii) Marquardt. The Gauss-
where the porosity of the filter bed,f, is evaluated by:
Nc& + N&(1 AL
1’
- td)
Table 2. Summarv of parameter estimation modes
(4)
The meaning of each of the above parameters is given in the Nomenclature. In the above model there are seven empirical parameters to be determined, namely r~,a, q,,, ap, /I, 82 and B’, which are reduced to five by taking the product of two, i.e. qa, t,y,, a,/?, j$ and B’.
Value of k,
Item Mode Mode Mode Mode
I 2a 2b 3
1 10 100 I
Supplementary cases: Mode 2c 3 Mode 2d 20 40 Mode 2e
Error calculation Absolute Absolute Absolute Relative
(from equations (7) and (8)) (from equations (7) and (8)) (from equations (7) and (8)) (from equation (9))
Absolute (from equations (7) and (8)) Absolute (from equations (7) and (8)) Absolute (from eauations (7) and (8))
Estimation of filtration coefficients
1085
Table 3. Sum of the square of errors and model parameters~ Sum o f the square o f errors Case No.
Mode
1) Fulvic acid 2a 2b 2) KCI = 4 0 m M 2a 2b 3) KC1 = 60 m M 2a 2b 3 4) KCI = 8 0 m M 2a 2b 5) v, = 0.21 x 10 2 m .,;-' 2a 2b 6) v, = 0.28 × 10 -2 m .,-~ 2a 2b 7) v, = 0.42 × 10 -2 m s -~ 2a 2b 3 2c 2d 2e
Model parameters
Effluent
Headloss
Total
//2
0.0935 0.0816 0.0588 0.0879 0.0215 0.0213
0.1928 0.2780 73.4684 0.4464 0.0115 0.0192 21.5737 0.0128 1.0035 1.3757 24.5051 1.3798 0.0762 0.1131 49.1670 0.1789 0.4430 2.3069 35.5344 8.5487 0.1030 1.3204 9.3508 0.1194 0.5373 4.1190 57.4511 1.4324 0.5413 10.3087 28.1234
0.2863 0.3596 73.5273 0.5343 0.0330 0.0404 21.5892 0.1209 1.1637 1.4132 24.5369 1.4173 0.1261 0.1621 49.2010 0.3557 0.5413 2.3645 35.5414 8.5821 0.1769 1.3339 9.3598 0.1916 0.6290 4.1749 57.4562 1.5136 0.6323 10.3423 28.1367
0.0169 0.0151 0.0119 0.0140 0.0147 0.0146 0.0128 0.0122 0.0120 0.0084 0.0074 0.0084 0.0052 0.0051 0.0041 0.0037 0.0032 0.11030 0.0020 0.0034 0.0025 0.0056 0.0059 0.0025 0.1 x 10 6 0.0023 0.0017 0.1571 × 10 -4 0.1 × 10 -6 0.0023 0.0019
0.0155 0.1081 0.1602 0.0375 0.0319 0.0376 0.0498 0.0491 0.0341 0.1767 0.0982 0.0576 0.0070 0.0335 0.0739 0.0135
0.0090 0.0722 0.0917 0.0559
0.0051 0.0812 0.0910 0.0336 0.0134
~/~ 0.1 x 10 -s O.l x 10 -s 0.1 x I0 -s 0.1 x 10 -5 0.2706 x 10 -3 0.1755 x 10 -3 0.4561 x 10 -3 0.1321 x 10 -3 0.1 x 10 -6 0.1 x 10 -2 0.5089 × 10 -3 I x 10 -4 0.1 x 10 -6 0.1 x 10 -6 0.1 x 10 6 0.1 x 10 -6 0.0961 0.1114 0.1343 0.1444 0.0860 0.1569 O. 1780 0.0877 0.0183 0.0483 0.0686 0.0206 0.0184 0.0559 0.0617
//' 0.0385 0.036~ 0.2348 0.0361 0.1111 x 10 -2 0.3507 x 10 -z 0.1245 0.1 x 10 -3 0.2213 0.1334 0.3405 0.1333 0.0395 0.0491 0.0417 0.040 0.1 x 10 -~ 0.1 x 10 -6 0.1022 x 10 -6 0.1 x 10 6 0.4408 x 10 -J O.Ol 0.0771 0.1 × 10 -6 0.1 × 10 6 0.1 x 10 -6 0.1 × IO -6 0.1 x 10 6 0.1 × 10 -6 0.1 x 10 6 0.1 x 10 6
aThe m i n i m u m sum o f the square o f errors, under each category and the model parameters selected are shown in bold.
Newton method is based on the linearization of the proposed model and uses the results of linear least squares in a succession of stages. This method, depending on the initial solutions, may converge very slowly, requiring a very large number of iterations. The solution may oscillate very widely, continually reversing direction (Kuester and Mize, 1973; Draper and Smith, 1981). The second metlhod, steepest descent, involves use of the sum of the: square of errors and use of an iterative process to find the minimum of this function. The derivative of this function would give the direction and the location of the next point. Although this method works well in the beginning, when the solution approaches the desired results, the process slows down, especially due to the zigzagging behaviour of the function (Draper and Smith, 1981). The Marquardt (1963) method represents a compromise between the Gauss-Newton and the steepest descent, and appears to combine the best features of both while avoiding their most serious limitations. It was proposed as an extension of the Gauss-Newton method to allow for convergence with relatively poor starting solutions for the unknown coefficients. Due to the fact that the filtration model parameters could vary over a wide range, the Marquardt algorithm is ideally suited for estimating them. As such, it is selected in this study. In this method, the Gauss--Newton normal equations are modified by adding a factor 2:
[A'A + 2._/]A~ = A t (]y1_. )
(6)
where: A is a coefficient; Y a dependent variable; I the identity matrix; 2 a factor for modifying normal equations; and the symbols _, t, ^ , A and • indicate that the variable is a vector, transposed, estimated, a deviation and obtained in the previous iteration, respectively. Here, when 2 approaches ~ , Marquardt's method is identical to the steepest descent. On the other hand, if 2 equals zero, the technique reduces to the Gauss-Newton method. A steepest descent procedure would be expected to converge for poor starting values but requires a lengthy solution time. The Gauss-Newton method, on the other hand, will converge rapidly for good starting estimates. These contrasting properties of the steepest descent and Gauss-Newton methods are used in the Marquardt technique. Initially the Marquardt method approaches steepest descent for identifying good starting solutions and then it imitates Gauss-Newton to achieve rapid convergence. Thus, in the Marquardt procedure, the initial values of 2 are large and will decrease towards zero as a solution is reached. As indicated earlier, in the case of the granular filtration model, there are two main equations (equations (1) and (3)), one for effluent concentration and the other for headloss. It is necessary to fit both equations simultaneously, as each one includes parameters which are evaluated from the other. However, the numerical values for effluent quality
1086
H.B. Dharmappa et al. \ y _ ]z.) = (y, _ ~,.) (for i = n
(co~co) and headloss (hi) are different by about 1-2
orders of magnitude. This causes problems in comparing the sum of the square of errors between effluent quality and headloss. Since the headloss is numerically bigger, it would have greater influence on minimisation of the sum of the square of errors. Thus direct use of the Marquardt algorithm would result in better fitting of the headloss values at the cost of effluent quality. To obviate this problem, the last term ( Y - ~ ' , ) in equation (6) is replaced by equations (7) and (8): Y - ~',) = k~(y, - ~,*) (for i = 1..... n)
(7)
2n)
where 2n indicates the total number of observations (including both effluent quality and headloss values), y, are the observed values for which the first n values correspond to the effluent quality and the last n values correspond to the headloss, and k, are the constant values with which the y, values are multiplied. Equations (7) and (8) imply that the constant factor k, is multiplied only with respect to the effluent quality values but not for the headloss. This is done principally to prop up the effluent quality values in
1.00
o.so-°g' °" a)Effluentqualityp o.,o
~
-
//7
.
o...
,///
-o=,o=i
o
0.40- /
--E--Modo
3
0.30.
0.20.
0.10.
~O
0.001
I 20
I 40
I 60 Timelmin)
I 80
I 100
120
20 b) Heedloss profile
18
E.14 14 .~ .
.
.
.
.5
.
=
= 0 Observed --U-Mode1
,li
_P 18
Mode 2a
~. Mode2b
"1" 6
--m.- Mode 3
0
0
(8)
I
I
I
I
[
20
40
60 Time, mln
80
100
120
Fig. 1. Model fitting for experimental data of Prasanthi et al. (1994) with fulvic acid (Case 1): (a) effluent quality profile; (b) headloss profile.
Estimation offiltrationcoefficients
0.8l~
1087
0.7
0.6 0.5 ~0.4 o
O Observed --I~-- Mode 1 Mode 2a
0.3 0.2
Mode 2b
Mode3
0.1 0
--
t 20
0
I 40
I 60
I 80
I 100
Time, min
I 120
I 140
I 160
180
b) Headloss profile
12
't
_T_7 0
20
40
60
80 100 Time, mln
120
140
160
11 180
Fig. 2. Model fitting for experimental data of Chang (1989) with va = 0.28 x 10-2m s -j (Case 6): (a) effluent quality profile; (b) headloss profile.
order to bring them up to that of headloss. The values of k, can be decided based on the magnitude of difference in numea6cal values of effluent quality and headloss. Direct application of equations (7) and (8) yields an absolute error between the observed and estimated parameters. In other work, Saez and Rittmann (1992) suggested that in a problem such as the one above, the calculation of relative error may give better results than absolute error. In order to calculate the relative error the right-hand side term in equation (6) is once again repl~tced by: Y - ~**') == (Y'-'v'*) (for i = 1 yi / Thus,
in total, three
2n)
(9)
modes of applying the
Marquardt algorithm for the estimation of filtration model parameters are derived: Mode 1, direct application of equation (6); Mode 2, use of equations (7) and (8); and Mode 3, use of equation (9). CASE STUDY All three modes of applying the Marquardt algorithm explained above are demonstrated using the experimental results of Prasanthi et al. (1994) and Chang (1989), and the predictions obtained are compared. In total seven experimental runs are considered, which have been carried out under different experimental conditions, as summarised in Table 1. The parameters are estimated under different predicting modes, as listed in Table 2. Mode 1 is the normal Marquardt algorithm with no modifications.
H. B. Dharmappa et al.
1088
In M o d e 2a the effluent y, values are multiplied by ki = 10, and Mode 2b includes multiplication with ki = 100. Finally, in Mode 3, the relative sum of the square of errors is calculated instead of absolute values. In addition, some supplementary predicting modes are selected with different k, values, as shown in Table 2.
RESULTS AND DISCUSSION Sum o f square o f errors The results from the estimation of model parameters are summarised in Table 3. Observing the sum of the square of errors for the effluent quality, it can be seen that values are a minimum under Mode
0.8
<>
<>
0.6.
0.5. (.1 0.4 (>Observed --O--Mode 1 --A--Mode 28 Mode 2b = Mode3
0.3 0.2 0.1 0
i
0
20
,i 40
i
60
I 8O
I
I
I
100
120
140
Time, min
160
18-r b) Headloss profile 16
14 I 12 E h 0
~
@ Obsocvod
--13--- Mode I -'-A-- Mode 2a ~-. Mode2b -..Xl-- IVlode3
8
0
I
I
I
20
40
60
I
80 Time, rain
I
I
I
100
120
140
1(10
Fig. 3. Model fitting for experimental data of Chang (1989) with v. = 0.42 x 10-:cm s-' (Case 7): (a) effluent quality profile; (b) headloss profile.
Estimation of filtration coefficients 0.9
1089
m
0.8
0:7
0
0.6
0
O
0.5
~
Observed
--.O-- Mode 2c •- - D - - Mode 2a
°0.4 0.3
X
Mode 2d Mode 2e
)K
Mode 2b
0.2
0.1
0
I 20
I 40
I
I
I
80
Tim 8°e.min
1oo
I 120
I 140
160
b) Headloss profile 16
14
12
E
0 u~l 0
.9o
O Observed --'0'--Mode 2c Mode 2a Mode 2d X Mode 2e I Mode 2b
"lO 8 "1"
6
OI 0
I 20
I 40
I 60
I 80
I 1O0
I 120
I 140
160
Time, rain Fig. 4. Model fitting for experimental data of Chang (1989) with v~ = 0.42 × 10-:m s-' (Case 7) for supplementary predicting modes: (a) effluent quality profiles; (b) headloss profile.
2b. Then they gradually increase under Modes 2a and l, with the sum of square of errors being higher under Mode I. This is due to the obvious reason that the significance of the sum of the square of errors for effluent quality increases with the ki values. Hence, the better fitting for effluent quality under Mode 2b. However, Mode 3 is unable to produce any consistent results. On most occasions the sum of the square of
errors obtained for effluent quality under Mode 3 is higher than values obtained from other modes of prediction. On the other hand, the sum of the square of errors for headloss and combined predictions are a minimum under Mode 1. They slightly increase under Mode 2a, while increasing drastically under Mode 2b. Here, again, there is no definite trend for the sum of
1090
H. B. Dharmappa et al.
the square of errors in Mode 3, where in most cases values are found to be higher than for Modes ! and 2a. Considering now the total sum of the square of errors obtained for combined predictions (Table 1), the values are a minimum under Mode ! and, as in the case of headloss, they increase under Modes 2a and 2b. Values of the total sum of the square of errors under Mode 2b are highest when predicted headloss values are high. Here, once again, Mode 3 did not show any specific trend; however, in most cases it gave higher values for the sum of the square of errors compared to those obtained under Modes l and 2a. The above results indicate that, while Mode 2b gives consistently better fitting for effluent quality, Mode 1 yields good fitting for headloss. In between these two extremes, Mode 2a appears to be a good compromise. It satisfies, to some extent, both the effluent and headloss predictions.
Mode 1. When the k, values are increased, the shape of the curve is changed to conform to the observed results. This change in shape of the effluent quality curve has significant impact on the model parameters as shown in Table 3, especially on ctpfl and f12. The values of f12 changed by about four orders of magnitude, while ~pfl values altered by a factor of three. As expected, as ki values are increased a better fit is given for the effluent quality data, while the headloss predictions exhibit an increasingly bad fit. Thus, it is necessary to select an appropriate mode of prediction which satisfies both effluent quality and headloss fittings. As shown in Fig. 4, Mode 2c, which exhibits a close fit with the headloss observations, does not provide a good fit for the effluent quality. On the other hand, Mode 2b provides the best fit for effluent quality observations while giving a bad fit for headloss. Thus a compromise solution could be the predictions under Modes 2a or 2d.
Nature o f fit
Model parameters
Selected predictions are presented in Figs 1-3. Figure l shows the predictions for experimental data obtained by Prasanthi et al. (1994) with fulvic acid. The data presented in the figure for the effluent quality give an S-shaped curve. The predictions under Modes l and 3 do not reflect the shape of the curve, whereas those from Mode 2b do, to some extent, conform to the shape with a minimum sum of the square of errors. However, the headloss is not predicted very well under Mode 2b. On the other hand, Mode 2a predicts headloss very well, while predicting effluent quality less closely than Mode 2b. However, Mode 2a makes an attempt to predict the S shape of the effluent quality profile. As such, Mode 2a can be considered to be a compromise between Modes 1 and 2b. Figure 2 shows the results obtained for Chang's (1989) data with velocity, va = l0 m h ~. The form of the observed data is closely followed by Mode 2a and 2b predictions. However, the predictions under Modes 1 and 3 give an altogether different shape for the effluent quality curve. This shifting of the effluent curve has a significant impact on the values of model parameters, as shown in Table 3. As in Fig. l, the headloss predictions in Fig. 2 under Mode 2b are significantly different from the observed values. Here, once again, the Mode 2a predictions appear to be a good compromise between Modes 1 and 2b. The above conclusions are further corroborated in Fig. 3. However, in this case, the difference between the predictions under Modes 2a and 2b is substantial. In order to obtain intermediate predictions, additional runs are carried out in which the k~ values are changed. Details of the supplementary mode of predictions are given in Table 2. The predicted results are presented in Fig. 4. It can be seen from this figure than the values of k~ determine the shape of the effluent quality curve. In Mode 2c, with k~ = 3, the prediction of effluent quality is similar to that of
The model parameters predicted are shown in Table 3 (model parameters obtained under Mode 2a are shown in bold and italics), and correlate well with the experimental conditions. The values of ~tpfl, the attachment coefficient, for Chang's (1989) data are found to be higher, reflecting the higher efficiency of particle removal as observed in the experimental study (Figs 2 and 3). On the other hand, the values of f12, the detachment coefficient, are found to be generally higher for the data of Prasanthi et al. (1994), reflecting the immediate breakthrough observed in their experimental studies (Fig. 1). The value of fl', the headloss coefficient (which is simply the fraction of particles deposited in the filter bed contributing to the additional surface area for headloss development), is found to be generally lower for Chang's (1989) data, despite showing a significant increase in headloss with respect to time. This may be attributed to greater removal of particles in Chang's experiments, whereby the actual fraction of particles contributing towards additional surface area for headloss development is reduced. Also, the variation in fl' can be attributed to the size of the particles being filtered. As reported in several earlier studies (e.g. Dharmappa et al., 1992), the smaller the particles are, the higher the headloss is, as they contribute more surface area towards headloss development. In Cases 1~,, as the size of the particles is very small (Table 1), the value offl' can be expected to be high. CONCLUSIONS In this paper a methodology comprising the Marquardt (1964) algorithm for determining the empirical parameters of the modified O'Melia and Ali model (Vigneswaran and Chang, 1986) is presented and is tested against literature data. Through the case study it is found that a multiplication factor for the
Estimation of filtration coefficients sum of the square of errors belonging to effluent quality could resulll in a more balanced estimation of parameters, while achieving a better fit for both effluent quality and headloss, simultaneously. In the case study it is found that the shape of the fitted effluent quality curve has significant impact on the values of the model parameters. If only the total sum of the square of errors is considered, there are several possible ways of fitting the observed data with respect to a filtration process. It is necessary to check the shape of the fitted curve as well. It should conform to the shape of the observed data. Towards this end, the proposed algorithm appears to be well suited to give the best curve fitted to a set of both observed effluent quality and headloss data. The value of the multiplication factor should be fixed based on the accuracy required under both effluent quality and headloss predictions, and the difference in the observed absolute values of effluent quality and headloss. The model parameters predicted reflected, to some extent, the experimental conditions adopted by Chang (1989) and ]?rasanthi et al. (1994). The values of ctpfl, the attachment coefficient, for Chang's data (1989) are found to be higher, reflecting the better efficiency of pardcle removal observed during experimental studies. On the other hand, the values of f12, the detachment coefficient, are found to be generally higher for the results of Prasanthi et al. (1994), reflecting the immediate breakthrough observed in their experimental work. The values of fl', the headloss coefficient, for Chang's (1989) exper-
1091
iments are found to be lower as there are more particles being removed and hence the actual fraction of particles contributing to the headloss development is minimum. The higher values of fl' obtained for the results of Prasanthi et al. (1994) can also be partly attributed to the smaller size of the particles filtered.
REFERENCES
Chang J. S. (t989) Mathematical modeling of transient behaviour of deep bed filtration. D.Eng. dissertation, Asian Institute of Technology, Bangkok. Darby J. L., Attanasio R. E. and Lawler D. F. (1992) Filtration of heterodisperse suspensions: Modeling of particle removal and headloss. Wat. Res. 26, 711-726. Dharmappa H. B., Verink J., Fujiwara O. and Vigneswaran S. (1992) Optimisation of granular bed filtration treating polydispersed suspension. Wat. Res. 26, 130%1318. Draper N. R. and Smith H. (1981) Applied Regression Analysis. 2nd edn. John Wiley & Sons, New York. Kuester J. L. and Mize J. H. (1973) Optimisation Techniques. McGraw-Hill Book Co., New York. Marquardt D. M. (1963) An algorithm for least-squares estimation of nonlinear parameters. J. Soc. Industrial Appl. Mathematics 11, 431-441. O'Melia C. R. and Ali W. (1978) The role of retained particles in deep bed filtration. Prog. Wat. Technol. 10, 16%182. Prasanthi H., Vigneswaran S., Waite T. D. and Ben Aim R. (1994) Filtration of submicron particles: effect of ionic strength and organic substances. Wat. Sci. Technol. 30, 149-158. Saez P. B. and Rittmann B. E. (1992) Model-parameter estimation using least squares. Wat. Res. 26, 711-726. Vigneswaran S. and Chang J. S. (1986) Mathematical modeling of the entire cycle of deep bed filtration. Wat. Air Soil Pollut. 29, 155-164.
Pergamon
Pll:S0043-1354(96)00329-6
War. Res. Vol. 31, No. 5, pp. 1083-1091, 1997 © 1997 Elsevier Science Ltd. All rights reserved Printed in Great Britain 0043-1354/97 $17.00 + 0.00
ESTIMATION OF EMPIRICAL COEFFICIENTS OF A G R A N U L A R FILTRATION MODEL H. B. D H A R M A P P A
I**, H. P R A S A N T H F *, M. U. K R I S H N A 3 a n d Y. X I A O )
~Department of Civil and Mining Engineering, University of Wollongong, Northfields Avenue, Wollongong, NSW 2522, Australia, 2School of Civil Engineering, University of Technology, Sydney, PO 123, Broadway NSW 2007, Australia and 3Environmental Engineering Program, Asian Institute of Technology, GPO Box 2754, Bangkok 10501, Thailand (Received February 1996; accepted in revised form November 1996)
Abstract--Inthis paper a methodology comprising a non-linear least-squares algorithm for determining the empirical parameters of a filtration model is presented and is tested against literature data. Through the case study it is found that incorporation of a multiplication factor for the sum of the square of errors associated with effluent quality could lead to a more balanced estimation of parameters, while achieving a better fit for both effluent quality and headloss. In the case study it is found that the shape of the fitted effluent quality curve has significant impact on the values of the model parameters. The proposed algorithm appears to be well suited for obtaining the best fit for both effluent quality and headloss, simultaneously. The value of the multiplication factor should be fixed based on the accuracy required under both effluent quality and headloss predictions. The model parameters predicted under each of the case studies, correlated well with the experimental observations. © 1997 Elsevier Science Ltd Key words--filtration, modelling, Marquardt algorithm, particle separation
NOMENCLATURE A= co = co = d~ = dp = f= = g = hf =
I= J = K= k~ =
L = Np = Arc = no = n, = S~ = $2 = t=
coefficient in equation (6) influent concentration (mg 1-~) effluent concentration (mg 1-~) average diameter of filter medium grain (m) diameter of particles (m) porosity in the filter medium at any given time porosity of clean bed acceleration due to gravity (m s -2) total headloss across filter bed at any given time (rn) identity matrix hydraulic gradient (m m -~) Koze~y's constant multiplication factor for the sum of the square of errors associated with effluent quality depth of filter medium (m) total number of particles retained in a unit area and depth AL of the filter bed total n umber of filter grains in a unit area and depth AL of the filter bed influent number concentration entering each layer (number ml -~) number of particles in the effluent (number ml -I) shape factor of particles shape factor of filter grain time (:0
*Author to whom all correspondence should be addressed [Fax: (+61) 42 2;1 3238]. The Fortran 77 listing of the algorithm is also available from this author.
v. = filtration rate (m s -)) Y = dependent variable in equation (6) Greek symbols ~t = attachment factor for filter grain to particle ~p = attachment factor for particle to particle fl = fraction of retained particles which act as particle collectors /32 = coefficient of detachment fl' = headloss coefficient AL, Ahf, At = small increments in the respective parameters Ed = porosity of the deposit on the filter grain rl = initial contact efficiency of single filter grain r/p = contact efficiency of a retained particle r/r = single collector efficiency of filter grain 2 = factor for modifying normal equations /~ = dynamic viscosity (kg m ~ s i) p = density o f water (kg m - ' ) Subscripts/symbols L t = time steps Symbols _, t, A, A and * in equations (6)-(8) indicate that the variable is a vector, transposed, estimated, a deviation and obtained in the previous iteration, respectively.
INTRODUCTION F i l t r a t i o n is o n e o f the m a i n s o l i d - l i q u i d s e p a r a t i o n processes widely used in a typical w a t e r t r e a t m e n t system. V a r i o u s theoretical f o r m u l a t i o n s have b e e n p r o p o s e d b y different r e s e a r c h e r s in o r d e r to p r e d i c t its p e r f o r m a n c e . A semi-empirical m o d e l p r o p o s e d by
1083
1084
H. B. Table
Dharmappa et al.
I. Experimental conditions Experimental conditions
Case No.# I) Fulvic acid
2) KCI = 40 mM
3) KCI = 60 mM
4) KCI = 80 mM
5) U‘r= 0.21 x 10-2 m s-1 6) U. = 0.28 x IO-’ m s-’ 7) “a = 0.42 x 10m2m s-’
Filtration of haematite particles with fulvic acid = 2 mg I-‘, KCI = 1 mM, glass beads as medium, d = 0.175 x lo-1 m, L = 3 x IO-* m, f= 0.4. dp = 0.096 x IO-~ m, co = 37-38 mg I-’ as Fez01, pH = 3, cr = 0.043 x lo-* m s-1 Filtration of haematite particles with ionic solution, KCI = 40mM. glass beads as medium, d=0.175 x lOm)m, L=3 x IO-‘m,f=0.4, d,=O.l8 x 10-6m, co=34.>36mg -I as FelOl, pH = 3, Ca= 0.043 x IO-2tn s-1 Filtration of haematite particles with ionic solution, KCI = 60 mM, glass beads as medium, -I as FeJOl, pH = 3, dc = 0.175 x IO-‘m, L = 3 x 10-zm, /= 0.4, dp = 0.46 x 10-6m, co = 34.5-36mg I)*= 0.043 x IO-’ m s-1 Filtration of haematite particles with ionic solution, KCI = 80 mM, glass beads as medium, -I as Fe201, pH = 3, d=O.l75 x IOm’m, L=3 x IO-*m,f=0.4, dp=0.99 x 10-6m, co=34.5-36mg Ua= 0.043 x IO-2 m s-1 Filtration with Kaolin clay suspensions, glass beads as medium, d = 350 x 10-6m (300 x 10m6400 x lO-6 m). f= 0.39, dp = I2 x lO-L m, L = 2 x IO-* m, co = 80 mg I-’ us = 0.21 x IO-‘m s-’ Filtration with Kaolin clay suspensions, glass beads as medium, d, = 350 x 10-Lm (300 x IO-” 400 x 1O-6m), f= 0.39, dp = 12 x lO-6 m, L = 2 x IO-* tn. co = 80 mg I-‘, ~1= 0.28 x IO-> m s-’ Filtration with Kaolin clay suspensions, glass beads as medium, d = 350 x IO-“tn (300 x lO-4400 x 1O-6m), f= 0.39, do = 12 x lO-6 m, L = 2 x IO-‘m, ca = 80 mg I-‘, u, = 0.42 x IO-> m s-’
#First four cases taken from Prasanthi et al. (1994) and the last three from Chang (1989).
O’Melia and Ali (1978) is one of the commonly used models in granular filtration. This model was initially developed for predicting the ripening stage of granular filtration. Later, Vigneswaran and Chang (1986) modified O’Melia and Ali’s model (1978) to extend it to the entire cycle of filtration. The effluent quality across an elemental depth AL can be calculated by:
where, the single collector removal efficiency, v~,, is given by:
Similarly, the headloss can be (Vigneswaran and Chang, 1986):
computed
by
Out of these the value of qa is calculated experimental observations using: 2 rla=
-3
(3)
dcln5 = L 0co ' cl.
(5)
The particle contact efficiency Q is normally assumed to be unity. This leaves three parameters, namely a$, p2 and /I’, to be determined. These parameters are usually determined through a trial and error method (e.g. Darby et al., 1992; Vigneswaran and Chang, 1986). However, Tare and Venkobachar (1985) reported an application of the Marquardt algorithm for determining the above model parameters. But the results of parameter estimation are not discussed in detail by the authors. The objectives of this study include: (i) to develop an algorithm comprising standard numerical techniques for determination of the above (a,/?, j12and fi’) model parameters; (ii) to demonstrate the use of proposed algorithm with experimental data from the literature; and (iii) to interpret the results obtained from the algorithm. ALGORITHM
1 + P’(~,l~~)(~,/~~)‘(s,/sl)* 2 1 + (N,/lv,)(~,/~,)‘(S,/sl)s 1
1 (l-5)
from
FOR PARAMETER
ESTIMATION
Equations (lH4) are principally non-linear. The parameters in these equations can be determined through non-linear regression techniques. One can use three types of technique: (i) Gauss-Newton, (ii) steepest descent and (iii) Marquardt. The Gauss-
where the porosity of the filter bed,f, is evaluated by:
Nc& + N&(1 AL
1’
- td)
Table 2. Summarv of parameter estimation modes
(4)
The meaning of each of the above parameters is given in the Nomenclature. In the above model there are seven empirical parameters to be determined, namely r~,a, q,,, ap, /I, 82 and B’, which are reduced to five by taking the product of two, i.e. qa, t,y,, a,/?, j$ and B’.
Value of k,
Item Mode Mode Mode Mode
I 2a 2b 3
1 10 100 I
Supplementary cases: Mode 2c 3 Mode 2d 20 40 Mode 2e
Error calculation Absolute Absolute Absolute Relative
(from equations (7) and (8)) (from equations (7) and (8)) (from equations (7) and (8)) (from equation (9))
Absolute (from equations (7) and (8)) Absolute (from equations (7) and (8)) Absolute (from eauations (7) and (8))
Estimation of filtration coefficients
1085
Table 3. Sum of the square of errors and model parameters~ Sum o f the square o f errors Case No.
Mode
1) Fulvic acid 2a 2b 2) KCI = 4 0 m M 2a 2b 3) KC1 = 60 m M 2a 2b 3 4) KCI = 8 0 m M 2a 2b 5) v, = 0.21 x 10 2 m .,;-' 2a 2b 6) v, = 0.28 × 10 -2 m .,-~ 2a 2b 7) v, = 0.42 × 10 -2 m s -~ 2a 2b 3 2c 2d 2e
Model parameters
Effluent
Headloss
Total
//2
0.0935 0.0816 0.0588 0.0879 0.0215 0.0213
0.1928 0.2780 73.4684 0.4464 0.0115 0.0192 21.5737 0.0128 1.0035 1.3757 24.5051 1.3798 0.0762 0.1131 49.1670 0.1789 0.4430 2.3069 35.5344 8.5487 0.1030 1.3204 9.3508 0.1194 0.5373 4.1190 57.4511 1.4324 0.5413 10.3087 28.1234
0.2863 0.3596 73.5273 0.5343 0.0330 0.0404 21.5892 0.1209 1.1637 1.4132 24.5369 1.4173 0.1261 0.1621 49.2010 0.3557 0.5413 2.3645 35.5414 8.5821 0.1769 1.3339 9.3598 0.1916 0.6290 4.1749 57.4562 1.5136 0.6323 10.3423 28.1367
0.0169 0.0151 0.0119 0.0140 0.0147 0.0146 0.0128 0.0122 0.0120 0.0084 0.0074 0.0084 0.0052 0.0051 0.0041 0.0037 0.0032 0.11030 0.0020 0.0034 0.0025 0.0056 0.0059 0.0025 0.1 x 10 6 0.0023 0.0017 0.1571 × 10 -4 0.1 × 10 -6 0.0023 0.0019
0.0155 0.1081 0.1602 0.0375 0.0319 0.0376 0.0498 0.0491 0.0341 0.1767 0.0982 0.0576 0.0070 0.0335 0.0739 0.0135
0.0090 0.0722 0.0917 0.0559
0.0051 0.0812 0.0910 0.0336 0.0134
~/~ 0.1 x 10 -s O.l x 10 -s 0.1 x I0 -s 0.1 x 10 -5 0.2706 x 10 -3 0.1755 x 10 -3 0.4561 x 10 -3 0.1321 x 10 -3 0.1 x 10 -6 0.1 x 10 -2 0.5089 × 10 -3 I x 10 -4 0.1 x 10 -6 0.1 x 10 -6 0.1 x 10 6 0.1 x 10 -6 0.0961 0.1114 0.1343 0.1444 0.0860 0.1569 O. 1780 0.0877 0.0183 0.0483 0.0686 0.0206 0.0184 0.0559 0.0617
//' 0.0385 0.036~ 0.2348 0.0361 0.1111 x 10 -2 0.3507 x 10 -z 0.1245 0.1 x 10 -3 0.2213 0.1334 0.3405 0.1333 0.0395 0.0491 0.0417 0.040 0.1 x 10 -~ 0.1 x 10 -6 0.1022 x 10 -6 0.1 x 10 6 0.4408 x 10 -J O.Ol 0.0771 0.1 × 10 -6 0.1 × 10 6 0.1 x 10 -6 0.1 × IO -6 0.1 x 10 6 0.1 × 10 -6 0.1 x 10 6 0.1 x 10 6
aThe m i n i m u m sum o f the square o f errors, under each category and the model parameters selected are shown in bold.
Newton method is based on the linearization of the proposed model and uses the results of linear least squares in a succession of stages. This method, depending on the initial solutions, may converge very slowly, requiring a very large number of iterations. The solution may oscillate very widely, continually reversing direction (Kuester and Mize, 1973; Draper and Smith, 1981). The second metlhod, steepest descent, involves use of the sum of the: square of errors and use of an iterative process to find the minimum of this function. The derivative of this function would give the direction and the location of the next point. Although this method works well in the beginning, when the solution approaches the desired results, the process slows down, especially due to the zigzagging behaviour of the function (Draper and Smith, 1981). The Marquardt (1963) method represents a compromise between the Gauss-Newton and the steepest descent, and appears to combine the best features of both while avoiding their most serious limitations. It was proposed as an extension of the Gauss-Newton method to allow for convergence with relatively poor starting solutions for the unknown coefficients. Due to the fact that the filtration model parameters could vary over a wide range, the Marquardt algorithm is ideally suited for estimating them. As such, it is selected in this study. In this method, the Gauss--Newton normal equations are modified by adding a factor 2:
[A'A + 2._/]A~ = A t (]y1_. )
(6)
where: A is a coefficient; Y a dependent variable; I the identity matrix; 2 a factor for modifying normal equations; and the symbols _, t, ^ , A and • indicate that the variable is a vector, transposed, estimated, a deviation and obtained in the previous iteration, respectively. Here, when 2 approaches ~ , Marquardt's method is identical to the steepest descent. On the other hand, if 2 equals zero, the technique reduces to the Gauss-Newton method. A steepest descent procedure would be expected to converge for poor starting values but requires a lengthy solution time. The Gauss-Newton method, on the other hand, will converge rapidly for good starting estimates. These contrasting properties of the steepest descent and Gauss-Newton methods are used in the Marquardt technique. Initially the Marquardt method approaches steepest descent for identifying good starting solutions and then it imitates Gauss-Newton to achieve rapid convergence. Thus, in the Marquardt procedure, the initial values of 2 are large and will decrease towards zero as a solution is reached. As indicated earlier, in the case of the granular filtration model, there are two main equations (equations (1) and (3)), one for effluent concentration and the other for headloss. It is necessary to fit both equations simultaneously, as each one includes parameters which are evaluated from the other. However, the numerical values for effluent quality
1086
H.B. Dharmappa et al. \ y _ ]z.) = (y, _ ~,.) (for i = n
(co~co) and headloss (hi) are different by about 1-2
orders of magnitude. This causes problems in comparing the sum of the square of errors between effluent quality and headloss. Since the headloss is numerically bigger, it would have greater influence on minimisation of the sum of the square of errors. Thus direct use of the Marquardt algorithm would result in better fitting of the headloss values at the cost of effluent quality. To obviate this problem, the last term ( Y - ~ ' , ) in equation (6) is replaced by equations (7) and (8): Y - ~',) = k~(y, - ~,*) (for i = 1..... n)
(7)
2n)
where 2n indicates the total number of observations (including both effluent quality and headloss values), y, are the observed values for which the first n values correspond to the effluent quality and the last n values correspond to the headloss, and k, are the constant values with which the y, values are multiplied. Equations (7) and (8) imply that the constant factor k, is multiplied only with respect to the effluent quality values but not for the headloss. This is done principally to prop up the effluent quality values in
1.00
o.so-°g' °" a)Effluentqualityp o.,o
~
-
//7
.
o...
,///
-o=,o=i
o
0.40- /
--E--Modo
3
0.30.
0.20.
0.10.
~O
0.001
I 20
I 40
I 60 Timelmin)
I 80
I 100
120
20 b) Heedloss profile
18
E.14 14 .~ .
.
.
.
.5
.
=
= 0 Observed --U-Mode1
,li
_P 18
Mode 2a
~. Mode2b
"1" 6
--m.- Mode 3
0
0
(8)
I
I
I
I
[
20
40
60 Time, mln
80
100
120
Fig. 1. Model fitting for experimental data of Prasanthi et al. (1994) with fulvic acid (Case 1): (a) effluent quality profile; (b) headloss profile.
Estimation offiltrationcoefficients
0.8l~
1087
0.7
0.6 0.5 ~0.4 o
O Observed --I~-- Mode 1 Mode 2a
0.3 0.2
Mode 2b
Mode3
0.1 0
--
t 20
0
I 40
I 60
I 80
I 100
Time, min
I 120
I 140
I 160
180
b) Headloss profile
12
't
_T_7 0
20
40
60
80 100 Time, mln
120
140
160
11 180
Fig. 2. Model fitting for experimental data of Chang (1989) with va = 0.28 x 10-2m s -j (Case 6): (a) effluent quality profile; (b) headloss profile.
order to bring them up to that of headloss. The values of k, can be decided based on the magnitude of difference in numea6cal values of effluent quality and headloss. Direct application of equations (7) and (8) yields an absolute error between the observed and estimated parameters. In other work, Saez and Rittmann (1992) suggested that in a problem such as the one above, the calculation of relative error may give better results than absolute error. In order to calculate the relative error the right-hand side term in equation (6) is once again repl~tced by: Y - ~**') == (Y'-'v'*) (for i = 1 yi / Thus,
in total, three
2n)
(9)
modes of applying the
Marquardt algorithm for the estimation of filtration model parameters are derived: Mode 1, direct application of equation (6); Mode 2, use of equations (7) and (8); and Mode 3, use of equation (9). CASE STUDY All three modes of applying the Marquardt algorithm explained above are demonstrated using the experimental results of Prasanthi et al. (1994) and Chang (1989), and the predictions obtained are compared. In total seven experimental runs are considered, which have been carried out under different experimental conditions, as summarised in Table 1. The parameters are estimated under different predicting modes, as listed in Table 2. Mode 1 is the normal Marquardt algorithm with no modifications.
H. B. Dharmappa et al.
1088
In M o d e 2a the effluent y, values are multiplied by ki = 10, and Mode 2b includes multiplication with ki = 100. Finally, in Mode 3, the relative sum of the square of errors is calculated instead of absolute values. In addition, some supplementary predicting modes are selected with different k, values, as shown in Table 2.
RESULTS AND DISCUSSION Sum o f square o f errors The results from the estimation of model parameters are summarised in Table 3. Observing the sum of the square of errors for the effluent quality, it can be seen that values are a minimum under Mode
0.8
<>
<>
0.6.
0.5. (.1 0.4 (>Observed --O--Mode 1 --A--Mode 28 Mode 2b = Mode3
0.3 0.2 0.1 0
i
0
20
,i 40
i
60
I 8O
I
I
I
100
120
140
Time, min
160
18-r b) Headloss profile 16
14 I 12 E h 0
~
@ Obsocvod
--13--- Mode I -'-A-- Mode 2a ~-. Mode2b -..Xl-- IVlode3
8
0
I
I
I
20
40
60
I
80 Time, rain
I
I
I
100
120
140
1(10
Fig. 3. Model fitting for experimental data of Chang (1989) with v. = 0.42 x 10-:cm s-' (Case 7): (a) effluent quality profile; (b) headloss profile.
Estimation of filtration coefficients 0.9
1089
m
0.8
0:7
0
0.6
0
O
0.5
~
Observed
--.O-- Mode 2c •- - D - - Mode 2a
°0.4 0.3
X
Mode 2d Mode 2e
)K
Mode 2b
0.2
0.1
0
I 20
I 40
I
I
I
80
Tim 8°e.min
1oo
I 120
I 140
160
b) Headloss profile 16
14
12
E
0 u~l 0
.9o
O Observed --'0'--Mode 2c Mode 2a Mode 2d X Mode 2e I Mode 2b
"lO 8 "1"
6
OI 0
I 20
I 40
I 60
I 80
I 1O0
I 120
I 140
160
Time, rain Fig. 4. Model fitting for experimental data of Chang (1989) with v~ = 0.42 × 10-:m s-' (Case 7) for supplementary predicting modes: (a) effluent quality profiles; (b) headloss profile.
2b. Then they gradually increase under Modes 2a and l, with the sum of square of errors being higher under Mode I. This is due to the obvious reason that the significance of the sum of the square of errors for effluent quality increases with the ki values. Hence, the better fitting for effluent quality under Mode 2b. However, Mode 3 is unable to produce any consistent results. On most occasions the sum of the square of
errors obtained for effluent quality under Mode 3 is higher than values obtained from other modes of prediction. On the other hand, the sum of the square of errors for headloss and combined predictions are a minimum under Mode 1. They slightly increase under Mode 2a, while increasing drastically under Mode 2b. Here, again, there is no definite trend for the sum of
1090
H. B. Dharmappa et al.
the square of errors in Mode 3, where in most cases values are found to be higher than for Modes ! and 2a. Considering now the total sum of the square of errors obtained for combined predictions (Table 1), the values are a minimum under Mode ! and, as in the case of headloss, they increase under Modes 2a and 2b. Values of the total sum of the square of errors under Mode 2b are highest when predicted headloss values are high. Here, once again, Mode 3 did not show any specific trend; however, in most cases it gave higher values for the sum of the square of errors compared to those obtained under Modes l and 2a. The above results indicate that, while Mode 2b gives consistently better fitting for effluent quality, Mode 1 yields good fitting for headloss. In between these two extremes, Mode 2a appears to be a good compromise. It satisfies, to some extent, both the effluent and headloss predictions.
Mode 1. When the k, values are increased, the shape of the curve is changed to conform to the observed results. This change in shape of the effluent quality curve has significant impact on the model parameters as shown in Table 3, especially on ctpfl and f12. The values of f12 changed by about four orders of magnitude, while ~pfl values altered by a factor of three. As expected, as ki values are increased a better fit is given for the effluent quality data, while the headloss predictions exhibit an increasingly bad fit. Thus, it is necessary to select an appropriate mode of prediction which satisfies both effluent quality and headloss fittings. As shown in Fig. 4, Mode 2c, which exhibits a close fit with the headloss observations, does not provide a good fit for the effluent quality. On the other hand, Mode 2b provides the best fit for effluent quality observations while giving a bad fit for headloss. Thus a compromise solution could be the predictions under Modes 2a or 2d.
Nature o f fit
Model parameters
Selected predictions are presented in Figs 1-3. Figure l shows the predictions for experimental data obtained by Prasanthi et al. (1994) with fulvic acid. The data presented in the figure for the effluent quality give an S-shaped curve. The predictions under Modes l and 3 do not reflect the shape of the curve, whereas those from Mode 2b do, to some extent, conform to the shape with a minimum sum of the square of errors. However, the headloss is not predicted very well under Mode 2b. On the other hand, Mode 2a predicts headloss very well, while predicting effluent quality less closely than Mode 2b. However, Mode 2a makes an attempt to predict the S shape of the effluent quality profile. As such, Mode 2a can be considered to be a compromise between Modes 1 and 2b. Figure 2 shows the results obtained for Chang's (1989) data with velocity, va = l0 m h ~. The form of the observed data is closely followed by Mode 2a and 2b predictions. However, the predictions under Modes 1 and 3 give an altogether different shape for the effluent quality curve. This shifting of the effluent curve has a significant impact on the values of model parameters, as shown in Table 3. As in Fig. l, the headloss predictions in Fig. 2 under Mode 2b are significantly different from the observed values. Here, once again, the Mode 2a predictions appear to be a good compromise between Modes 1 and 2b. The above conclusions are further corroborated in Fig. 3. However, in this case, the difference between the predictions under Modes 2a and 2b is substantial. In order to obtain intermediate predictions, additional runs are carried out in which the k~ values are changed. Details of the supplementary mode of predictions are given in Table 2. The predicted results are presented in Fig. 4. It can be seen from this figure than the values of k~ determine the shape of the effluent quality curve. In Mode 2c, with k~ = 3, the prediction of effluent quality is similar to that of
The model parameters predicted are shown in Table 3 (model parameters obtained under Mode 2a are shown in bold and italics), and correlate well with the experimental conditions. The values of ~tpfl, the attachment coefficient, for Chang's (1989) data are found to be higher, reflecting the higher efficiency of particle removal as observed in the experimental study (Figs 2 and 3). On the other hand, the values of f12, the detachment coefficient, are found to be generally higher for the data of Prasanthi et al. (1994), reflecting the immediate breakthrough observed in their experimental studies (Fig. 1). The value of fl', the headloss coefficient (which is simply the fraction of particles deposited in the filter bed contributing to the additional surface area for headloss development), is found to be generally lower for Chang's (1989) data, despite showing a significant increase in headloss with respect to time. This may be attributed to greater removal of particles in Chang's experiments, whereby the actual fraction of particles contributing towards additional surface area for headloss development is reduced. Also, the variation in fl' can be attributed to the size of the particles being filtered. As reported in several earlier studies (e.g. Dharmappa et al., 1992), the smaller the particles are, the higher the headloss is, as they contribute more surface area towards headloss development. In Cases 1~,, as the size of the particles is very small (Table 1), the value offl' can be expected to be high. CONCLUSIONS In this paper a methodology comprising the Marquardt (1964) algorithm for determining the empirical parameters of the modified O'Melia and Ali model (Vigneswaran and Chang, 1986) is presented and is tested against literature data. Through the case study it is found that a multiplication factor for the
Estimation of filtration coefficients sum of the square of errors belonging to effluent quality could resulll in a more balanced estimation of parameters, while achieving a better fit for both effluent quality and headloss, simultaneously. In the case study it is found that the shape of the fitted effluent quality curve has significant impact on the values of the model parameters. If only the total sum of the square of errors is considered, there are several possible ways of fitting the observed data with respect to a filtration process. It is necessary to check the shape of the fitted curve as well. It should conform to the shape of the observed data. Towards this end, the proposed algorithm appears to be well suited to give the best curve fitted to a set of both observed effluent quality and headloss data. The value of the multiplication factor should be fixed based on the accuracy required under both effluent quality and headloss predictions, and the difference in the observed absolute values of effluent quality and headloss. The model parameters predicted reflected, to some extent, the experimental conditions adopted by Chang (1989) and ]?rasanthi et al. (1994). The values of ctpfl, the attachment coefficient, for Chang's data (1989) are found to be higher, reflecting the better efficiency of pardcle removal observed during experimental studies. On the other hand, the values of f12, the detachment coefficient, are found to be generally higher for the results of Prasanthi et al. (1994), reflecting the immediate breakthrough observed in their experimental work. The values of fl', the headloss coefficient, for Chang's (1989) exper-
1091
iments are found to be lower as there are more particles being removed and hence the actual fraction of particles contributing to the headloss development is minimum. The higher values of fl' obtained for the results of Prasanthi et al. (1994) can also be partly attributed to the smaller size of the particles filtered.
REFERENCES
Chang J. S. (t989) Mathematical modeling of transient behaviour of deep bed filtration. D.Eng. dissertation, Asian Institute of Technology, Bangkok. Darby J. L., Attanasio R. E. and Lawler D. F. (1992) Filtration of heterodisperse suspensions: Modeling of particle removal and headloss. Wat. Res. 26, 711-726. Dharmappa H. B., Verink J., Fujiwara O. and Vigneswaran S. (1992) Optimisation of granular bed filtration treating polydispersed suspension. Wat. Res. 26, 130%1318. Draper N. R. and Smith H. (1981) Applied Regression Analysis. 2nd edn. John Wiley & Sons, New York. Kuester J. L. and Mize J. H. (1973) Optimisation Techniques. McGraw-Hill Book Co., New York. Marquardt D. M. (1963) An algorithm for least-squares estimation of nonlinear parameters. J. Soc. Industrial Appl. Mathematics 11, 431-441. O'Melia C. R. and Ali W. (1978) The role of retained particles in deep bed filtration. Prog. Wat. Technol. 10, 16%182. Prasanthi H., Vigneswaran S., Waite T. D. and Ben Aim R. (1994) Filtration of submicron particles: effect of ionic strength and organic substances. Wat. Sci. Technol. 30, 149-158. Saez P. B. and Rittmann B. E. (1992) Model-parameter estimation using least squares. Wat. Res. 26, 711-726. Vigneswaran S. and Chang J. S. (1986) Mathematical modeling of the entire cycle of deep bed filtration. Wat. Air Soil Pollut. 29, 155-164.