123
Ind. Eng. Chem. Fundam. 1982, 2 1 , 123-131 io =
Nomenclature
dimensionless local rate of energy absorption cA = molar concentration of the reactant cA0= initial molar concentration of the reactant 3 = effective diffusivity e’’’ = local rate of energy absorption e”’,, = r@cAJ0= characteristic scale for e”’ & = rate of energy absorption within the slab E = energy absorbed within the slab I = intensity of the radiation Io = a/r = intensity of the diffusive emission K = local extinction coefficient L = geometric width of the slab n = reaction order R = hemispherical reflectivity t = time T = hemispherical transmissivity x = geometric coordinate Greek Letters a = power emitted per unit surface by the boundary x = 0 6 = molar extinction coefficient y = cA/.cAO = dimensionless reactant concentration yav= dimensionless average reactant concentration 7 = *~?cAO’-~IO/L= dimensionless time p = come of the angle between the direction of the radiation intensity and the positive x axis - x L = dimensionless geometric coordinate pD = reflectivity of the boundary x = L pD = diffuse component of p p s = specular component of p T = optical coordinate 7 = optical coordinate at t = 0 (itu =
f p$ .+
optical thickness
= optical thickness at t = 0 IC. = a ) / r L c ~ ~ ” - =~ 1dimensionless ~4 number ?o
4 = quantum yield L i t e r a t u r e Cited
Bandlnl, E.; Spamlgloll, C.; Santarelll, F. Chem. Eng. Sci. 1977, 32, 89. Beach, H. L.; Ozlslk, M. N.; Sleweft, C. E. Int. J . Heat Mass Transfer 1971, 14, 1551. Bhattacharya, A.; Deshpande, P. K. Chem. Eng. Sci. 1979, 34, 145. Boffi, V. C.; Santarelll, F.; Stramlgloll, C., J . Quant. Spectrosc. Radiat. Transfer 1977, 18, 189. Cerda. J.; Irazoqul, H. A.; Cassano, A. E. AIChEJ. 1973, 19, 963. Chandrasehkar, S. “Radlatlve Transfer”; Oxford University Press: London, 1950. Felder, R. M.; Hill, F. B. Chem. Eng. Sci. 1989, 24, 385. Irazoqui, H. A.; Cerda, J.; Cassano, A. E. AIChEJ. 1973, 19, 460. &lagelll,F.; Santarelll, F. Chem. Eng. Sci. 1978, 33, 611. Ozlslk, M. N. ”Radiative Transfer”; Why: New York, 1973. Pasquall, G.; Santarelll, F. Chem. Eng. Commun. 1978, 2, 271. Roger, M.; Vlllermux, J. Chem. Eng. J . 1979, 17, 219. Shendalmn, L. H.; HIII, F. 8. Chem. Eng. J . 1971, 2 , 261. Schechter, R. S.; Wlssler, E. H. Appi. Sci. Res., 1980, A9, 334. Spadonl, G.; Stramigloll, C.; Santarelll, F. Chem. Eng. Sci. 198Oa, 35, 925. Spadonl, G.; Stramlgloll, C.; Santarelll, F. Chem. f n g . Commun. 1980b, 4 , 643. Splga, G.; Santarelli, F.; Stramigloll, C. Int. J . Heat Mass Transfer, 1980, 2 3 , 841. Stramigioll, C.; Santarelll, F.; Foraboschi, F. P. Appl. Sci. Res. Feb 1977, 33, 23. Stramlgbll, C.; Spadoni, G.; Santarelli, F. “The Role of a Reflecting Boundary In Improvlng the Operation of an Annular Photoreactor”; submitted for publlcatlon In Ind. Eng Chem, Fundam.
.
Received for review December 30, 1980 Accepted October 26, 1981
This work was supported by C.N.R.-Roma (Grant No. 79.01273.03/115.3013) and by Minister0 P.1.-Roma.
Simulation of Multicomponent Adsorption Beds. Model Analysis and Numerical Solution Masslmo Morbldelll, Albert0 Servlda, Gluseppe Storll, and Serglo Carrh Istituto di Chimica Fisica, Elettrochimica e Metallurgia, Politecnico di Milano, Piazza Leonard0 da Vinci 32, 20 133 Milano, Italy
The dynamics of multicomponent adsorption performance in a fixed bed adsorber is examined numerically. For this purpose two different models of inter- and intraparticle mass transfer are considered. The particle internal poroslty may be either neglected or taken into account. The model includes axial dispersion, variable linear velocity, and local equilibrium, as well as nonequilibrium, between the pore fluid and the solid phases. A generalized nonlinear multicomponent adsorption isotherm is also included. A new numerical procedure, based upon the Lax-Wendroff method, is presented. I t is applied to simulate the dynamics of fixed bed adsorbers in response to a step or a pulse change in feed composition. The influence of the operating conditions and the relevant physicochemical parameters on the efficiency of the separation processes through fixed bed adsorbers are numerically investigated.
Introduction
Multicomponent adsorption is often encountered in industrial applications to separate or purify liquid and gaseous mixtures. A mathematical model of fixed bed adsorbers is very useful to obtain a proper design and to define the optimal operating conditions for this equipment. Different models for the adsorption process of gaseous or liquid components have been proposed in the literature. There also exist various numerical techniques developed for the solution of the fixed bed adsorber models. Analytical solutions, on the other hand, have been obtained only for the single-solute case with linear equilib0198-4313/82/1021-0123$01.25/0
rium at the surface (Lapidus and Amundson, 1952;Rosen, 1952) or with Langmuir adsorption kinetics (Thomas, 1944). However, for multicomponent fixed bed adsorbers the solution can be obtained only by numerical means. Cooney and Lighfoot (1965)proposed an approximate solution of this problem assuming the presence of a constant pattern along the unit. Based on this assumption, the original system of partial differential equations (PDE) is reduced to a system of ordinary differential equations (ODE). However, this hypothesis eventually turns out to be quite restrictive, since it requires constant feed composition and space velocity along the adsorber, favorable 0 1982 American Chemical Society
124
Ind. Eng. Chem. Fundam., Vol. 21, No. 2, 1982
equilibrium between the gaseous and the adsorbed phases, and finally long adsorber beds. Von Rosenberg et al. (1977) developed a different numerical procedure based on the centered difference method as a solution starter, and on a backward difference equation later in the process. This procedure is very efficient for single-solute problems, while it appears to be more complex for multicomponent mixtures. In fact in this case it requires the solution of an algebraic system of nonlinear equations at each point of the finite difference grid. This method does not contemplate axial dispersion and variation in space velocity. By adopting this same approximation, Hsieh et al. (1977) developed a model of multicomponent liquid phase adsorption in fixed beds, taking into account both separately and in combination the effects of liquid and/or solid phase mass transfer resistance with a different kind of governing adsorption isotherm. In this work a third-order method of characteristic algorithm, previously developed by Vanier (1970) for the numerical solution of systems of semilinear hyperbolic equations, is used. None of the methods hereby mentioned takes mass dispersion into account in the axial direction. This restricts applicability of such methods to cases where step function feed composition is used. Namely in the case of pulse feed, commonly used in the elution technique, the behavior of the adsorber bed is largely affected by axial dispersion. The statistical moments method (Grubner, 1968) has been widely used for the description of gas-solid chromatography beds in which the feed always takes the pulse form. This method can on principle be used for the simulation of fixed bed adsorbers. I t does not, however, provide any real solution to the problem, but an approximate curve represented in general by a series form (for example, the Gramm-Charlier series), which shows the same statistical moments of the exact solution. Furthermore, this method, which reduces the system of PDE to a more simple system of ODE, can be applied only to linear adsorption isotherms. In a recent paper, Liapis and Rippin (1978) have applied the orthogonal collocation method for the numerical solution of the model with axial dispersion. In this model the entire concentration profile in the particle is evaluated at each axial position and at each time value. While this method appears to be more powerful than of all the previous ones, its application is very time consuming since a two-dimensional collocation procedure is required. The aim of this work is to develop a numerical procedure for the solution of the mathematical model of an isothermal fixed bed adsorber in which the effects of axial dispersion, inter- and intraphase mass transfer resistance, and variable space velocity along the unit are taken into account. While the first two effects have already been explored by other authors, the last one has not received the same amount of attention. However, when mixtures with high concentrations of the adsorbable components are fed to the adsorber and the adsorption process leads to a decrease in the global mass flow rate, the space velocity cannot be assumed to be constant. This is the case of industrial fixed bed adsorbers in which the separation of some components is achieved for example, the separation process (McKetta and Kobe, 1960) of n-paraffins from mixtures of isoparaffins and aromatics. In this process the paraffm concentrations in the feed stream are considerable and the regeneration of the bed is performed by decreasing the pressure from 60 to 1 psi so that the particle pores are left almost completely empty. When the new operation is started, a significant amount of fluid is substracted from
Figure 1. Scheme of the particle pore structure (Vermeulen et al., 1973): (1)external mass transfer; (2) pore diffusion; (3) surface adsorption or equilibrium; (4) solid diffusion.
the gaseous stream by the adsorbent bed, and then the space velocity changes along the unit axis. Particular attention has been devoted in developing the numerical procedure to the pulse-response problem. In fact, analysis appears to be more efficient than the step-response to the evaluation of the physicochemical parameters involved in the model. Interparticle and Intraparticle Mass Transport The adsorption of a multicomponent mixture on a solid particle can generally be regarded as a mass transport in which the rate expression shall include the transfer processes both inside and outside the adsorbent particle itself. Many different models have been presented in the literature. We shall here review only two of the most significant. The first, as recently used by Hsieh et al. (1977),assumes the linear driving force approximation when describing the mass transport both inside and outside the adsorbent and assumes also the equilibrium at the fluid-particle interface. The rate of adsorption of the ith component is thus given bY
where C and F are the concentrations in the fluid bulk and the solid phases, respectively, and surface values are represented by the index s. The diffusion process thus represented allows great simplification of the numerical solution of the model. The description of the competitive effects which are typical of multicomponent adsorption is left only with the adsorption isotherm. The most commonly used isotherm is of the Langmuir type
This equation, coupled with eq 1,constitutes the complete model of the adsorption process. Other types of adsorption isotherms can be used to replace eq 3 and the most general expression is one proposed by Fritz and Schliinder (1974) ui(C:)d,o
rp =
N
1
+ Cjbj(Cj”)d” 1
(3)
where di0and dijare adjustable parameters. In particular for dio = djj = 1, eq 3 reduces to the Langmuir adsorption isotherm in eq 2, and for bj = 0 to the Freundlich equation. The second approach to modeling the fluid particle system is the one applied by Liapis and Rippin (1977). In this approach, previously adopted by Kucera (1965) in gas-solid chromatography problems, the solid particle is regarded as shown in Figure 1. Three different phases
Ind. Eng. Chem. Fundam., Vol. 21, No. 2, 1982 125
are present: the fluid-bdk phase, the pore fluid phase, and the solid phase. The first two are connected by the mass transfer through the external film. Once the adsorbable molecule has reached the particle surface, it diffuses with *simultaneous adsorption through the pore liquid phase within the particle. This process can be described by a material balance of the generic adsorbable component i, in the adsorbent particle
- complete
___
model
lumped mide! ( k s ! ' lumped model [ k,)'
with the boundary conditions 0251
0
aCpi/ar = 0;
r =0
(6)
and the initial condition C, = 0 at t = 0. Even in this case the evaluation of ri,which is assumed to be in equilibrium with Cpithroughout the particle, is performed by means of an adsorption isotherm of the type reported by eq 2 or 3. Different modifications of this model can be obtained whenever the diffusion in the particle is not represented by normal diffusion in the pore liquid phase (pore-diffusion), as in eq 4, but by diffusion in the solid phase (solid-diffusion) or in the pore surface (surface-diffusion). All of these models have been reported in detail by Liapis and Rippin (1977). The model given by eq 4-6 is obviously more accurate than the previous one (eq 1);however, the solution requires greater numerical effort. That is, at each time value a second-order differential equation must be solved with separate boundary conditions, instead of an algebraic equation. In order to simplify subsequent calculations, a lumping procedure can be applied to eq 4. Integrating both sides of this equation in the particle volume, we obtain
where Cpi and Ti are the volume average values of the concentration in the pore liquid and in the solid phase, respectively. The first term on the right-hand side of eq 7 is obtained from the Laplacian term in eq 4 by using the boundary condition 5. This term can be expressed in terms of the average composition Cpi, by means of the linear driving force approximation kfi[Ci - Cpi(Rp)I = kli(Ci - Cpi)
(8)
where kli is the global mass transfer coefficient, which can be evaluated from the external and internal mass transfer coefficients, kfi and ksi, respectively, as 1 -1= - +1(9) h i
Kfi
ksiep
Two problems arise from this procedure. The first is the evaluation of the term Ti of eq 7. The approximation generally used is
where f(C,,J indicates the adsorption isotherm represented by the right-hand side of eq 2 or eq 3. Equation 10 is obviously correct for linear adsorption isotherms, but it involves some approximation for nonlinear isotherms. The second difficulty lies with the evaluation of the internal
1
I
20
I
I
LO
I
I
I
I
60 90 Real time I . [mini
I
I
EO
I
!
;
Figure 2. Time evolution of the bulk concentration of an adsorbable component in a finite bath (data from Table 11; component: 2; t = 0.8; Co = 0.408 M.
mass transfer coefficient, kSi. It is defined, similarly as kfi in eq 5, by the relationship
The value of ksi is then time dependent and can be evaluated exactly only from eq 11. This implies knowledge of Cpi(R ) and dCPi/drlrSRP and therefore the solution of the comprete problem (eq 4-6). Normally a time average value of kSi is used, and in particular the one proposed by Glueckauf (1955) Dpi
k,i = 5-
RP
In order to assess the accuracy of the adopted lumping procedure, we consider the problem of a finite bath of the adsorbable components in which a certain amount of these components is introduced. The material balance of the ith component leads to the equation
where t is the void fraction of the bath and Ci is the concentration in the fluid phase. Equation 13 can be integrated coupled either with eq 4-6 or with eq 7 and 8. For illustrative purposes Figure 2 shows the external fluid concentration vs. time curve, through both the complete and the lumped model with two different constant values of the mass transfer coefficient: (kJ1 and (Ksi)". The former has been evaluated by eq 12 and the latter through a fitting procedure of the complete model curve. These average values have been compared in Figure 3 with the exact Itsi evaluated at each time by means of eq 11. As a matter of fact, final conclusions about the comparison of the two models cannot be drawn, since such a comparison depends on the adsorption isotherm type as well as on the mass transfer coefficients, kfi and the equilibrium ratios, Ki. However, the two external concentration vs. time profiles in Figure 2 with &,.)I and (kBi)" appear to be comparable. It is worthwhile to note that (kSi)Ihas been evaluated by eq 12 which has been derived by Glueckauf for the pure diffusion problem, i.e., without reaction, while here we assume equilibrium between the solid and the pore liquid phases, that is infinite adsorption rate. This result is rather surprising because the rate of reaction is known
~
Ind. Eng. Chem. Fundam., Vol. 21, No. 2, 1982
126
1.0
where the superscript indicates the feed stream conditions, and the (3N) initial conditions
I\ 1
0
.
ci = C?(z); cpi=
t =0
t =0
C,t(z);
z
E (0,l)
ri = r:(z);
1
I
20
~
I
I
60
10
I
80
I
1
100
I
1
1
113
Real time I Lminl
Figure 3. Plot of k, vs. time: (-) k, evaluated by eq 11 for the example of Figure 2.
to play an important role on the value of ke Namely, the presence of the consumption of the component in the particle, due either to a chemical reaction or to an adsorption reaction process, increases the concentration gradient at the particle surface and thus the value of kSi (Figure 3) through eq 11. The preceding analysis points out the approximations involved in the lumping procedure. Application of eq 12 appears justified, despite the significant difference between the value of the effective mass transfer coefficients (kG)*and (k,)". This lumping procedure, which is widely used in process design calculations, has been adopted in our work due to the significant simplifications in the model equations. In particular, the fluid particles mass transfer will be described through eq 1 or eq 7 and 8, coupled with the adsorption isotherm 2 for the surface equilibrium. While it is established that adsorption on the surface is seldom rate limiting, the nonequilibrium model has also been considered as a generalization. In it a kinetic adsorption is assumed between the fluid and the solid phases and the equilibrium relationship 2 is then replaced in the particle model by the equation
where k& is the adsorption-rate constant of component i. An approach of this type, although if a linear adsorption isotherm was employed, has been recently used by Hattori et al. (1980) to describe the adsorption of aromatic hydrocarbons on zeolites.
Modeling of the Fluid-Phase Flow In a Fixed Bed The model developed above describes the adsorption of solutes from a flowing fluid stream onto a fixed bed of particles. The material balance for the ith component in the fluid phase, including axial dispersion and variable space velocity along the bed axis, is
with the (2N + 1) boundary conditions
aci -- 0 ; _ az
u=uF;
z=L,t>O z=o
(16)
t =o (17) The index 0 indicates the saturation conditions for the bed at the beginning of the process. In the formulation of eq 15 the following assumptions have been adopted: (1)isothermal conditions, (2) no radial gradient, (3) a constant axial dispersion coefficient, and (4) spherical and uniformly sized adsorbent particles. Assuming N components in the mixture, there are 3N + 1 unknowns in the model; the external fluid phase Ci(z,t), the pore fluid phase Cpi(z,t), and the solid phase concentration Ti(z,t) for each component, and the space velocity u(z,t). The corresponding equations are, apart from the N equations 15, the N material balances on the particle (eq 1or eq 7 and 8) and the N surface equilibrium relationships (eq 2 or 3). If local equilibrium between solid and fluid phases at each point in the pore is not assumed, the latter equations are replaced by eq 14. The last equation, necessary for the evaluation of the unknown u, is-a stoichiometric equation N
Pm =
CjCj
(18)
1
where pm is the mixture molar density, which can be evaluated through the Amagat equation
where pj, is the molar density of the pure jth component. Combining eq 18 and 19 gives N
CiCi/Pi = 1
(20)
1
and with this equation the mathematical model is complete. It is thus constituted by a parabolic system of 3N 1 partial differential equations. Numerical Solution The numerical procedure developed for the solution of the mathematical model derives from the method previously applied by Lax and Wendroff for conservation-type equations (Amundson and Aris, 1973). The application of the procedure will be developed in detail for the case of no local equilibrium between solid and fluid phases (eq 14); use will be made of the lumped pore diffusion particle model, given by eq 7 and 8. The model equations in dimensionless form are
+
Ind. Eng. Chem. Fundam., Vol. 21, No. 2, 1982
$IO8
with the boundary conditions
a$i
-=o; as
$.Os I ?
pi
pi'';
(6s)2
(34)
(25)
s=o
and the initial conditions $i
+ 6s) - 2$i(O,s) + qi(6,s - 6s)
where 6s is the discretization interval along the bed axis. The second derivatives with respect to 6 are obtained by differentiation of eq 21 and eq 23 and 24. From eq 23
s=1,6>0
u=l;
$i(O,s
=
127
and eq 24
yi = y;O;
s
E (O,l), 6 > 0
(26)
where $i =
ci/co;pi = Cpi/co; yi = iyrm;
= u/uF; s = z / L ; 7 = L +F; 6 = t/7; Pe = uFdp/(DLe); pi* = pi/C0; Xi = klir(3/Rp); = [rmp,(i- ep)/col;K ~ = * cai (27) bi = kci7co; u
In the definition of the dimensionless concentrations, the total fluid phase has been used in order to obtain 0 I $i I 1 and 0 I pi I 1, and to enable the model to describe both adsorption and desorption operations. If local equilibrium is assumed, iz, and eq 24 is replaced by eq 3. The detailed description of the modification which occurs in the numerical procedure to solve the model in this case is reported in the Appendix. The evaluation of the concentration profiles at each time is performed by means of Taylor expansions truncated after the second term
-
The numerical procedure herewith developed consists of the explicit evaluation of all the derivatives involved in eq 20-30. Using eq 22 the concentration of the Nth component can be evaluated from N-1
$N = P N * ( ~-
C j$j/pj*) 1
(31)
Substituting eq 31 in eq 21 with i = N and then adding all eq 21 with i = 1, 2, ..., N - 1 multiplied by l/pj* it follows
which will be used to evaluate the dimensionless space velocity along the bed. The derivatives in eq 28-30 will now be evaluated at each value of s referring to the system constituted by: eq 23 and 24 with i = 1, N; eq 21 with i = 1, N - 1; eq 31 and 32. The derivatives a$i/a6 for i = 1,N - 1 are evaluated from eq 21, dpi/a6 and dyi/a6 for i = 1, N from eq 23 and eq 24, respectively. The derivatives with respect to s at a given value 0 are estimated by means of a centered difference method. For example
+ 6s) - +i(e,s - 6s) 26s
(33)
Finally from eq 21
(37)
where the space velocity variation with time has been neglected, and &/as is evaluated by eq 32. The terms a3$,/a6as2and a2$i/a6ds require further manipulations, however. Differentiating eq 21 with respect to s
where all the terms can be estimated at each 6 value through the above-mentioned discretization procedure along s. A similar procedure can be applied to the calculation of the term d3$i/a6ds2 by further differentiation of eq 38 with respect to s. In brief, given a concentration profile $ ,: p/', and y/' at 6 = 0, the new profiles at 6 = 66 are evaluated through eq 28-30, in which the derivatives are calculated as above outlined. It may be well worth pointing out that the value of $N is not calculated through the Taylor series expansion in eq 28, but directly from eq 31. Once the new concentrations profiles are found, the new space velocity profile along s is also evaluated using a Taylor series with respect to s truncated after the second term
v(6,s
+ 6s) = u(6,s) + el as
o,s
6s
+
$Io,#y (39)
where both the derivatives at (0,s) can be estimated by means of eq 32. By iterating the alternate use of eq 28-30 and eq 39 it is possible to build up the time evolution of the concentration and space velocity profiles along the adsorber. On the whole, this procedure consists of a centered discretization along the s axis and a marching integration in time. Convergence Characteristic The proposed method is explicit in the time variable and therefore, as with all methods of this type, it exhibits instabilities for particular values of the integration intervals 6s and 66. An analysis of the convergence characteristics has been performed by truncating the Taylor series (eq 28-30) after the first and the third term. In the latter case, the third derivatives of the concentrations yi, $it and vi
128
Ind. Eng. Chem. Fundam., Vol. 21, No. 2, 1982
Table I. Comparison of the Numerical Solutions Obtained through Taylor Expansion (Eq 28-30) Truncated at the Second and Third Term (Data from Table 11) third order
second order -
time, e
,j11.1 g 1F
1.785 2.379 2.974 3.569 4.164 4.7 59 5.354 5.948 6.543 7.138
3.5321 X 2.0070 X lo-* 9.2998 X 10.' 3.7053 X 10-I 1.1043 1.6465 1.6109 1.4921 1.3377 1.1837
+ J+zF
+INIF 3.8528 X 1.4268 X 5.2153 X 1.8422 X 5.7306 X 1.2903 X 2.3175 X 3.8390 X 5.7959 x 7.7276 X
3.5305 X 2.0063 X lo-' 9.2964 X lo-* 3.7043 X lo-' 1.1042 1.6467 1.6110 1.4921 1.3371 1.1837
lo-'
10.' 10.' 10.' 10.' 10-1 lo-'
with respect to 0, can still be evaluated explicitly by further differentiation of eq 35-37. The three solutions thus obtained for an illustrative example have been compared. The first-order approximation is unacceptable because of the small values of 60 that are required to overcome instabilities. the second- and third-order solutions are reported in Table I. The reported differences appear to be of little significance, and because of the greater computation time required by the third-order approximation, the second order has been selected. No rigorous relationship between 6s and 60 can be given to ensure stability of the method. For single solute cases in which no axial dispersion is present (DL= 0), the parabolic system becomes hyperbolic with characteristic equations dO/ds = l/v;
s = constant
(40)
From the theory of hyperbolic first-order partial differential equations, it follows that a necessary condition for stability is 68/6s 5 l / v
(41)
when DL # 0, we can obtain the relationship L 68 I(6s)2PedP
by analogy with such equations
aY dY - +aat
ax
a2Y
- b-
ax2
+cy = 0
(43)
which are stable for 6t I ( 6 ~ ) ~ / [+2~ b( S X ) ~whenever ] we use an explicit finite-difference scheme similar to the one developed here with Taylor series expansion (eq 28-30) truncated at the first term. Although the criteria given by eq 41 and 43 are not rigorous, they can be very useful to guide the choice of the values 60 and 6s to perform calculations. With a suitable selection of these values the method appeared to be efficient (it requires about 1 min of CPU on a UNIVAC 1108 to calculate three-component adsorption on a 0.5-m fixed-bed for 2 h of operation). The accuracy and efficiency of the method has also been compared with other solution techniques previously published in the literature. In particular, the problem reported by Hsieh et al. (1977) in their Table 6 has been solved. In this case the fluid-particle mass transfer model used is the one given by eq 1 assuming solid diffusion as the determining step and local equilibrium at the surface. The solution obtained along with the ones given by Hsieh et al. (1977) and Cooney and Lighfoot (1965) are shown in Figure 4. The noticeable differences in the figure are caused by the constant pattern approximation adopted by Cooney and Lighfoot and by the very large integration step of Hsieh et al. Actually those methods require a very small
3.8508 X 1.4263 X 5.2130 X 1.8422 X 5.7288 X 1.2901 X 2.3172 X 3.8385 X 5.7957 x 7.7276 X
7ea
PPI
lo-' 10.' 10.' 10.'
lo-' lo-' 10.'
b
Figure 4. Comparison of the numerical solutions obtained by different methods (data from Table 6 of Hsieh et al., 1977).
Table 11. Parameter Values for Numerical Solutions
component m-xylene
p-xylene
Kj*
hj
Pj*
uj*
CjF
32.97 188.4
7 7
1.004 1.0
439.6 879.2
3.93 3.92
C, = 7.85 M
T', = 1.75 x 10-3
d p / L = 1.4 X E p = 0.2 p = 0.2528 T = 1400 S
E
moug
= 0.42
L = 100 cm PE = 1.7 tiF = 3.0 X
cm/s
amount of computation time and thus must be regarded as extremely efficient, although approximate "short-cut" type methods. This same problem has also been solved by the orthogonal collocation method (Villadsen and Michelsen, 1978). The obtained results have not been reported because they are indistinguishable from our solution. Furthermore, the computation time required by this procedure was comparable to the one required by the technique developed here. Examples of Applications In some industrial applications the adsorbent bed is regenerated by heating or pressure decrease instead of desorption. In such cases, after regeneration, the particle pores are empty and therefore in the subsequent adsorption processes the global mass flow rate varies along the bed. An example of this is the separation process (McKetta and Kobe, 1960) of n-paraffins from mixtures containing isoparaffins and aromatics. The net decrease of the global mass flow rate implies a variation of the space velocity along the adsorber axis. For illustrative purposes the separation process of a m-xylene-p-xylene mixture in which the adsorbent bed is initially assumed empty is examined. All the data used in the calculation have been evaluated elsewhere (Santacesaria et al., 1982) and are
Ind. Eng. Chem. Fundam., Vol. 21, No. 2, 1982
q
L-
5
15-
-
$ 0150
: 0
I
3
10-
e
.-\
;i
' I
I
___ -
conrlanl velocity variable velocily
i"z II
I
as-
0071-
1
! I
0 121
3
129
I
\
\
/
1
I
3
Gtmenrionierr time
Gimensianlerr time 4
4
Figure 5. Effect of variable linear velocity on the breakthrough curves of a binary liquid mixture.
Ql
Figure 7. Effect of axial dispersion on the pulse-response curve for an inert component. I
--- P,-. --P.,!7
IO
I Oimenuonlesr lime 4
0
ia
16
hnenrionlrrr time
20
9
Figure 6. Effect of variable velocity on the exiting flow rate: Q = (+l')vuiable/
(+l)mnst.nt.
summarized in Table 11. Figure 5 shows the step-response of a liquid phase multicomponent adsorption bed assuming constant and variable velocity. In Figure 6 the ratio s2 of the flow rates at the exiting section of component 1(Table 11)evaluated both neglecting and accounting for the space velocity variation along the bed, is shown as a function of time. The value of s2 at each time is then obtained from the ratio of the product $ J ~ V for the two curves reported in Figure 5. It is worthwhile to note that the influence of the constant space velocity approximation on the components flow rate is significant in the time interval in which the separation of the component should be performed. This can be confirmed by comparison of Figures 5 and 6. For gaseous mixtures the variation of space velocity is also due to the pressure drop along the equipment. Although this problem has not been examined here in detail, it can easily be solved by means of the same procedure by simply adding the momentum balance to the system. For dilute mixtures, as reported by Zwiebel(1969), the pressure-velocity relationship can be described by the Ergun equation (Ergun, 1952). Another important problem is the simulation of units with rapidly changing feed composition. An example of this mode of operation are the pulse-response tests frequently used in laboratory runs to identify the physicochemical parameters of the model. Another example is the separation of the xylene isomers (Seko, 1979; Seko et al., 1980; Thornton, 1970) in which, on principle, a sequence of pulses is fed to the adsorber in order to increase its efficiency. The effect of axial dispersion in this case is large.
Figure 8. Effect of axial dispersion on the pulse-response curve for a binary liquid mixture.
Figure 7 reports the well known effect of the Peclet number on the peak shape of an inert component. It may be of some importance to remark that the presence of the internal particle porosity also produces a delay on the peak exit for a nonadsorbable component. In particular from Figure 7 it is observed that the maximum peak occurs at 0 = 1.31, instead of 0 = 1which would be obtained for ep = 0. More important still is the influence of axial dispersion when multicomponent adsorption is performed. Figure 8 shows the influence of the Peclet number on the selectivity of a binary liquid mixture separation (the data are reported in Table 11). The simulation of these features requires a mathematical model which includes both axial dispersion and the ability to take into account pulse types of feed stream. For this reason all the previously presented numerical procedures (Hsieh et al., 1977; Cooney and Lighfoot, 1965) are inadequate. Finally, some numerical results are reported in order to illustrate the behavior of the model for different values of the operating conditions and the physical parameters. The simulation of the separation of a liquid m- (1)-p-xylene (2) mixture is considered. AU the data have been evaluated elsewhere (Santacesaria et al., 1982) and are reported in Table 111. The third component in the table is toluene, which is used, only in the pulse-response runs, as a desorbent. The model examined is the one constituted by eq 21-26. Figure 9 shows a comparison of two pulse-response runs, in which only the space velocity is changed. The velocity affects not only the shape of the breakthrough curves and their position on the time scale but also the efficiency of the whole separation process. As velocity increases, the adsorber capacity increases too, while the
130
Ind. Eng. Chem. Fundam., Vol. 21, No. 2, 1982
Table 111. Parameter Values for Numerical Solutions
Ki*
component m-xylene (1) p-xylene (2) toluene (3)
38.18 218.16 72.72
hi 7.259 7.259 7.259
C, = 9.09 M d p / L= 3.589 X EP = 0.2 p = 0.2183 7 = 2419.5 s
Pi* 0.8669 0.8669 1.o
I', = 1.75 E = 0.42
x
mol/g
X
cm/s
CP,
csta 0.393 0.392 0.0
Ui*
879.73 1759.46 1319.59
3.93 3.92 0.09
L = 39 cm Pe = 1.7 uF = 6.77
Cst = step-feed concentration.
Cp,
= pulse-feed concentration.
Dimensionless time
Figure 9. Effect of the flow rate on the separation of a binary liquid mixture.
9
Figure 11. Sensitivity to X, and K;* of the pulse-response curves.
increased by 20% so that the selectivity of the separation remains unchanged. However, the shape of the pulse and therefore the efficiency of the separation changes. We can then conclude that the degree of separation in multicomponent firred bed adsorption depends not only on the ratio of the equilibrium constants, that is the selectivity, but also on the absolute values of the equilibrium constants themselves.
Acknowledgment We gratefully acknowledge the Italian Consiglio Nazionale delle Ricerche (Progetto Finalizzato Chimica Fine e Secondaria) for financial support. 6
1
0
meniionieii
time
4
Figure 10. Sensitivity to Xi and Ki*of the breakthrough curves.
degree of separation decreases. All this is particularly important in the optimal design and scale-up of this type of equipment. Figures 10 and 11analyze the model sensitivity to the physicochemical parameters for the step and pulse-response respectively. The effects of both the kinetic hi and the equilibrium Ki* parameters have been examined separately. In general it can be observed that the variations in the results obtained are comparable to the assigned incrementa of the parameters, which means that this model does not exhibit parametric sensitivity in the sense introduced by Bilous and Amundson (1956) for chemical reactors. Actually, the pulse-response runs appear to be more sensitive than the stepresponse ones. This fact could make the use of the first type of laboratory experiment commendable to allow a more accurate evaluation of the parameters involved in the model. Another interesting feature of the model is observed in the comparison of curves 1 and 3 shown in Figure 11. The equilibrium constant values of both components have been
Appendix The equations of the model, assuming local equilibrium between the solid and the pore fluid phases and constant space velocity, are as follows.
(-42) Ki* vi Ti =
N
From eq A3 it is readily seen that ayi N ayi avj -=E,-ae Ilavj ae
where ayi/dqj can be evaluated analytically through eq A3. Defining the following matrices P 1 [\ciil;
f
1[vi];
gE
[ril
Ind. Eng. Chern. Fundam., Vol. 21, No. 2, 1982 131
L
= [lij] where lij = 0 for i # j ; Xi for i = j A = [aij]where aij = ayi/dpj D = aA/a6[dij] where dij = duij/a6
I = [ i k j ] where ikj = 0 for k
z j ; 1 for k
= j (A5)
the system of eq Al-A2 and eq A4 can be rewritten a6
‘
.as
Pe\L
Substituting eq A8 in eq A7, and solving with respect to afla6 gives af/ae = [L(P - O][C,I
+ PA]-’
(A9)
which allows the explicit evaluation of all derivatives a(pi/ae. The first derivatives of y i and rli can be readily evaluated by eq A4 and eq A6, respectively. In order to apply Taylor series expansion (eq 28-30), it is now necessary to evaluate the second derivatives of all the involved variables with respect to 6. Differentiating eq AS
In the same way from eq A6, it is possible to evaluate a2p/ae2,which is identical with the nonequilibrium case previously developed. Finally, differentiating eq A7 and substituting eq A10, it follows that
Both first and second derivatives of rli, pi, and yi are then explicitly evaluated, and the numerical procedure can now be performed in the same way as previously described. Nomenclature C = external fluid phase concentration, M C, = reference concentration, M C, = pore fluid phase concentration, M DL = axial effective diffusivity, cm2/s D, = intrapellet effective diffusivity, cm2/s d, = particle diameter, cm k, = kinetic constant, L/mol s K = equilibrium constant, L/mol K* = dimensionless equilibrium constant kl = overall mass transfer coefficient, cm/s kf = film mass transfer coefficient, cm/s k, = internal mass transfer coefficient, cm/s L = bed length, cm N = number of component in the mixture P e = Peclet number for mass transfer; P e = uFL/DLt
R, = particle radius, cm s = dimensionless axial coordinate u = superficial fluid velocity, cm/s V, = particle volume Greek Letters r = solid phase concentration, mol/g r, = solid phase concentration at saturation conditions, mol/g y = solid phase dimensionless concentration t = bed void fraction 6, = particle porosity 6 = dimensionless time X = dimensionless overall mass transfer p = defined by eq 27 v = dimensionless velocity p = molar density, M p* = dimensionless molar density ps = compact solid density, g/L u = dimensionless kinetic constant T = residence time, s cp = pore fluid phase dimensionless concentration = external fluid phase dimensionless concentration Subscripts i = component ith m = mixture Superscripts F = feed 0 = initial conditions s = surface conditions _ -- volume average value Literature Cited
+
Amundson, N. R.; Arls, R. “Mathematical Methods in Chemical Engineering: Fkst-order Partlal Mfferentlal Equations wkh Applications”; vol. 2;Prentlce Hall: EnglewoobCllffs, NJ, 1973. Blious, C. H.; Amundson, N. R. A I C M J. 1958, 2,117. Cooney, D. 0.; Lighfoot, E. N. Ind. Eng. Chem. Fundem. 1965, 4 , 233. Ergun, S. Chem. Eng. Rog. 1952. 48, 89. Frltz, W.; SchlMder, E. U. Chem. Eng. Scl. 1974, 29,1279. Glueckauf, E. Trans. Faraday Soc.1955, 57, 540. GrClbner, 0. A&. Chromatogr. 1988, 6 , 173. Ha!orl, T.; Akizukl, K; Kamogawa, K Murakaml. Y. I n Rws, L. V. C. Proceedlngs of 5th International Conference on Zeolites. in Naples”; Heyden: London, 1980. Hsleh, J. S. C.; Twlan. R. M.; Tien, C. AIChE J. 1977, 23,263. Kucera, E. J. Chromatogr. 1985, 79,237. Lapldus, L.; Amundson, N. R. J. Phys. Chem. 1952, 56, 984. Llapls, A. I.; Rlppln, D. W. T. Chem. Eng. Sci. 1977, 22,619. Llapls, A. I., Rippln, D. W. T. Chem. Eng. Sci. 1978. 33, 593. McKetta, J. J.; Kobe, K. A. “Advances In Petroleum Chemistry and Refining”; WUey-Interscience: New York, 1960. Rosen, J. B. J. Chem. Phys. 1952, 20,387. Santacesarla, E.; Morbldelll, M.; Servida, A.; Storti, G.; Carri, S. Ind. Eng. Chem. Roc. Des. Dev., 1982, In press. Seko, M. 011 Gas J . 1979, 77, 81. Seko, M.; Mlyake, T.; Inada, K. W o c a r b o n Process. 1980, 59(1), 133. Thomas, H. C. J. Am. Chem. Soc. 1944, 66, 1664. Thornton, D. P. &d+ocarbon Process. 1970, 49(1 I), 151. Vank, C. R., Ph.D. Dlssertetlon, Syracuse Unlverstty, Syracuse, NY. 1970. Vermeulen, T.; Klein, 0.: Hlester, N. K. “Adsorption and Ion Exchange” In “Chemical Engineers' Handbook”; Perry, R. H. Chllton, C. H., Ed.; 5th ed; McQraw-Hill: New York, 1973;Chapter 16. Villadsen, J.; Mlchelsen, N. L. “Solution of Differentlal Equatlon Models by Polynomlal Approxlmation”; Prentlce-Hall: Englewood Cliffs, 1978. Von Rosenberg. D. U.; Chambers, R. P.; Swam, 0.A. Ind. Eng. Chem. Fundam. 1977, 76. 154. Zwlebel, I. Ind. Eng. Chem. Fundam. 1989, 8 , 803.
Receiued for review December 30, 1980 Accepted December 14, 1981
