- Email: [email protected]
Moment theory of breakthrough curves for fixed-bed adsorbers and reactors
The Chemical Engineering Journal, 16 (1978) 211- 222 @ Elsevier Sequoia S.A., Lausanne - Printed in the Netherlands
Moment Theory of Breakthrough Curves for Fixed-bed Adsorbers and Reactors M.-S. RAZAVI, Department (Received
B. J. MCCOY and R. G. CARBONELL
of Chemical Engineering,
University
of California,
Davis, Calif. 95616 (U.S.A.)
14 March 1977; in final form 8 May 1978)
Abstract
sary. Optimization of the process depends on the relation of the breakpoint to the parameters of the system. For linear processes the Goldstein J function problem [ 1,2] is a widely used model which leads to a closed-form solution for the BTC. Although this model does not explicitly contain axial dispersion, intraparticle diffusion or details of adsorption, semi-empirical methods for combining such “resistances” have been proposed by Vermeulen et al. [ 21. Masamune and Smith [3] and Schneider and Smith [4] extended Rosen’s treatment [ 51 and derived an integral expression for the BTC for a model containing intraparticle diffusion and reversible adsorption but neglecting axial dispersion and volume capacity of porous particles. The more complete model [6,7] used by Smith and coworkers [4, 81 to evaluate kinetic and transport parameters by moment analysis of pulse response experiments, reduces to a single second order differential equation under a pseudo steady state condition. This same differential equation, which describes a variety of column processes, can be solved exactly, as we shall show, and it yields an expression for the BTC which We compare with the present theory. The approach studied in the present paper to obtain the BTC was introduced by Gelbin and coworkers [9] and is essentially the inverse of a technique used by Ammons et al. [lo, 111 to apply moment analysis to BTCs. Briefly, the basic idea is to use the convolution theorem to integrate the response to a 6 -function input and to obtain an expression for the BTC. The pulse response, characterized by its moments, is given by a series expansion in terms of Hermite polynomials. It is in the truncation of this series that the approximation enters. For most industrial-sized columns this approximation is satisfactory. In
For linear rate expressions in isothermal fixed-bed adsorbers or reactors, an approximate expression for the breakthrough curve (BTC) is obtained by applying the convolution theorem and in tegmting a Hermite polynomial series expansion for the pulse response .of the column. The approximate solution is compared with exact solutions such as that of the Goldstein J function model and a dispersed-flow column with homogeneous first order reaction. If the series is truncated so that only two moments appear, i.e. the first absolute moment pi and the second central moment j.i2, then the approximate model gives BTCs which are in good agreement with the exact solutions when u;/ (2l.i2Y2 > 2.5. The agreement is improved for longer columns and lower fluid velocities. For a model describing mass transfer to porous particles and intraparticle diffusion as well as adsorption, the approximate solution compares well with experimental data except for very short columns. The semi-empirical method of combining mass-transfer resistances is evaluated and is found to be quite accurate. An economic study ofpumping cost and product value is developed to optimize flow rates in adsorption processes.
1. INTRODUCTION
Prediction of the breakthrough curve (BTC), i.e. the output response to a step input, is necessary for analysis of unsteady processes in fixed-bed systems, e.g. adsorption or ionexchange columns, tubular reactors or heat regenerators. Knowing the BTC, we can define the breakpoint as a point on the BTC when regeneration of the bed becomes neces211
212
essence, the approach sacrifices the accuracy of BTCs for very short columns so that complex linear transport and kinetic processes can be described in longer columns. The final expression for the BTC is a relatively simple algebraic equation, so that no numerical integration is required. Our primary objective in the present study is to evaluate the accuracy of this procedure.
2. THEORY
The general forms of the governing equations for the models considered here, expressed in terms of the concentration c(t,z) in the fluid phase as a function of time and position along the column and the concentration ci( t,, r) in the solids in a certain column cross section as a function of time and radial position in a particle, are
ac -- + at
DC&;_
6?(C,CiR)
=
(1)
0
A(S) = al(S) + a2(S)a,(S)
-
a$
+ bi’(C,~*,CJ
= 0
(2)
The functions 61 and A’ are assumed to be linear. The initial and boundary conditions for c( t.2) are
(9)
The dimensionless Laplace-transformed boundary conditions are Q&z=
0) = C,(S)
(10)
E(S, 2 + -) is finite
(11)
For finite length columns the Danckwertz boundary conditions will provide a more accurate accounting of axial dispersion. However, considering that we shall limit the application to long columns and that we shall be interested in processes where mass transfer or chemical reaction effects make significant contributions to pulse broadening, we apply eqns. (10) and (11). Then eqn. (6) has a solution which, after the convolution theorem is applied, becomes C(T,Z)=JCB(T-r)L-‘exp(hZ)dr 0
(12)
where h is the negative root of X”/Pe- h + A(S) = 0. If C,(T) = cg(t,z
;;_
-23
= 0) = 6(T)
(13)
i.e. if the input is a Dirac delta function, then the integral in eqn. (12) gives C,(T,Z)
= L-l exp(A2)
(14)
c(t = 0,z) = q(z)
(3)
c(t,z = 0) = cg(t)
(4)
Thus the Laplace inversion of the exponential is simply the output response of the system to a delta function input. Now the BTC is the response to a step function input; thus if
c(t,z + -) is finite
(5)
C,(T) = C,(T,Z
= 0) = U(T)
(15)
where
Solving eqn. (2) with appropriate initial and boundary conditions allows cia to be eliminated from eqn. (1), which may be made dimensionless for convenience. The resulting equation in Laplace transform space has the form
then eqn. (12) will yield the expression for the BTC:
1 d2c .-- ----Pe dZ2
&,(T,Z)
dc dZ
+ A(S)6 = 0
(6)
where the Peclet number Pe = ud/D is based on the column diameter d. The expression A(S) may be evaluated as follows: if R(((;;, cia) = al(S)c
+ es(S)cia
(7)
and Cz = as(S)C then
(8)
U(T
..7)=
1
T>7
0
T<7
1 = fC,(r,Z) 0
d7
(16)
(17)
Therefore, an integral over the pulse response is what is needed to calculate the BTC. In many ca6es it is impossible to invert exp(X2) directly to C, (T.2); this is the limitation to obtaining the exact analytical expression of the BTC. However, C, (T, 2) can be approximated by a series expansion in the orthogonal Hermite polynomials H,(X); see
213
for example Kubin [6] , Carbonell and McCoy [12] and Radeke et al. [9] . Hermite polynomials are chosen since they have the weighting function exp(-X2); thus higher order terms in the series describe perturbations to a Gaussian distribution. The expansion can be written g A,(Z)e-X’Hn(X) n=O
Cb(T,Z)= M,(Z)
(18)
(l)), the front part of the peak should be nearly Gaussian, and the model should be adequate for calculating breakpoints. If the expansion (18) is substituted into eqn. (17), we obtain after integration, using properties of H,(X),
G(T,n
(
in terms of (2&,@/2
t-cc; =
-
2(2Ms)l’2 x
(=f2P2
exp(-X2)H,,_1(X)
T-M; x=
M;
.= ?4Mo erfX+erf
--
(1%
(2P2P2
and the dimensionless moments of concentration MO(~) = mo(e)/co
(20)
M; =p;v/L
(21)
M2 = p2v2/L2
(22)
Using the orthogonality property of H,,(X) we can calculate the expansion coefficients in eqn. (18) as
A0 = (2aM2)-1’2
(23)
Al=0
(24)
A2 = 0
(25)
(26)
-exp
(-$$)
If the expansion is truncated, this result has the unsatisfactory property that for long times (T+ - ) the concentration does not rise to its input value. The reason for this behavior is that part of C6(T,2) h in the region T < 0 owingto its representation as a perturbed Gaussian and is not included in the integ ration from 0 to = of eqn. (17). This deficiency is overcome by replacing the lower limit of integration in eqn. (17) by .--00; the result of the integration is [9]
C,(T,Z)=?hM,(Z)X X 1 +erfX
1 The first term in the expansion is Gaussian, and the higher order terms contribute distortions and asymmetries. Conditions for the uniform convergence of the series (18), called a Gram-Charlier series, are discussed by Ord [13]. Any function of order exp(--X2) can be represented by the series. Although we cannot prove that chromatographic output responses to Sfunction inputs are always nearly Gaussian, we recognize that in many experimental situations Gaussian outputs are observed. In fact, the assumption of a Gaussian shape is the basis for much analysis of chromatographic data [14] . We also should note that the worst deviation from the near-Gaussian shape is due to tailing, which can be caused by slow desorption. However, for determining breakpoints we need to know accurately the shape of the front of the peak, not of the tail. In the absence of channelling, which we cannot describe using the axial dispersion model (eqn.
Hn-r \--(2.$1,211,(27)
-2(2M2)l12
2 Ane-X2Hn-1(X) n=3 (28)
By replacing t in the interval (O,=) with t in (--O”, -), we include in the integral portions of the concentration curve which may not be zero for t Q 0. For this reason eqn. (27), because it omits this portion, never reaches its proper value as t + Q). However, eqn. (28) has the disadvantage that for small flow velocities or very short columns the concentration at T = 0 is not zero as it should be according to our initial condition. Instead, the BTC at T = 0 has the value C,(T=
O,Z)=
%M,(Z) l-erf 1
MM; (W2)1/2
I
(2g)
if we truncate before n = 3. Keeping the term involving the third moment improves the agreement but not substantially. Of course, if many terms were included in the expansion arbitrarily satisfactory behavior at T = 0 could be guaranteed, but handling such large expansions is.impracticable for the analytical solu-
214
tions we want. However, in eqn. (29) if p’,/ is greater than 2.5, C,(T= 0,Z) = 0, (2P P as required, and eqn. (28) with n < 3, and of course n > 3, will be accurate representations of the BTC. We might also add that eqn. (27) is accurate if p;/(2p2)112 > 2.5. The validity of the derived approximate expression eqn. (28) for predicting the BTC can best be determined by comparison with exact solutions of different mathematical models. We shall also compare our expressions with experimental data. The advantage of our method is that an analytical expression for the BTC, eqn. (28), can be obtained for any linear model by expressing its coefficients in terms of the moments obtained from that particular model. We realize that third and higher moments may be difficult to compute from experimental data; however, expressions for these moments in terms of the parameters of the system are known for many of the models discussed in this paper.
a discontinuity at t = z/u, which is due to the first order differential equations of the model. ,The solution representing the concentration in the solids,
t
I
0
Ci(t,Z)
=
1 --J{mk(t-z/u),k’z/u)
t > z/u
(37) is continuous because J(O,y ) = 1. Since the solution (J function) is represented by an integral, for simplikity of application it has been presented in the form of graphs [2] of C/c0 uersus mk(t -z/u) for different values of the parameter k’du. The mathematical behavior of the J function has been thoroughly discussed by Goldstein [l]. Table 1 shows the relation between the various parameters of different formulations.
TABLE 1
*
Correspondence of nomenclature in Figs. 2i6, 216(a) and 217 of Hougen and Watson [ 151, Fig. 16-20 of Vermeulen [ 21 and the present paper
3. SIMPLE MODEL FOR BTCs
The following equations define what we may call the Goldstein J function model: (30)
Hougen and Watson
Vermeulen
This paper
a2 br Y/Y 0
N NT X = J(N, NT)
k’z/v mk(t--z/u) c/co
NW0
y.
CilCO
In this table, b = mk and 7 = t -z/v.
(31) with initial condition,
c(t = 0,z) = 0
initial condition,
c*(t = 0,z) = 0
boundary
condition,
(32)
(33) c(t,z = 0) = ceU(t) (34)
To express our approximate solution in terms of the parameters of the system, we require expressions for moments, which are obtained from derivatives of the Laplacetransformed concentration:
Equations (30) and (31) obviously are a special case of eqns. (1) and (2). With the linear equilibrium relation c* = mci, the solution can be written
m, = (-l)Yillil
-=
mO(d=co
mL(t--z/u)} where the Jfunction is defined
d”c ds”
(38)
For the fluid phase, the moments for the Goldstein J function problem are (39)
J{k'z/u,
J&y)
= 1 - ey feSxIo
by
{2(xy)‘/2) dx
(36)
&(Z) = ;,! 742
0
Here IO is the modified Bessel function of the first hind. Since J(x, 0) = eqx, the solution has
32
= -2.z - k’ u m2k2
(41)
2k’
rl(3@)= G-msk3
(42)
For material adsorbed moments are ms,&)
=
onto the solid phase the (43)
co/m
(44)
2z k’ +2__
E12,s=
;
-
m2k2
(45)
m2k2
32 2k’ E13,S= ‘y ;$p
1
0.005
(46)
+m
0.001 1
2
5
10
20
50
100
ZOO
5OO
1OOO
mk(l-z/v)=NT
We have compared the approximate solution, eqn. (28) when n > 3, with the exact solution as presented by Hougen and Watson [15] in their Figs. 216(a), 216 and 217, and by Vermeulen [2] in Fig. 16-20. We note that the curves corresponding to values of k’/ u < 5 in Hougen and Watson’s [15] Figs. 216(a) and 216 are slightly in error [ 161. The error has been transferred from earlier studies [17, 181. However, Brinkley [19, 201 supports the accuracy of Vermeulen’s [2] Fig. 16-20. It is interesting to note that eqn. 16-169 of Vermeulen [ 21, empirically constructed as an approximation to the J function solution, is identical to our eqn.(28) with n < 3 when we use the moments of the Goldstein J function model, eqns. (39) - (42). We have compared our approximate solutions with exact graphical representations of the J function in Figs. 1 - 3. Figure 1 shows the effect of the lower limit of integration, 0 in eqn. (27) and --00 in eqn. (28), when n < 3 ’ 0.999 0.995 0.93
0.2 0.1 0.M O.O2
-J-Function
0.005 O.Mn 1
2
I
1’
5
10
I
1
I
II
2O
50
loo
ZOO
mk (t-z/v)
i I
1
6a
IOW
= NT
Fig. 1. Comparison of the approximate expressions, eqns. (27) and (28) for n < 3, with the J function solution given by Vermeulen [ 21 in his Fig. 16-20. Note: for k’z/u = 36, eqns. (27) and (28) coincide.
Fig. 2. Comparison of the approximate expressions, eqn. (28) for n < 3 and n = 3, with the J function solu. tion given by Vermeulen [ 21 in his Fig. 16-20. Note for k’zh = 9 and 35 and mk(t -z/v) 5 5, the approximate solution for n = 3 coincides with the J function.
1, / 0.2
IJ
01 1 0
-
J-Function
LVem~~len)
-..---___-_
J-Funclh n&3 “:3
(Brinkley)
’
1
1
I
I
I
I
0.2
0.4
0.6
0.8
1.0
1.2
1.4
t,
min
Fig. 3. Comparison of the approximate expressions in real time, eqn. (28) for n < 3 and n = 3, with the J function solution given by Vermeulen [ 21 and by Brinkley [20] ; k’du = 1, mk = 5 min-l, z/u = 0.1 min.
for selected values of k’zlu = N. For N 5 35 the agreement between the approximate and exact solutions is quite satisfactory. At N = 35, the estimate of c/c, varies from the exact value by less than 2% at the two extremes of the BTC. Even for N = 9, the estimate of the midpoint of the BTC is good. The tmsatisfactory aspect of using a lower limit of zero in the convolution integral is shown for N = 1 and 9, where the BTC does not rise to c/c, = 1. For N > 35 the effect of the different lower limits of integration, 0 or -00, is negligible. Figures 2 and 3 show how including the third moment in eqn. (28) improves the comparison with the exact curves. Except for N = 1 and a portion of N = 9, the values from eqn . (28) with n > 3 overlap the exact solution. However, even for N = k’du = 1 the agreement is fair, and this is also demon-
216
strated in Fig. 3, where a BTC is plotted for mk = 5 min-’ and z/u = 0.1 min. The tabulations of Brinkley [20] are used to extend the BTC to small values of time in Fig. 3. Equation (28) with n = 3 matches the exact BTC fairly well, except that it does not show the discontinuity at t = z/u = 0.1 min. The presence of unrealistic inflection points due to the third order term in the series is obvious. For the parameters of Fig. 3, p ;/ (2&+‘2 = 0.75, which constitutes a severe ,t.est of the theory. We have also compared eqn. (28) with moments for the adsorbed gas on the solid, eqns. (43) - (46). While the agreement between exact and approximate curves is fair for n < 3, the agreement is excellent for n = 3 for all values of k’z/v. From this comparison of the approximate BTC expressions with the exact Jfunction solution, we conclude that -00 rather than zero should be used as the lower limit of integration. The agreement is improved for longer columns and lower flow rates, i.e. for larger values of k’zh. For most industrial applications k’z/u is larger than 10, otherwise the separation is unlikely to be economical. Short columns with inlet pipe diameters much smaller than the column diameter will experience radial concentration gradients and radial dispersion, and will not fit the basic model of eqn. (1) at all.
4. A MODEL WITH A FIRST.ORDER TriANSFER
RATE OF
IN THE FLUID PHASE
Here we consider a model whose exact solution can easily be formulated and compared with the approximate solution. The model with a first order rate term approximately describes a variety of processes in tubular reactors with axial dispersion, e.g. a first order homogeneous reaction, a heterogeneous packed-column reactor where the rate of mass transfer from the flowing stream to the particles can be written linearly [21], the model of Suzuki and Smith [2$] for reaction and adsorption inside porous particles of a’ packed column when a pseudo steady state prevails inside the particles [ 221, and adsorption or chemical reaction on the surface of non-porous particles in the column under
conditions that a pseudo steady state applies at the surface [22]. Because of the loss of the reacting species, moments may not be easily measured by pulse response experiments, and must be estimated for a definite reactor model such as the simple one proposed hem. In dimensionless variables we have 1 a2c ac _._ _--_-__
Pe az2
ac
az
aT
-Kc=0
(47)
with initial condition,
C(T = 0,Z)
boundary
condition,
C(T, 2 = 0) = U(T)
boundary
condition,
C(T, 2 + -)
The exact solution
= 0
(43)
is finite (50)
is
C(T,Z)=+(,exp[y\l-(l+4$2\] 2 --~.
x erfc
(49)
X
T(l
+ 4~ /Pe)1’2
+
(4T/Pe)l12 +exp[iZPe\l+
X
2 +
erfc
(1 +k)r2]X
T(l +
4~ /Pe)1’2
(4T/Pe)1/2-
(51)
which has been obtained by Laplace transformation by Cho [23] to describe transport ’ processes in soil. The approximate solution (n < 3) to this model can be obtained by substituting the following expressions of the moments into eqn. (28): Ms(Z)=exp[i-ZPe]l--(l+E)Liql
(52)
z (53)
M; @) = (1 + 4K /pe)l/2
(54)
M3(Z)=Pe
-b/2
122
1+; (
1
(55)
217
This gives C(T,Z)=,-iexp[iZPe/l--(I
2 - T(l + 4s /Pe)lj2
+z)“2/]
L
Pe
\~
The second term inside the brackets is of the order of the third moment. Next we compare the two solutions, eqns. (51) and (56), and assess the accuracy of the approximate solution. For the limiting case when Pe -+ 00, i.e. when axial dispersion has no significant effect on the concentration distribution in the column, the two solutions become identical and equal to fexp( C(T, 2) =
+
--KZ)U(T---2)
wz
(57)
T=Z
eXp(-KZ)
For moderate values of Pe and 2, the graphs corresponding to exact and approximate solutions are plotted in Figs. 4 - 6. The comparison of the curves indicates that for high values of Pe and 2 the approximation is quite good. In Fig. 4 the exact expression for the BTC is plotted for 2 = 10 and K = 0.1 with a range of Peclet numbers. Decreasing Pe, which implies more axial dispersion, flattens the curves and rounds the sharp comer of the curve corresponding to Pe + 00. Only the approximate solution at Pe = 1 is presented (as dots) in Fig. 4, since it is the least accurate and is in error by only 3% at most. The error decreases for large values of Pe, but for smaller values of 2, i.e. z < lOd, the error would be greater.
(56) Figure 5 demonstrates the shape of the BTC in different locations of the bed for Pe = 4OandK = 0.1, and shows C(T + =, 2) decreasing with increasing 2. The approximate solution, shown by solid dots, is in good agreement with the exact solution even when the column is short, 2 = 1. The broken lines represent the exact BTCs when D = 0, eqn. (57). The locus of comers is where T = 2 and thus from eqn. (57) is given by exp(-O.lT). When both Pe and 2 are small, as in Fig. 6, the approximate solution obviously deviates from the exact. The approximate solution including the third moment has three inflection points, while the solution excluding the third moment has one inflection point. Although neither approximation fits the exact curve very well, the case represented by these values of Pe and 2 is hardly representative of an industrial process. High values of Pe and 2 are equivalent to high values of the parameter ~;/(2p~)l’~ = M;(2M2)1’2 since in this model
Pi (2/J2)1’2
= J- (PeZ)l12 2
4K (1 + --)
“* (58)
Large values of Pe and 2 also imply that the third moment is negligible. For packed columns, gases and liquids usually have values of Pe ranging from 20 to 40; see for example Vermeulen [ 21. Also, we usually have 2 % 1 since Z = z/d.
0.4
s ”
0.2
-Exact l
0
Approximate
0.1
0 6
a
10 T=tv/d
Fig. 4. BTCk for 2 = z/d = 10 and
‘*
14
16
K = 0.1 and different values of Pe = vd/D. Exact, eqn. (51 .), and approximate .tions, eqn. (56). nearly coincide except for Pe = 1 where approximate values are shown as dots.
218
+_-__ and the coupling
boundary
condition
(6%
= k!f (C -- CiR)
The moments for a 6 -function beginning at t = 0 are
P
‘, ;:
inlet condition
p; = ;(1+6,) 0
12
3
4
5
6
18
9
1
I
I
10
11
12
and
T=tvld
22
Fig. 5. BTC!s for various positions in the column with Pe = 40 and K = 0.1. The approximate solution is represented as dots and eqn. (57) is shown as broken lines.
P2
=--y
6 1 +
D(l + 6,)’ ;s
I
where (65) and
x(&+&)1 -3
-2
I
I
I
I
I
I
I
-1
0
1
2
3
4
5
6
T=tr/d
Fig. 6. BTCs for a short column, z/d = 1, with Pe = 1 and K = 0.1: -, exact solution, eqn. (51); - - -, approximate solution, eqn. (56) when R < 3 (/.I3 = 0); ------ 9approximate solution, eqn. (10) when n = 3.
5. A MODEL DESCRIBING DETAILS OF ADSORPTION INSIDE PARTICLES
The model developed by Kubin [6] and Kucera [7] and used extensively by Smith and coworkers [4,8] takes into account the effects of (a) axial dispersion, (b) fluid to particle mass transfer, (c) intraparticle diffusion and (d) reversible adsorption. The chief restriction is that the transport and adsorption processes be linear. The governing differential equations are D~~--+~-
0
hi
pP
acads
P
at
- -at -__
-
(59) =0
(60)
(61)
(66)
Ammons et al. [ 10,111 obtained rate parameters for adsorption of methyl mercuric chloride on activated carbon by the moment technique. Because the capacity of the carbon for the adsorbate was so large (large value of KX), pulse experiments were not feasible. Ammons carried out the experiment using a step function input and related the calculations of the moments to the analysis of the experimental BTCs. Values of the first moment p; and the second central moment p2 are presented by Ammons [ll] and may be used in the approximate theory to calculate BTCs. Since Ammons did not calculate pa, we use eqns. (27) and (28) with n < 3 to plot the BTCs of Figs. 7 and 8. For the BTCs in Fig. 7, the ratio pi/ is 2.66 for BTC 4 and 1.88 for BTC 7, (2P2 P2 and the agreement between the approximate theory, eqn. (28), and the experimental data is good, with better agreement for larger values of the ratio. In Fig. 8, P;/(@~)~‘~ equals 1.11. Since this ratio is substantially less than 2.5, the agreement is less satisfactory. In this figure we have also plotted eqn. (27) with n < 3, the solution obtained when the lower limit of integration is zero. Although this curve is a better approximation for small times, it is poorer for large times.
219
0
tx10-3,
SK
Fig. 7. Comparison of the approximate expression for BTC (dots), eqn. (28) for n < 3, with experimental results (line), BTCs 4 and 7 of Ammons [ 111.
processes can be combined without the aid of a detailed theory such as that described by eqns. (59) - (62). The overall transfer rate constant kOp represents the combined resistance and can be used in a mathematical model such as the Goldstein J function problem [l] . To assess the accuracy of this technique, we examine expressions of moments that can be obtained from the combination of resistances in the Goldstein J function problem and compare them with results of the detailed model based on equations describing the actual transport phenomena, i.e. eqns. (59) - (62). According to Vermeulen 123 we make the following substitutions in eqns. (30) and (31): FE’= ko, Aa&
(67)
and k = kopap
(68)
where A = I!&,(1
- Ly)
(69)
The overall transfer constant is expressed by 0 0
4
8
12
16
20
24
t x 10-3, SC
Fig. 8. Comparison of the approximate expressions (n < 3) eqn. (27) (-------) and eqn. (28) (- - -) with the experimental BTC 10 (-) of Ammons [ 111.
Another theory for the BTC which involves a numerical integration has been presented [3 - 51. This approach is based on the approximations D = 0 and flaci/a t = 0 in eqns. (59) and (60), namely that axial dispersion and volume capacity of the porous particles are negligible. Neglecting the volume capacity alters the first moment and shifts the BTC to smaller values of t; the shift is less than 6% for the experiments discussed by Schneider and Smith [4] . Neglecting axial dispersion reduces the second moment by as much as 70% for their experiments and thereby sharpens the BTCs rather unrealistically. Since for the conditions of experiment under consideration the approximations in our theory have much less effect than these discrepancies, we forego a detailed comparison (see ref. 22 for details). Vermeulen [2] suggests that mass transfer resistances involved in packedcolumn
1 1 _=_-++._+ kop kcp
1
a, -
kcf
k
(70)
where (71)
k,, = kp + k,,,lA and K = k:,,/K:,
(72)
When axial dispersion is included in the suggested manner 1 -= kcf
1 -+-
1
(73)
k:,
kf
in terms of kb = au 2/oap
(74)
we obtain k,
=
K:p,(l
---a) --a
+
+
-1
=XP~
3kf
5 + R2KXp, k’ads 15Di
+
KXpp(l --a) -$
t
(75)
From eqns. (40) and (41) the moments for the Goldstein J function problem become /.& =f
1 --a 1 + -.-a-
(76)
220
and
tb =p;
z 1 --a PJ;~ cl2 = -.v
R2Ki2p;
_+---+
-a
(
R2p;KX2 +--+3kfR
p;KX2
aD
CY2
u2
(77)
For comparing these estimated moments with the calculated moments according to the model (eqns. (59) - (62)) we note that p,Ki/P is the critical parameter. If pp K;//? S 1 then the expressions for moments are identical; however, if pP KX/fl Q 1 the expressions are entirely different. Since the former inequality usually holds in adsorption processes, the combinations of resistances technique must be considered quite accurate. To illustrate this, consider the adsorption of ethane from helium onto silica gel with the parameters measured by Schneider and Smith [4]. For c(\, eqn. (76) gives 84.61 s while eqn. (63) gives 86.69 s, a difference of 2.4%. For p2, eqn. (77) gives 414.1 s2 while eqn. (64) gives 445.3 s2, a difference of 7%. In this example KX = 14.6 ml g-’ silica gel so that ppKXI p = 33.94 which is sufficiently greater than unity to give excellent agreement between the estimated and correct values of the moments. For adsorption of propane and n-butane the equilibrium adsorption constants are larger than for ethane and the agreement is even better. The third moments, however, are poorly estimated, and thus the estimation procedure will not accurately predict the skewness of pulses or BTCs. However, since for most industrial processes pP Ka/fla 1, we may conclude that the estimation method is quite accurate for the analysis and design of adsorption in packed columns.
6. THE BREAKPOINT
AND OPTIMIZATION
(78)
in terms of a constant a that is related to cb by
1503
&ads
-adp2
OF
COLUMN OPERATION
A most important design parameter in the adsorption process is the time t,,, corresponding to the breakpoint, that fixes how long the column should be in operation. The time t,, occurs when an arbitrarily small but signif* icant effluent concentration cb appear%at the end of the column. Since the first moment locates the center of the BTC and the second moment the width, for a symmetric BTC the breakpoint can be defined as
(79) which follows from eqn. (28), it < 3. Thus if cb is known, a can be calculated from eqn. (79) and hence tb from eqn. (78). To maximize the utilization of a column before it is regenerated, we wish to maximize tb. From eqn. (78), we see that the parameters that increase P; and/or decrease p2 will increase tb. Knowing how these moments depend on the parameters, e.g. as in eqns. (63) - (66), will enable us to optimize tb. By comparing the effect of the rate parameters on p; and p2, we can conclude that increasing Q, k&s, kf and decreasing R2 and D will increase tb. To determine the optimum velocity, or rate of throughput, we express the net income of an adsorption column of given length as the difference between the value of the product and the pumping cost. We assume that other costs are independent of the velocity. The value of the product, 7 tbuo”d2/4, is given in terms of the price y per unit volume of product and the throughput after time tb. The pumping cost WV; follows from the equation for the pressure drop in a packed column. The net income function is then
(80) which is to be maximized with respect to v = V&J. As an illustration, we use the expressions for moments, eqns. (63) - (66), along with the parameters for the adsorption of ethane from helium onto silica gel given by Schneider and Smith [4]. To give a more realistic commercial meaning to the product price we considered a column cross section 100 times larger than that used in the laboratory apparatus. Values of o and cuy are chosen as 21.0 $ crnm2 s2 and 0.05 $ cm-‘, respectively . The net income function f for various values of a, along with the pumping cost, are plotted in Fig. 9. The optimum velocity increases for increasing a, i.e. for increasing purity, while the net income decreases. The positions of the optimum velocity for differ-
221
(see also ref. 24). The moment expressions and present theory constitute a foundation for optimizing BTCs for column processes.
NOMENCLATURE i
a UP
0
2
4
6
0
10
12
14
v, cm/sac
Fig. 9. Plot of income us. velocity for various purities of effluent. The broken line shows the values of optimum velocity.
ent values of a are shown by the broken line. For values of a 2 11.5 there would be a negative net income owing to the high purity required of cb. It is interesting that the range of optimum velocity is less than 5 cm s-l in this case of adsorption of ethane on silica gel.
c(t, 2) cad&, r) c,(t, z= 0) q(t = 0,z) Ci(t,
CO c
7. CONCLUDING REMARKS
A supposition of the work described herein is that exact analytical, or more usually numerical, solutions for BTCs provide in many cases more information than is needed for the design and optimization of industrial column processes. For such processes governed by linear rate expressions, and therefore linear differential equations, the moment method, combined with the Duhamel theorem, presents a solution that is relatively easily obtained and simple in form (an algebraic expression). For many industrial columns, where @C(;/(~P#~ > 2.5 is generally satisfied, only terms up to the second central moment are needed to represent the BTC accurately. The first moment II> locates the center of the BTC and the second central moment p2 determines the spread of the BTC. For shorter columns with larger velocities, the third moment can be computed and this will extend the accuracy of the BTC representation by including asymmetries. For example, Kubin [ 61 has presented expressions for the third central moment for the adsorption model used by Schneider and Smith [4]
6
D Di d
f &I
k, k’
r, 2)
coefficient in expression for tb, eqn. (78) surface area of 30 -0)/R, spherical sorbent particles per unit volume of bed, 1 cm-’ concentration in interparticle space, mol cmm3 concentration of adsorbed species based on unit void volume of particles, mol cmm3 concentration input, mol cm-l initial concentration profile in interparticle space, mol cme3 concentration in intraparticle pore space, mol cmm3 characteristic constant concentration, mol cmW3 c/co, dimensionless form of preceding concentrations dimensionless concentration response to a unit step function dimensionless concentration response to a Dirac delta function input Laplace transform of C effective axial dispersion coefficient cm2 s-l effeciive intraparticle diffusion coefficient, cm2 s-l diameter of cylindrical column, cm net income, $ Hermite polynomial overall transfer rate constants in Goldstein J function problem, S-l
k’ads
adsorption
rate constant,
cm3
g-ls-1 kd kP
k POW k
desorption rate constant, 1 s-l particle phase diffusion rate constant, 1 s-l rate constant of pore diffusion, Is-l mass transfer coefficient at sphere surface, cm s-l k:,,/k,, adsorption equilibrium
222
L m m,(z)
M,(Z) M;(Z) M,(Z)
Pe r R 63
61’
tb
T u UC i
constant, cm3 g-l length of column c*/cI, equilibrium constant c(t,z) t “dt, nth moment of J 0 concentration m&)/co, zeroth dimensionless moment p; u/L, first dimensionless moment ~12 u2/L2, second central moment, dimensionless ud/D, Peclet number for axial dispersion radial coordinate in the spherical particles, cm radius of the spherical particles, cm di(c, CiR), term in eqn. (1) representing mass transfer from the interparticle fluid to the surface of particles 61’(Cads, Ci), telXl in eqn. (2) representing chemical reaction or adsorption within the pores of particles Laplace transform parameter, 1 s- ’ dimensionless Laplace transform parameter time dorresponding to the breakpoint, s td/u, dimensionless time interstitial velocity of fluid, cm s-l cw, superficial velocity, cm s-l length coordinate of the bed, cm z/d, dimensionless length coordinate
Greek symbols a inter-particle void fraction of the packed column porosity, intraparticle void fraction 0 I.4
1 m,
-
s ct dt, first absolute moment 0 0
‘0 Jc(t
moo moment
- p;)ndt,
nth central
t P PP
r/R, dimensionless radial coordinate fluid density, g cm- 3 apparent particle density, g cm- 3
REFERENCES 1 S. Goldstein, Proc. R. Sot. London, Ser. A, 219 (1953) 151. 2 T. Vermeulen, G. Klein and N. K. Hiester, in R. H. Perry and C. H. Chilton (eds.), Chemical Engineer’s Handbook, 5th edn., McGraw-Hill, New York, 1973, Sect. 16. 3 S. Masamune and J. M. Smith,A.Z.Ch.E. J., 11 (1965) 34. 4 P. Schneider and J. M. Smith,A.Z.Ch.E. J., 14 (1968) 762. 5 J. B. Rosen, J. Chem. Phys., 20 (1952) 387. 6 M. Kubin, Coil. Czech. Chem. Commun., 30 (1964) 1104, 2900. 7 E. Kucera, J. Chromotogr., 19 (1965) 237. 8 T. Furusawa, M. Suzuki and J. M. Smith, Catal. Rev., 13 (1976) 43. 9 K. Rade ke, K. Wiedemann and D. Gelbin, Chem. Tech. (Leipzig), 28 (1976) 476. 10 R. D. Ammons, N. A. Dougharty and J. M. Smith, Znd. Eng. Chem. Fundam., 16 (1977) 263. 11 R. D. Ammons, Adsorption of methyl mercuric chloride on activated carbon, M.S. Thesis, Univ. Calif., Davis, 1975. 12 R. G. Carbonell and B. J. McCoy, Chem. Eng. J., 9 (1975) 115. 13 J. K. Grd, Families of Frequency Distributions, Hafner Publ. Co., New York, 1972, p. 28. 14 J. C. Giddings, Dynamics of Chromatography, Part I: Principles and Theory, Marcel Dekker, New York, 1965. 15 0. A. Hougen and K. M. Watson, Chemical Process Principles, Wiley, New York, 1946 ; p. 1080. 16 A. Klinkenberg, Znd. Eng. Chem., 46 (1954) 2285. 17 T. E. W. Schumann, J. Fmnklin Inst., 208 (1929) 405. 18 C. C. Fumas, Tmns. A.Z.Ch.E., 29 (1930) 142. 19 S. R. Brinkley and R. F. Brinkley, Math: Tables Aids Comp., 2 (1947) 221. 20 S. R. Brinkley, Math: Tables Aids Camp., 6 (1952) 40. 21 M. Suzuki and J. M. Smith, Chem. Eng. Sci., 26 (1971) 221. 22 M.-S. Razavi, A method for predicting breakthrough curves in packed column processes, M.S. Thesis, Univ. Calif., Davis, 1978. 23 C. M. Cho, Can. J. Soil&i., 51 (1971) 339. 24 B. J. McCoy and R. G. Carbonell, A.Z.Ch.E. J., 24 (1978) 159.
Moment Theory of Breakthrough Curves for Fixed-bed Adsorbers and Reactors M.-S. RAZAVI, Department (Received
B. J. MCCOY and R. G. CARBONELL
of Chemical Engineering,
University
of California,
Davis, Calif. 95616 (U.S.A.)
14 March 1977; in final form 8 May 1978)
Abstract
sary. Optimization of the process depends on the relation of the breakpoint to the parameters of the system. For linear processes the Goldstein J function problem [ 1,2] is a widely used model which leads to a closed-form solution for the BTC. Although this model does not explicitly contain axial dispersion, intraparticle diffusion or details of adsorption, semi-empirical methods for combining such “resistances” have been proposed by Vermeulen et al. [ 21. Masamune and Smith [3] and Schneider and Smith [4] extended Rosen’s treatment [ 51 and derived an integral expression for the BTC for a model containing intraparticle diffusion and reversible adsorption but neglecting axial dispersion and volume capacity of porous particles. The more complete model [6,7] used by Smith and coworkers [4, 81 to evaluate kinetic and transport parameters by moment analysis of pulse response experiments, reduces to a single second order differential equation under a pseudo steady state condition. This same differential equation, which describes a variety of column processes, can be solved exactly, as we shall show, and it yields an expression for the BTC which We compare with the present theory. The approach studied in the present paper to obtain the BTC was introduced by Gelbin and coworkers [9] and is essentially the inverse of a technique used by Ammons et al. [lo, 111 to apply moment analysis to BTCs. Briefly, the basic idea is to use the convolution theorem to integrate the response to a 6 -function input and to obtain an expression for the BTC. The pulse response, characterized by its moments, is given by a series expansion in terms of Hermite polynomials. It is in the truncation of this series that the approximation enters. For most industrial-sized columns this approximation is satisfactory. In
For linear rate expressions in isothermal fixed-bed adsorbers or reactors, an approximate expression for the breakthrough curve (BTC) is obtained by applying the convolution theorem and in tegmting a Hermite polynomial series expansion for the pulse response .of the column. The approximate solution is compared with exact solutions such as that of the Goldstein J function model and a dispersed-flow column with homogeneous first order reaction. If the series is truncated so that only two moments appear, i.e. the first absolute moment pi and the second central moment j.i2, then the approximate model gives BTCs which are in good agreement with the exact solutions when u;/ (2l.i2Y2 > 2.5. The agreement is improved for longer columns and lower fluid velocities. For a model describing mass transfer to porous particles and intraparticle diffusion as well as adsorption, the approximate solution compares well with experimental data except for very short columns. The semi-empirical method of combining mass-transfer resistances is evaluated and is found to be quite accurate. An economic study ofpumping cost and product value is developed to optimize flow rates in adsorption processes.
1. INTRODUCTION
Prediction of the breakthrough curve (BTC), i.e. the output response to a step input, is necessary for analysis of unsteady processes in fixed-bed systems, e.g. adsorption or ionexchange columns, tubular reactors or heat regenerators. Knowing the BTC, we can define the breakpoint as a point on the BTC when regeneration of the bed becomes neces211
212
essence, the approach sacrifices the accuracy of BTCs for very short columns so that complex linear transport and kinetic processes can be described in longer columns. The final expression for the BTC is a relatively simple algebraic equation, so that no numerical integration is required. Our primary objective in the present study is to evaluate the accuracy of this procedure.
2. THEORY
The general forms of the governing equations for the models considered here, expressed in terms of the concentration c(t,z) in the fluid phase as a function of time and position along the column and the concentration ci( t,, r) in the solids in a certain column cross section as a function of time and radial position in a particle, are
ac -- + at
DC&;_
6?(C,CiR)
=
(1)
0
A(S) = al(S) + a2(S)a,(S)
-
a$
+ bi’(C,~*,CJ
= 0
(2)
The functions 61 and A’ are assumed to be linear. The initial and boundary conditions for c( t.2) are
(9)
The dimensionless Laplace-transformed boundary conditions are Q&z=
0) = C,(S)
(10)
E(S, 2 + -) is finite
(11)
For finite length columns the Danckwertz boundary conditions will provide a more accurate accounting of axial dispersion. However, considering that we shall limit the application to long columns and that we shall be interested in processes where mass transfer or chemical reaction effects make significant contributions to pulse broadening, we apply eqns. (10) and (11). Then eqn. (6) has a solution which, after the convolution theorem is applied, becomes C(T,Z)=JCB(T-r)L-‘exp(hZ)dr 0
(12)
where h is the negative root of X”/Pe- h + A(S) = 0. If C,(T) = cg(t,z
;;_
-23
= 0) = 6(T)
(13)
i.e. if the input is a Dirac delta function, then the integral in eqn. (12) gives C,(T,Z)
= L-l exp(A2)
(14)
c(t = 0,z) = q(z)
(3)
c(t,z = 0) = cg(t)
(4)
Thus the Laplace inversion of the exponential is simply the output response of the system to a delta function input. Now the BTC is the response to a step function input; thus if
c(t,z + -) is finite
(5)
C,(T) = C,(T,Z
= 0) = U(T)
(15)
where
Solving eqn. (2) with appropriate initial and boundary conditions allows cia to be eliminated from eqn. (1), which may be made dimensionless for convenience. The resulting equation in Laplace transform space has the form
then eqn. (12) will yield the expression for the BTC:
1 d2c .-- ----Pe dZ2
&,(T,Z)
dc dZ
+ A(S)6 = 0
(6)
where the Peclet number Pe = ud/D is based on the column diameter d. The expression A(S) may be evaluated as follows: if R(((;;, cia) = al(S)c
+ es(S)cia
(7)
and Cz = as(S)C then
(8)
U(T
..7)=
1
T>7
0
T<7
1 = fC,(r,Z) 0
d7
(16)
(17)
Therefore, an integral over the pulse response is what is needed to calculate the BTC. In many ca6es it is impossible to invert exp(X2) directly to C, (T.2); this is the limitation to obtaining the exact analytical expression of the BTC. However, C, (T, 2) can be approximated by a series expansion in the orthogonal Hermite polynomials H,(X); see
213
for example Kubin [6] , Carbonell and McCoy [12] and Radeke et al. [9] . Hermite polynomials are chosen since they have the weighting function exp(-X2); thus higher order terms in the series describe perturbations to a Gaussian distribution. The expansion can be written g A,(Z)e-X’Hn(X) n=O
Cb(T,Z)= M,(Z)
(18)
(l)), the front part of the peak should be nearly Gaussian, and the model should be adequate for calculating breakpoints. If the expansion (18) is substituted into eqn. (17), we obtain after integration, using properties of H,(X),
G(T,n
(
in terms of (2&,@/2
t-cc; =
-
2(2Ms)l’2 x
(=f2P2
exp(-X2)H,,_1(X)
T-M; x=
M;
.= ?4Mo erfX+erf
--
(1%
(2P2P2
and the dimensionless moments of concentration MO(~) = mo(e)/co
(20)
M; =p;v/L
(21)
M2 = p2v2/L2
(22)
Using the orthogonality property of H,,(X) we can calculate the expansion coefficients in eqn. (18) as
A0 = (2aM2)-1’2
(23)
Al=0
(24)
A2 = 0
(25)
(26)
-exp
(-$$)
If the expansion is truncated, this result has the unsatisfactory property that for long times (T+ - ) the concentration does not rise to its input value. The reason for this behavior is that part of C6(T,2) h in the region T < 0 owingto its representation as a perturbed Gaussian and is not included in the integ ration from 0 to = of eqn. (17). This deficiency is overcome by replacing the lower limit of integration in eqn. (17) by .--00; the result of the integration is [9]
C,(T,Z)=?hM,(Z)X X 1 +erfX
1 The first term in the expansion is Gaussian, and the higher order terms contribute distortions and asymmetries. Conditions for the uniform convergence of the series (18), called a Gram-Charlier series, are discussed by Ord [13]. Any function of order exp(--X2) can be represented by the series. Although we cannot prove that chromatographic output responses to Sfunction inputs are always nearly Gaussian, we recognize that in many experimental situations Gaussian outputs are observed. In fact, the assumption of a Gaussian shape is the basis for much analysis of chromatographic data [14] . We also should note that the worst deviation from the near-Gaussian shape is due to tailing, which can be caused by slow desorption. However, for determining breakpoints we need to know accurately the shape of the front of the peak, not of the tail. In the absence of channelling, which we cannot describe using the axial dispersion model (eqn.
Hn-r \--(2.$1,211,(27)
-2(2M2)l12
2 Ane-X2Hn-1(X) n=3 (28)
By replacing t in the interval (O,=) with t in (--O”, -), we include in the integral portions of the concentration curve which may not be zero for t Q 0. For this reason eqn. (27), because it omits this portion, never reaches its proper value as t + Q). However, eqn. (28) has the disadvantage that for small flow velocities or very short columns the concentration at T = 0 is not zero as it should be according to our initial condition. Instead, the BTC at T = 0 has the value C,(T=
O,Z)=
%M,(Z) l-erf 1
MM; (W2)1/2
I
(2g)
if we truncate before n = 3. Keeping the term involving the third moment improves the agreement but not substantially. Of course, if many terms were included in the expansion arbitrarily satisfactory behavior at T = 0 could be guaranteed, but handling such large expansions is.impracticable for the analytical solu-
214
tions we want. However, in eqn. (29) if p’,/ is greater than 2.5, C,(T= 0,Z) = 0, (2P P as required, and eqn. (28) with n < 3, and of course n > 3, will be accurate representations of the BTC. We might also add that eqn. (27) is accurate if p;/(2p2)112 > 2.5. The validity of the derived approximate expression eqn. (28) for predicting the BTC can best be determined by comparison with exact solutions of different mathematical models. We shall also compare our expressions with experimental data. The advantage of our method is that an analytical expression for the BTC, eqn. (28), can be obtained for any linear model by expressing its coefficients in terms of the moments obtained from that particular model. We realize that third and higher moments may be difficult to compute from experimental data; however, expressions for these moments in terms of the parameters of the system are known for many of the models discussed in this paper.
a discontinuity at t = z/u, which is due to the first order differential equations of the model. ,The solution representing the concentration in the solids,
t
I
0
Ci(t,Z)
=
1 --J{mk(t-z/u),k’z/u)
t > z/u
(37) is continuous because J(O,y ) = 1. Since the solution (J function) is represented by an integral, for simplikity of application it has been presented in the form of graphs [2] of C/c0 uersus mk(t -z/u) for different values of the parameter k’du. The mathematical behavior of the J function has been thoroughly discussed by Goldstein [l]. Table 1 shows the relation between the various parameters of different formulations.
TABLE 1
*
Correspondence of nomenclature in Figs. 2i6, 216(a) and 217 of Hougen and Watson [ 151, Fig. 16-20 of Vermeulen [ 21 and the present paper
3. SIMPLE MODEL FOR BTCs
The following equations define what we may call the Goldstein J function model: (30)
Hougen and Watson
Vermeulen
This paper
a2 br Y/Y 0
N NT X = J(N, NT)
k’z/v mk(t--z/u) c/co
NW0
y.
CilCO
In this table, b = mk and 7 = t -z/v.
(31) with initial condition,
c(t = 0,z) = 0
initial condition,
c*(t = 0,z) = 0
boundary
condition,
(32)
(33) c(t,z = 0) = ceU(t) (34)
To express our approximate solution in terms of the parameters of the system, we require expressions for moments, which are obtained from derivatives of the Laplacetransformed concentration:
Equations (30) and (31) obviously are a special case of eqns. (1) and (2). With the linear equilibrium relation c* = mci, the solution can be written
m, = (-l)Yillil
-=
mO(d=co
mL(t--z/u)} where the Jfunction is defined
d”c ds”
(38)
For the fluid phase, the moments for the Goldstein J function problem are (39)
J{k'z/u,
J&y)
= 1 - ey feSxIo
by
{2(xy)‘/2) dx
(36)
&(Z) = ;,! 742
0
Here IO is the modified Bessel function of the first hind. Since J(x, 0) = eqx, the solution has
32
= -2.z - k’ u m2k2
(41)
2k’
rl(3@)= G-msk3
(42)
For material adsorbed moments are ms,&)
=
onto the solid phase the (43)
co/m
(44)
2z k’ +2__
E12,s=
;
-
m2k2
(45)
m2k2
32 2k’ E13,S= ‘y ;$p
1
0.005
(46)
+m
0.001 1
2
5
10
20
50
100
ZOO
5OO
1OOO
mk(l-z/v)=NT
We have compared the approximate solution, eqn. (28) when n > 3, with the exact solution as presented by Hougen and Watson [15] in their Figs. 216(a), 216 and 217, and by Vermeulen [2] in Fig. 16-20. We note that the curves corresponding to values of k’/ u < 5 in Hougen and Watson’s [15] Figs. 216(a) and 216 are slightly in error [ 161. The error has been transferred from earlier studies [17, 181. However, Brinkley [19, 201 supports the accuracy of Vermeulen’s [2] Fig. 16-20. It is interesting to note that eqn. 16-169 of Vermeulen [ 21, empirically constructed as an approximation to the J function solution, is identical to our eqn.(28) with n < 3 when we use the moments of the Goldstein J function model, eqns. (39) - (42). We have compared our approximate solutions with exact graphical representations of the J function in Figs. 1 - 3. Figure 1 shows the effect of the lower limit of integration, 0 in eqn. (27) and --00 in eqn. (28), when n < 3 ’ 0.999 0.995 0.93
0.2 0.1 0.M O.O2
-J-Function
0.005 O.Mn 1
2
I
1’
5
10
I
1
I
II
2O
50
loo
ZOO
mk (t-z/v)
i I
1
6a
IOW
= NT
Fig. 1. Comparison of the approximate expressions, eqns. (27) and (28) for n < 3, with the J function solution given by Vermeulen [ 21 in his Fig. 16-20. Note: for k’z/u = 36, eqns. (27) and (28) coincide.
Fig. 2. Comparison of the approximate expressions, eqn. (28) for n < 3 and n = 3, with the J function solu. tion given by Vermeulen [ 21 in his Fig. 16-20. Note for k’zh = 9 and 35 and mk(t -z/v) 5 5, the approximate solution for n = 3 coincides with the J function.
1, / 0.2
IJ
01 1 0
-
J-Function
LVem~~len)
-..---___-_
J-Funclh n&3 “:3
(Brinkley)
’
1
1
I
I
I
I
0.2
0.4
0.6
0.8
1.0
1.2
1.4
t,
min
Fig. 3. Comparison of the approximate expressions in real time, eqn. (28) for n < 3 and n = 3, with the J function solution given by Vermeulen [ 21 and by Brinkley [20] ; k’du = 1, mk = 5 min-l, z/u = 0.1 min.
for selected values of k’zlu = N. For N 5 35 the agreement between the approximate and exact solutions is quite satisfactory. At N = 35, the estimate of c/c, varies from the exact value by less than 2% at the two extremes of the BTC. Even for N = 9, the estimate of the midpoint of the BTC is good. The tmsatisfactory aspect of using a lower limit of zero in the convolution integral is shown for N = 1 and 9, where the BTC does not rise to c/c, = 1. For N > 35 the effect of the different lower limits of integration, 0 or -00, is negligible. Figures 2 and 3 show how including the third moment in eqn. (28) improves the comparison with the exact curves. Except for N = 1 and a portion of N = 9, the values from eqn . (28) with n > 3 overlap the exact solution. However, even for N = k’du = 1 the agreement is fair, and this is also demon-
216
strated in Fig. 3, where a BTC is plotted for mk = 5 min-’ and z/u = 0.1 min. The tabulations of Brinkley [20] are used to extend the BTC to small values of time in Fig. 3. Equation (28) with n = 3 matches the exact BTC fairly well, except that it does not show the discontinuity at t = z/u = 0.1 min. The presence of unrealistic inflection points due to the third order term in the series is obvious. For the parameters of Fig. 3, p ;/ (2&+‘2 = 0.75, which constitutes a severe ,t.est of the theory. We have also compared eqn. (28) with moments for the adsorbed gas on the solid, eqns. (43) - (46). While the agreement between exact and approximate curves is fair for n < 3, the agreement is excellent for n = 3 for all values of k’z/v. From this comparison of the approximate BTC expressions with the exact Jfunction solution, we conclude that -00 rather than zero should be used as the lower limit of integration. The agreement is improved for longer columns and lower flow rates, i.e. for larger values of k’zh. For most industrial applications k’z/u is larger than 10, otherwise the separation is unlikely to be economical. Short columns with inlet pipe diameters much smaller than the column diameter will experience radial concentration gradients and radial dispersion, and will not fit the basic model of eqn. (1) at all.
4. A MODEL WITH A FIRST.ORDER TriANSFER
RATE OF
IN THE FLUID PHASE
Here we consider a model whose exact solution can easily be formulated and compared with the approximate solution. The model with a first order rate term approximately describes a variety of processes in tubular reactors with axial dispersion, e.g. a first order homogeneous reaction, a heterogeneous packed-column reactor where the rate of mass transfer from the flowing stream to the particles can be written linearly [21], the model of Suzuki and Smith [2$] for reaction and adsorption inside porous particles of a’ packed column when a pseudo steady state prevails inside the particles [ 221, and adsorption or chemical reaction on the surface of non-porous particles in the column under
conditions that a pseudo steady state applies at the surface [22]. Because of the loss of the reacting species, moments may not be easily measured by pulse response experiments, and must be estimated for a definite reactor model such as the simple one proposed hem. In dimensionless variables we have 1 a2c ac _._ _--_-__
Pe az2
ac
az
aT
-Kc=0
(47)
with initial condition,
C(T = 0,Z)
boundary
condition,
C(T, 2 = 0) = U(T)
boundary
condition,
C(T, 2 + -)
The exact solution
= 0
(43)
is finite (50)
is
C(T,Z)=+(,exp[y\l-(l+4$2\] 2 --~.
x erfc
(49)
X
T(l
+ 4~ /Pe)1’2
+
(4T/Pe)l12 +exp[iZPe\l+
X
2 +
erfc
(1 +k)r2]X
T(l +
4~ /Pe)1’2
(4T/Pe)1/2-
(51)
which has been obtained by Laplace transformation by Cho [23] to describe transport ’ processes in soil. The approximate solution (n < 3) to this model can be obtained by substituting the following expressions of the moments into eqn. (28): Ms(Z)=exp[i-ZPe]l--(l+E)Liql
(52)
z (53)
M; @) = (1 + 4K /pe)l/2
(54)
M3(Z)=Pe
-b/2
122
1+; (
1
(55)
217
This gives C(T,Z)=,-iexp[iZPe/l--(I
2 - T(l + 4s /Pe)lj2
+z)“2/]
L
Pe
\~
The second term inside the brackets is of the order of the third moment. Next we compare the two solutions, eqns. (51) and (56), and assess the accuracy of the approximate solution. For the limiting case when Pe -+ 00, i.e. when axial dispersion has no significant effect on the concentration distribution in the column, the two solutions become identical and equal to fexp( C(T, 2) =
+
--KZ)U(T---2)
wz
(57)
T=Z
eXp(-KZ)
For moderate values of Pe and 2, the graphs corresponding to exact and approximate solutions are plotted in Figs. 4 - 6. The comparison of the curves indicates that for high values of Pe and 2 the approximation is quite good. In Fig. 4 the exact expression for the BTC is plotted for 2 = 10 and K = 0.1 with a range of Peclet numbers. Decreasing Pe, which implies more axial dispersion, flattens the curves and rounds the sharp comer of the curve corresponding to Pe + 00. Only the approximate solution at Pe = 1 is presented (as dots) in Fig. 4, since it is the least accurate and is in error by only 3% at most. The error decreases for large values of Pe, but for smaller values of 2, i.e. z < lOd, the error would be greater.
(56) Figure 5 demonstrates the shape of the BTC in different locations of the bed for Pe = 4OandK = 0.1, and shows C(T + =, 2) decreasing with increasing 2. The approximate solution, shown by solid dots, is in good agreement with the exact solution even when the column is short, 2 = 1. The broken lines represent the exact BTCs when D = 0, eqn. (57). The locus of comers is where T = 2 and thus from eqn. (57) is given by exp(-O.lT). When both Pe and 2 are small, as in Fig. 6, the approximate solution obviously deviates from the exact. The approximate solution including the third moment has three inflection points, while the solution excluding the third moment has one inflection point. Although neither approximation fits the exact curve very well, the case represented by these values of Pe and 2 is hardly representative of an industrial process. High values of Pe and 2 are equivalent to high values of the parameter ~;/(2p~)l’~ = M;(2M2)1’2 since in this model
Pi (2/J2)1’2
= J- (PeZ)l12 2
4K (1 + --)
“* (58)
Large values of Pe and 2 also imply that the third moment is negligible. For packed columns, gases and liquids usually have values of Pe ranging from 20 to 40; see for example Vermeulen [ 21. Also, we usually have 2 % 1 since Z = z/d.
0.4
s ”
0.2
-Exact l
0
Approximate
0.1
0 6
a
10 T=tv/d
Fig. 4. BTCk for 2 = z/d = 10 and
‘*
14
16
K = 0.1 and different values of Pe = vd/D. Exact, eqn. (51 .), and approximate .tions, eqn. (56). nearly coincide except for Pe = 1 where approximate values are shown as dots.
218
+_-__ and the coupling
boundary
condition
(6%
= k!f (C -- CiR)
The moments for a 6 -function beginning at t = 0 are
P
‘, ;:
inlet condition
p; = ;(1+6,) 0
12
3
4
5
6
18
9
1
I
I
10
11
12
and
T=tvld
22
Fig. 5. BTC!s for various positions in the column with Pe = 40 and K = 0.1. The approximate solution is represented as dots and eqn. (57) is shown as broken lines.
P2
=--y
6 1 +
D(l + 6,)’ ;s
I
where (65) and
x(&+&)1 -3
-2
I
I
I
I
I
I
I
-1
0
1
2
3
4
5
6
T=tr/d
Fig. 6. BTCs for a short column, z/d = 1, with Pe = 1 and K = 0.1: -, exact solution, eqn. (51); - - -, approximate solution, eqn. (56) when R < 3 (/.I3 = 0); ------ 9approximate solution, eqn. (10) when n = 3.
5. A MODEL DESCRIBING DETAILS OF ADSORPTION INSIDE PARTICLES
The model developed by Kubin [6] and Kucera [7] and used extensively by Smith and coworkers [4,8] takes into account the effects of (a) axial dispersion, (b) fluid to particle mass transfer, (c) intraparticle diffusion and (d) reversible adsorption. The chief restriction is that the transport and adsorption processes be linear. The governing differential equations are D~~--+~-
0
hi
pP
acads
P
at
- -at -__
-
(59) =0
(60)
(61)
(66)
Ammons et al. [ 10,111 obtained rate parameters for adsorption of methyl mercuric chloride on activated carbon by the moment technique. Because the capacity of the carbon for the adsorbate was so large (large value of KX), pulse experiments were not feasible. Ammons carried out the experiment using a step function input and related the calculations of the moments to the analysis of the experimental BTCs. Values of the first moment p; and the second central moment p2 are presented by Ammons [ll] and may be used in the approximate theory to calculate BTCs. Since Ammons did not calculate pa, we use eqns. (27) and (28) with n < 3 to plot the BTCs of Figs. 7 and 8. For the BTCs in Fig. 7, the ratio pi/ is 2.66 for BTC 4 and 1.88 for BTC 7, (2P2 P2 and the agreement between the approximate theory, eqn. (28), and the experimental data is good, with better agreement for larger values of the ratio. In Fig. 8, P;/(@~)~‘~ equals 1.11. Since this ratio is substantially less than 2.5, the agreement is less satisfactory. In this figure we have also plotted eqn. (27) with n < 3, the solution obtained when the lower limit of integration is zero. Although this curve is a better approximation for small times, it is poorer for large times.
219
0
tx10-3,
SK
Fig. 7. Comparison of the approximate expression for BTC (dots), eqn. (28) for n < 3, with experimental results (line), BTCs 4 and 7 of Ammons [ 111.
processes can be combined without the aid of a detailed theory such as that described by eqns. (59) - (62). The overall transfer rate constant kOp represents the combined resistance and can be used in a mathematical model such as the Goldstein J function problem [l] . To assess the accuracy of this technique, we examine expressions of moments that can be obtained from the combination of resistances in the Goldstein J function problem and compare them with results of the detailed model based on equations describing the actual transport phenomena, i.e. eqns. (59) - (62). According to Vermeulen 123 we make the following substitutions in eqns. (30) and (31): FE’= ko, Aa&
(67)
and k = kopap
(68)
where A = I!&,(1
- Ly)
(69)
The overall transfer constant is expressed by 0 0
4
8
12
16
20
24
t x 10-3, SC
Fig. 8. Comparison of the approximate expressions (n < 3) eqn. (27) (-------) and eqn. (28) (- - -) with the experimental BTC 10 (-) of Ammons [ 111.
Another theory for the BTC which involves a numerical integration has been presented [3 - 51. This approach is based on the approximations D = 0 and flaci/a t = 0 in eqns. (59) and (60), namely that axial dispersion and volume capacity of the porous particles are negligible. Neglecting the volume capacity alters the first moment and shifts the BTC to smaller values of t; the shift is less than 6% for the experiments discussed by Schneider and Smith [4] . Neglecting axial dispersion reduces the second moment by as much as 70% for their experiments and thereby sharpens the BTCs rather unrealistically. Since for the conditions of experiment under consideration the approximations in our theory have much less effect than these discrepancies, we forego a detailed comparison (see ref. 22 for details). Vermeulen [2] suggests that mass transfer resistances involved in packedcolumn
1 1 _=_-++._+ kop kcp
1
a, -
kcf
k
(70)
where (71)
k,, = kp + k,,,lA and K = k:,,/K:,
(72)
When axial dispersion is included in the suggested manner 1 -= kcf
1 -+-
1
(73)
k:,
kf
in terms of kb = au 2/oap
(74)
we obtain k,
=
K:p,(l
---a) --a
+
+
-1
=XP~
3kf
5 + R2KXp, k’ads 15Di
+
KXpp(l --a) -$
t
(75)
From eqns. (40) and (41) the moments for the Goldstein J function problem become /.& =f
1 --a 1 + -.-a-
(76)
220
and
tb =p;
z 1 --a PJ;~ cl2 = -.v
R2Ki2p;
_+---+
-a
(
R2p;KX2 +--+3kfR
p;KX2
aD
CY2
u2
(77)
For comparing these estimated moments with the calculated moments according to the model (eqns. (59) - (62)) we note that p,Ki/P is the critical parameter. If pp K;//? S 1 then the expressions for moments are identical; however, if pP KX/fl Q 1 the expressions are entirely different. Since the former inequality usually holds in adsorption processes, the combinations of resistances technique must be considered quite accurate. To illustrate this, consider the adsorption of ethane from helium onto silica gel with the parameters measured by Schneider and Smith [4]. For c(\, eqn. (76) gives 84.61 s while eqn. (63) gives 86.69 s, a difference of 2.4%. For p2, eqn. (77) gives 414.1 s2 while eqn. (64) gives 445.3 s2, a difference of 7%. In this example KX = 14.6 ml g-’ silica gel so that ppKXI p = 33.94 which is sufficiently greater than unity to give excellent agreement between the estimated and correct values of the moments. For adsorption of propane and n-butane the equilibrium adsorption constants are larger than for ethane and the agreement is even better. The third moments, however, are poorly estimated, and thus the estimation procedure will not accurately predict the skewness of pulses or BTCs. However, since for most industrial processes pP Ka/fla 1, we may conclude that the estimation method is quite accurate for the analysis and design of adsorption in packed columns.
6. THE BREAKPOINT
AND OPTIMIZATION
(78)
in terms of a constant a that is related to cb by
1503
&ads
-adp2
OF
COLUMN OPERATION
A most important design parameter in the adsorption process is the time t,,, corresponding to the breakpoint, that fixes how long the column should be in operation. The time t,, occurs when an arbitrarily small but signif* icant effluent concentration cb appear%at the end of the column. Since the first moment locates the center of the BTC and the second moment the width, for a symmetric BTC the breakpoint can be defined as
(79) which follows from eqn. (28), it < 3. Thus if cb is known, a can be calculated from eqn. (79) and hence tb from eqn. (78). To maximize the utilization of a column before it is regenerated, we wish to maximize tb. From eqn. (78), we see that the parameters that increase P; and/or decrease p2 will increase tb. Knowing how these moments depend on the parameters, e.g. as in eqns. (63) - (66), will enable us to optimize tb. By comparing the effect of the rate parameters on p; and p2, we can conclude that increasing Q, k&s, kf and decreasing R2 and D will increase tb. To determine the optimum velocity, or rate of throughput, we express the net income of an adsorption column of given length as the difference between the value of the product and the pumping cost. We assume that other costs are independent of the velocity. The value of the product, 7 tbuo”d2/4, is given in terms of the price y per unit volume of product and the throughput after time tb. The pumping cost WV; follows from the equation for the pressure drop in a packed column. The net income function is then
(80) which is to be maximized with respect to v = V&J. As an illustration, we use the expressions for moments, eqns. (63) - (66), along with the parameters for the adsorption of ethane from helium onto silica gel given by Schneider and Smith [4]. To give a more realistic commercial meaning to the product price we considered a column cross section 100 times larger than that used in the laboratory apparatus. Values of o and cuy are chosen as 21.0 $ crnm2 s2 and 0.05 $ cm-‘, respectively . The net income function f for various values of a, along with the pumping cost, are plotted in Fig. 9. The optimum velocity increases for increasing a, i.e. for increasing purity, while the net income decreases. The positions of the optimum velocity for differ-
221
(see also ref. 24). The moment expressions and present theory constitute a foundation for optimizing BTCs for column processes.
NOMENCLATURE i
a UP
0
2
4
6
0
10
12
14
v, cm/sac
Fig. 9. Plot of income us. velocity for various purities of effluent. The broken line shows the values of optimum velocity.
ent values of a are shown by the broken line. For values of a 2 11.5 there would be a negative net income owing to the high purity required of cb. It is interesting that the range of optimum velocity is less than 5 cm s-l in this case of adsorption of ethane on silica gel.
c(t, 2) cad&, r) c,(t, z= 0) q(t = 0,z) Ci(t,
CO c
7. CONCLUDING REMARKS
A supposition of the work described herein is that exact analytical, or more usually numerical, solutions for BTCs provide in many cases more information than is needed for the design and optimization of industrial column processes. For such processes governed by linear rate expressions, and therefore linear differential equations, the moment method, combined with the Duhamel theorem, presents a solution that is relatively easily obtained and simple in form (an algebraic expression). For many industrial columns, where @C(;/(~P#~ > 2.5 is generally satisfied, only terms up to the second central moment are needed to represent the BTC accurately. The first moment II> locates the center of the BTC and the second central moment p2 determines the spread of the BTC. For shorter columns with larger velocities, the third moment can be computed and this will extend the accuracy of the BTC representation by including asymmetries. For example, Kubin [ 61 has presented expressions for the third central moment for the adsorption model used by Schneider and Smith [4]
6
D Di d
f &I
k, k’
r, 2)
coefficient in expression for tb, eqn. (78) surface area of 30 -0)/R, spherical sorbent particles per unit volume of bed, 1 cm-’ concentration in interparticle space, mol cmm3 concentration of adsorbed species based on unit void volume of particles, mol cmm3 concentration input, mol cm-l initial concentration profile in interparticle space, mol cme3 concentration in intraparticle pore space, mol cmm3 characteristic constant concentration, mol cmW3 c/co, dimensionless form of preceding concentrations dimensionless concentration response to a unit step function dimensionless concentration response to a Dirac delta function input Laplace transform of C effective axial dispersion coefficient cm2 s-l effeciive intraparticle diffusion coefficient, cm2 s-l diameter of cylindrical column, cm net income, $ Hermite polynomial overall transfer rate constants in Goldstein J function problem, S-l
k’ads
adsorption
rate constant,
cm3
g-ls-1 kd kP
k POW k
desorption rate constant, 1 s-l particle phase diffusion rate constant, 1 s-l rate constant of pore diffusion, Is-l mass transfer coefficient at sphere surface, cm s-l k:,,/k,, adsorption equilibrium
222
L m m,(z)
M,(Z) M;(Z) M,(Z)
Pe r R 63
61’
tb
T u UC i
constant, cm3 g-l length of column c*/cI, equilibrium constant c(t,z) t “dt, nth moment of J 0 concentration m&)/co, zeroth dimensionless moment p; u/L, first dimensionless moment ~12 u2/L2, second central moment, dimensionless ud/D, Peclet number for axial dispersion radial coordinate in the spherical particles, cm radius of the spherical particles, cm di(c, CiR), term in eqn. (1) representing mass transfer from the interparticle fluid to the surface of particles 61’(Cads, Ci), telXl in eqn. (2) representing chemical reaction or adsorption within the pores of particles Laplace transform parameter, 1 s- ’ dimensionless Laplace transform parameter time dorresponding to the breakpoint, s td/u, dimensionless time interstitial velocity of fluid, cm s-l cw, superficial velocity, cm s-l length coordinate of the bed, cm z/d, dimensionless length coordinate
Greek symbols a inter-particle void fraction of the packed column porosity, intraparticle void fraction 0 I.4
1 m,
-
s ct dt, first absolute moment 0 0
‘0 Jc(t
moo moment
- p;)ndt,
nth central
t P PP
r/R, dimensionless radial coordinate fluid density, g cm- 3 apparent particle density, g cm- 3
REFERENCES 1 S. Goldstein, Proc. R. Sot. London, Ser. A, 219 (1953) 151. 2 T. Vermeulen, G. Klein and N. K. Hiester, in R. H. Perry and C. H. Chilton (eds.), Chemical Engineer’s Handbook, 5th edn., McGraw-Hill, New York, 1973, Sect. 16. 3 S. Masamune and J. M. Smith,A.Z.Ch.E. J., 11 (1965) 34. 4 P. Schneider and J. M. Smith,A.Z.Ch.E. J., 14 (1968) 762. 5 J. B. Rosen, J. Chem. Phys., 20 (1952) 387. 6 M. Kubin, Coil. Czech. Chem. Commun., 30 (1964) 1104, 2900. 7 E. Kucera, J. Chromotogr., 19 (1965) 237. 8 T. Furusawa, M. Suzuki and J. M. Smith, Catal. Rev., 13 (1976) 43. 9 K. Rade ke, K. Wiedemann and D. Gelbin, Chem. Tech. (Leipzig), 28 (1976) 476. 10 R. D. Ammons, N. A. Dougharty and J. M. Smith, Znd. Eng. Chem. Fundam., 16 (1977) 263. 11 R. D. Ammons, Adsorption of methyl mercuric chloride on activated carbon, M.S. Thesis, Univ. Calif., Davis, 1975. 12 R. G. Carbonell and B. J. McCoy, Chem. Eng. J., 9 (1975) 115. 13 J. K. Grd, Families of Frequency Distributions, Hafner Publ. Co., New York, 1972, p. 28. 14 J. C. Giddings, Dynamics of Chromatography, Part I: Principles and Theory, Marcel Dekker, New York, 1965. 15 0. A. Hougen and K. M. Watson, Chemical Process Principles, Wiley, New York, 1946 ; p. 1080. 16 A. Klinkenberg, Znd. Eng. Chem., 46 (1954) 2285. 17 T. E. W. Schumann, J. Fmnklin Inst., 208 (1929) 405. 18 C. C. Fumas, Tmns. A.Z.Ch.E., 29 (1930) 142. 19 S. R. Brinkley and R. F. Brinkley, Math: Tables Aids Comp., 2 (1947) 221. 20 S. R. Brinkley, Math: Tables Aids Camp., 6 (1952) 40. 21 M. Suzuki and J. M. Smith, Chem. Eng. Sci., 26 (1971) 221. 22 M.-S. Razavi, A method for predicting breakthrough curves in packed column processes, M.S. Thesis, Univ. Calif., Davis, 1978. 23 C. M. Cho, Can. J. Soil&i., 51 (1971) 339. 24 B. J. McCoy and R. G. Carbonell, A.Z.Ch.E. J., 24 (1978) 159.