Download Analysis of Mass Transport in Packed Beds: Dispersion Term and Laplace Transform Solutions and more Exams Chemistry in PDF only on Docsity!
MODEL 4
Chromatographic Separations
Mathematical Aspects:
(1) The Laplace transformation. (2) Analysis of transient response by the moment technique.
Chromatography In the separation process known as chromatography, a gas or a liquid containing low concentrations of a mixture of chemical species, passes through a system (e.g. a packed bed) containing a solid phase that interacts with the species in the liquid to provide separation. The mechanism employed is often reversible adsorption on the solid surface. Due to the different affinities of the various species in the mixture to the solid surface (i.e. different quantitative liquid/solid or gas/solid partitioning), some species will be retained longer on the solid phase than others, leading to a later elution. The process is illustrated in Figure 1.
Figure 1. Chromatographic separation: at t=0, a slug containing a mixture A+B is injected into the main flow and the concentrations of A and B are measured in the effluent as a function of time (elution curves). If the two compounds have different affinities for the solid phase (for example, one adsorbs more strongly than the other on the solid surface), elution will occur at different times and the compounds will be separated. In the example shown, B elutes later than A.
The objective of this model is to analyze the evolution of the concentration profiles of a chemical species in a packed bed. Since the chemical species to be separated are usually present in dilute concentrations, transport of each species can be analyzed separately. The analysis will be based on an equation representing the mole balance of a chemical species at a point in the bed. To establish the governing equation, we will use as dependent variable the average species concentration over an averaging volume that is small with respect to the size of the column but
that is sufficiently large to smooth out concentration variations at the pore scale. The averaging volume is illustrated in Figure 2.
Figure 2. Example of an averaging volume in a porous medium.
Let cp be the point (pore-level) concentration of a species in the fluid. The volume-averaged concentration is defined by
V f
p f
c dV V
c (1)
The use of this average concentration allows us to formulate the problem describing the transport of the species without having to calculate the rapidly varying cp. The analysis will be performed under the following assumptions and conditions: (1) Average fluid motion occurs along the main direction of motion (z) only. This means that, even though the velocity field is three dimensional at the pore level, the average velocity vector has only a nonzero z component. The average velocity will be characterized by means of the superficial velocity, v, calculated as
A B
Q
v = (2)
a k (c Kcˆ ) z
c D z
c v t
c 2 v m s
2 − − ∂
ε (8)
This equation has two dependent variables (c and cˆ^ s). To complete the formulation of the
problem, we need a mole balance of species in the solid phase:
k (c Kcˆ ) t
cˆ^ m s
s (^) = − ∂
This mole balance assumes that the only mechanism that affects the adsorbed concentration of species is mass transfer from the fluid. Equations (8) and (9) must be solved subject to appropriate initial and boundary conditions. In a typical chromatographic process, the bed is initially free of solute,
c=0, t=0 (10)
cˆs = 0 , t=0 (11)
and the inlet concentration is a prescribed function of time:
c=ci(t), z=0 (12)
A second boundary condition for c is needed. The second condition is not straight forward sometimes and has been historically a matter of debate in the literature. Under certain conditions, it can be assumed that changes in concentration with position at the exit of the bed are small and state mathematically
z
c
∂
, z=L (13)
where L is the length of packing. Furthermore, if the bed is long enough, condition (13) can be applied letting L→∞, which somewhat simplifies the solution. The PDEs (8) and (9) constitute a system of two linear, coupled, homogeneous differential equations. Note that the only non homogeneous condition is equation (12). To solve this problem, we will apply the Laplace transformation method. In what follows, we will introduce the Laplace transform as a tool for the solution of differential equations.
The Laplace transform The Laplace transform is a particular case of an integral transform, which is an operator that, when applied to a function, generates a different function with a new independent variable. For a function f(t), the general integral transform takes the form
b
a
f (s) K(s,t)f(t) dt (14)
where f (s) is the transform of f(t), a and b are integration limits in the domain in which f is
defined, and K(s,t) is the kernel of the transformation. The functional form of the kernel and the values of a and b define a particular transformation. One of the most useful transforms is the Laplace transform, in which the kernel is an exponential function:
∞ = =^ − 0
f (s) L f(t) e st^ f(t)dt (15)
where the operation L {f(t)} means "Laplace transform of f(t)." The Laplace transform is a useful tool for the solution of linear differential equations. As we will see later in this chapter, the usefulness of the Laplace transform is that it converts ordinary differential equations into algebraic expressions, and partial differential equations into ordinary differential equations. A flow chart of the application of the Laplace transform is shown in Figure 3. The success of the application of the Laplace transform to solve a specific problem is assured only if the three steps in these procedures can be performed analytically. Exception to this are methods that can be used to invert a Laplace domain solution numerically, but these are not used very much in practice, and methods that do not seek to find the final solution but only certain properties of it, such as the method of moments that will be analyzed later in this chapter.
Figure 3. Typical steps followed in the application of the Laplace transformation to the solution of ODEs or PDEs.
The following example gives an introduction to the use of Laplace transform to solve ODEs.
∞ −
∞ −
0
i 1 1 R
st 0
st (^1) (c c) kc dt t
dt e dt
dc e (19)
The left-hand side of this equation can be resolved applying integration by parts, as follows,
∞ − ∞ −
∞ − (^) = + 0
st 10 st 1 0
st (^1) dt e c s e cdt dt
dc e (20)
We restrict the values of the Laplace variable to s>0 and recognize that the last integral in this equation is the transform of c 1 , and use the initial condition (18) to get
0 1 0
st (^1) dt c sc dt
dc
∫e^ =− +
∞ − (^) (21)
Substituting this equation into equation (19) and operating with the right-hand side yields
1 R R
i 1 0 k c t
t
c s c c
Similarly, application of the Laplace transform to equation (17) leads to
2 R R
1 2 0 k c t
t
c s c c
Since ci(t) is known (and, consequently, ci can be determined), equation (22) can be solved to
find c 1 :
+α
s
c t
s
c c i R
0 1 (24)
where
k t
R
α = + (25)
Also, from equation (23) we find
+α
s
c t
s
c c 1 R
0 2 (26)
The next step is to find the inverse transform of the functions in equations (24) and (26). To do this, we need to know the inlet concentration as a function of time, ci(t). Let us take the simple case in which ci is a constant for t≥0. We then have (using equation 15)
s
c e s
c (s) e cdt c i 0 i st 0
i st i =
∞ − − ∞
Equation (24) now can be written as
s(s )
t
c s
c c R
0 i 1 +α
+α
Because of the linear nature of the Laplace operator, and using the notation
L −^1 {f (s)} =f(t) (29)
equation (28) implies (taking the inverse)
s(s )
t
c s
c (t) c^1 R
1 i 1 0 L^ L (30)
The definition of the Laplace transform (equation 15) can be used to show that
{ } +α
−α = s
L e t (31)
To find the second inverse transform in equation (30), we make use of partial fraction decomposition. We seek to find constants A and B such that
B
s
A
s(s )
which implies
1=A(S+α)+Bs (33)
Since this equation must be valid for any s, we can equate coefficients of terms having the same power of s:
s^0 : αA=1 ⇒ A=1/α (34)
s^1 : A+B=0 ⇒ B = − 1 /α (35)
Equation (40) now can be written as follows
α + α
α +α
α
+α
+α
R
i R^2
0 0 2 (s )
(s )
1 s s
t
c (s )
t
c s
c c (44)
Using equations (31) and (41) (or directly equation 15), it is possible to show that
t 2
(^1) ( 1 t)e (s )
− s^ −α = −α
L (45)
Now we can find the inverse transform of equation (44), using equations (31), (38), (41) and (45). The result is, after manipulations,
[ t]
(^22) R
t i R
2 0 1 (^1 t)e t
c e t
t c (t) c 1 −^ α − +α −α α
This example shows how a system of two ODEs can be transformed to a system of two algebraic equations in the Laplace domain, whose solution can be found analytically. A key aspect of the procedure is the knowledge of the inverse Laplace transform of the solution in the Laplace domain. In practice, tables of Laplace transforms available in the literature are used to find the inverse transform. In addition, there are certain properties of the Laplace transformation that aid in the solution of problems such as the one considered in this section. Some of these properties are presented below.
Properties of the Laplace transformation All the properties presented in this section can be proven from the definition of Laplace transform (equation 15). Some of these properties were derived and/or used in the previous section. In what follows, a and b are constants, and n is an integer, n>0.
(1) Linearity:
L {af (t)+bg(t)} =af(s)+bg(s) (47)
(2) Transform of derivatives:
sf(s) f( 0 ) dt
df = −
L (48)
s f(s) [s f( 0 ) s f'( 0 ) f ( 0 )] dt
d f n n 1 n 2 (n 1 ) n
n = − −^ + − + −
L L (49)
(3) Transform of integrals:
ξ ξ
a
0
t
a
f(s) f( ) d s
L f( )d (50)
(4) Shifting properties:
L {e −atf^ (t)} =f(s+a) (51)
L {H (t−a)f(t−a)} =e−asf(s), a>0 (52)
where H is the step function (sometimes called Heaviside function), defined by
1,t a
0 ,t a H( t a) (53)
(5) Other properties:
{ } (^) n
n n n ds
d f L t f(t) =(− 1 ) (54)
{ } f(s/a)
a
L f (at) = , a>0 (55)
(6) Convolution theorem: The definition of the convolution product of two functions, f(t) and g(t) is
t
0
f *g(t) f(t )g() d (56)
The convolution theorem states that
L {f *g(t)} =f(s)g(s) (57)
This theorem is particularly useful to obtain the inverse Laplace transform of the product of two functions in the Laplace domain of two functions whose individual inverse transforms are known. For example, consider the inverse transform of the function
(s )^3
h (s) +α
Here, we can use
vy=vz=0, and that vx=vx(t,y). Since the fluid motion is caused exclusively by the plate motion, there will be no pressure gradients in the x direction, and the x component of the Navier-Stokes equations simplifies to
2
x^2 x y
v t
v ∂
=ν ∂
where ν is the kinematic viscosity (ν=μ/ρ). The initial and boundary conditions are
vx=0, t=0 (65)
vx=v 0 , y=0 (66)
vx=0, y→∞ (67)
For the case of constant v 0 , a simple change of variables will make this problem identical as the leaching from an infinite solid treated at the end Chapter II. The solution was found in that case using similarity transformation. Here, we will find the solution by Laplace transform considering the more general case v 0 =v 0 (t) which, in general, cannot be solved by similarity transformation. We start by taking the Laplace transform of equation (64). For the left-hand side of the equation, use of equations (48) and (65) yields
x x x x (^) sv v ( 0 ) sv t
v = − =
L (68)
For the right-hand side of the equation, we note that the Laplace operator is interchangeable with the derivative with respect to x, since the operator involves integration with respect to time. This means that
2
(^2) x 2
(^2) x dy
d v y
v
L (69)
The transform v (^) x is a function of y and s but, since there are no derivatives with respect to s in
the formulation, we may consider s to be a parameter of the function vx (y). This is the reason for the use of the total derivative in equation (69). Upon substitution of equations (68) and (69), equation (64) becomes
v 0
s dy
d v 2 x
(^2) x = ν
To find the boundary conditions for this ODE, we take the Laplace transform of the boundary conditions (66) and (67) to get
v (^) x = v 0 , y=0 (71)
vx = 0 , y→∞ (72)
The solution of equation (70) is
ν
ν
s Bexp y s vx Aexp y (73)
Application of condition (72) leads to B=0, and condition (71) yields A= v 0 and thus
ν
s v (^) x v 0 exp y (74)
Let us first consider the case of constant v 0. In this case,
s
v v 0 = 0 (75)
and
ν
s exp y s
v v (^) x 0 (76)
To find vx(t,y) we need the inverse transform of this expression. From Appendix D, Rice and Do, transform number 61, we have
s
e 2 t
erfc
−α s
α L (77)
where erfc is the complementary error function, defined by
erfc(x)=1−erf(x) (78)
We see that the inverse transform of equation (76) yields
ν
2 t
y v (^) x (t,y) v 0 erfc (79)
which is the solution for the velocity for the case of constant v 0. Note how it is possible to
express this solution in terms of a similarity variable of the form y / t.
V
n c (^) A 0 = A, t=0 (84)
where V is the volume of reactive mixture on the reactor. From this moment on, no more reactant is added. A mole balance of A in the reactor yields (equation I-9 with cAi=0)
c 0 t
k dt
dc A R
A =
Using the initial condition (84), the solution is
R
A A 0 R t
t c c exp 1 kt (86)
Suppose that the process will be carried out until 90% of the added reactant is either consumed or leaves the reactor; that is, until cA=0.90cA0. This means that the process will be carried out until a time tf such that (equation 86)
R
f R t
t
- 1 exp 1 kt ⇒ R
R f 1 kt
- 3 t t
The loading of the reactor with reactant at the beginning of the process is not instantaneous in practice. Let t 0 be the time necessary to loan the nA moles of reactant at the beginning of the process. We can see that the loading process can be considered instantaneous for all practical purposes if t 0 «tf. That is, equation (86) will be valid only if this constraint is satisfied. Consider now the case in which the process is cyclic: at time tf, a load of 0.9nA moles of reactant is added to the reactor (over a period of time t 0 «tf) to compensate for the decay in reactant concentration. This means that at time tf the concentration of A in the reactor will be again cA0 and the process will start over. The analysis of this case is simple: we can "shift" the time scale so that t=0 coincides with the beginning of each cycle, and equation (86) then will represent the concentration evolution within each cycle. Alternatively, we can formulate the problem so that we incorporate each load mathematically into the mole balance. This would be more convenient to treat problems with a more complex loading pattern. One way to do this is to add a term in the mole balance (85) to incorporate addition of A:
c L(t ) t
k dt
dc A R
A =
Each term in this balance represents moles of A per unit time and unit reactor volume so that L(t) must represents the moles of A added at time t per unit time and reactor volume. Since the addition is instantaneous, L(t) must be infinite at the loading time t=tf (and all subsequent cyclic loading times) since a finite number of moles nA are added in a very short time interval (∆t→0). That is
f
f 0 ,t t
,t t L( t) (89)
But L(t) must reflect the fact that nA moles are added so that, over the loading interval ∆t→0, we must have that
V
n L( t)dt A
tf t/ 2
tf t/ 2
∫^ =
− ∆
A mathematical representation for L(t) is accomplished by means of the Dirac delta function, δ(t), which has the properties (see more formal definition below)
δ = 0 ,t 0
,t 0 (t ) (91)
and
(t)dt (t)dt 1
t/ 2
t/ 2
∫ δ^ = ∫δ^ =
∆
−∆
∞
−∞
Using this definition, we can incorporate instantaneous loads into the mole balance as follows. Let nAn be the number of moles of A loaded into the reactor at time tn, with n=1,2,… If these loads are instantaneous then we can say that the nth value of L(t) is
(t t ) V
n L (^) n (t)= Anδ − n (93)
Note that this satisfies equation (90) due to property (92). The sum of all loads will give the total addition of A, so that the mole balance, equation (88), can be expressed as
(t t ) V
n c t
k dt
dc n alln
An A R
A (^) = δ −
In summary, multiple additions of reactant can be treated two ways: (1) use the solution (96) to represent concentration evolution between loads, and then shift the time axes and apply the appropriate initial condition to solve for each period between loads, or (2) solve directly equation (94). The latter treatment is more convenient when using Laplace transform to find a solution, as will be illustrated in the next section. Equations (91) and (92) represent properties of the Dirac delta function, but they are far from being a mathematically rigorous definition. A formal definition follows from the use of sequences of functions that are continuous and differentiable at each point. Consider the following set of functions:
where a−∆t/2≤ξ≤a+∆t/2. The integral on the right-hand side of this equation is unity. Letting ∆t→0 yields the following property for the delta function:
∫ (t a)f(t)dt f(a)
∞
−∞
δ − = (102)
valid for any function f(t).
t
ψn
n=
Figure 6. The first three functions of the sequence given by equation (95).
Next, notice that
−∞ 1 ,t a
0 ,t a ( a)d
t (103)
which is the definition of the step function (equation 53). We conclude that
( a)d H(t a)
t
∫δτ− τ=^ −
−∞
or
(t a ) dt
dH (t a ) =δ −
Finally, the Laplace transform of the delta function is
∞ δ − = −^ δ − =^ − 0
L (t a) e st^ (t a)dt e as, a≥ 0 (106)
Note that, for a=0,
L {δ ( t)} = 1 (107)
Now we can solve equation (94). Taking the Laplace transform of equation (94) and using equation (106) leads to
exp( st ) V
n c t
s c k n alln
An A R
A = −
which leads to
R
n alln
An A s k 1 /t
exp( st ) V
n c
Taking the inverse transform of this equation and taking advantage of the linearity of the operation yields
R
n all n
A An^1 s k 1 /t
exp( st ) V
n c L (110)
Consider the property (52). We can see that
L −^1 {exp( −stn)f(s)} =f(t−tn)H(t−tn) (111)
In this case (equation 110), we have
s k 1 /t R
f (s)
which means that (see equation 31)
f (t)= exp[− ( k+ 1 /tR)t] (113)
