Effect of radial velocity profiles on axial dispersion in packed beds: asymptotic behaviour

Abstract

Axial dispersion of a solute in a flow through packed beds arises from molecular diffusion, velocity variations at the interstitial-pore scale and systematic radial velocity profiles. Currently, there is not a widely accepted procedure for estimating the last contribution. One possible reason is the fact that radial velocity profiles in packed beds are intimately related on the radial porosity profile, which is not a deterministic property, and velocity variations are relatively mild and confined to a region close to the vessel walls. A widespread notion is that only in the zone about a particle radius from the wall is the fluid velocity significantly larger (wall channelling) than in the remaining of the bed core. Based on literature results and effective-model evaluations, the impacts of porosity profiles and related velocity profiles are explored in this manuscript. It is found that the dispersion of a solute can be increased several times in magnitude solely due to the wall channelling when the detailed porosity profile is considered. This can be explained in terms of small velocity differences between a second region (up to around five particle diameters) and the innermost bed region. Guidelines are also discussed for predicting likely levels of dispersion. A further aspect concerns the simplification of the treatment in the wall-zone.

Appendices

Appendix 1: Taylor-Aris formulation

The mass conservation Eq. (1) for a species in a packed bed is re-written here:

$$\psi (\rho )\frac{\partial C}{{\partial t}} = \frac{1}{{r_{t}^{2} }}\frac{\partial }{\rho \,\partial \rho }\left( {\rho D_{rad} (\rho )\,\frac{\partial C}{{\partial \rho }}} \right) + D_{ch} (\rho )\frac{{\partial^{2} C}}{{\partial \,z^{2} }} - \psi (\rho )\,v(\rho )\frac{\partial C}{{\partial z}};\quad \frac{\partial C}{{\partial \rho }} = 0\quad \rho \, = \,0,{1}$$

(40)

Equation (40) will be used for a bed infinitely extended in the axial direction with a tracer initially present (t = 0) according to an axi-symmetric distribution C_I (ρ, z), and without being fed at t > 0. For convenience, a change of variables from z to ζ it is made, with ζ = z-$\overline{v}$ t, where $\overline{v}$ = $(2/\varepsilon )\int_{0}^{1} {\psi (\rho )\,v(\rho )\,\rho \,d\rho }$ is the average interstitial velocity and ζ is the mobile axial coordinate. Equation (40) then becomes:

$$\psi (\rho )\frac{\partial C}{{\partial t}} = \frac{1}{{r_{t}^{2} }}\frac{\partial }{\rho \,\partial \rho }\left( {\rho D_{rad} (\rho )\,\frac{\partial C}{{\partial \rho }}} \right) + D_{ch} (\rho )\frac{{\partial^{2} C}}{{\partial \,\zeta^{2} }} - \psi (\rho )\,\overline{v}\,\chi (\rho )\,\frac{\partial C}{{\partial \zeta }};\quad \frac{\partial C}{{\partial \rho }}\, = \,0,\rho \, = \,0,{1}$$

(41)

where the relative excess velocity $\chi = (v - \overline{v})/\overline{v}$ verifies:

$$\int_{0}^{1} {\psi (\rho )\chi (\rho )\,\rho \,d\rho } = 0$$

(42)

Note that in Eq. (41) ($\partial C/\partial t$) is the partial derivate of C at constant ζ.

The normalized local spatial moment of order n is defined as:

$$c_{n} (t,\rho )\, = \,\tfrac{1}{{M_{I} }}\int_{ - \infty }^{\infty } {\zeta^{n} \,C(t,\zeta ,\rho )\,d\zeta } ;$$

(43)

where M_I = N/(ε π $r_{t}^{2}$) and N is the initial moles of tracer, which will be kept inside the bed at all times. The cross-section averaged spatial moment of order n is defined as:

$$m_{n} (t)\, = \,\left( {2/\varepsilon } \right)\int_{0}^{1} {\psi (\rho )\;c_{n} (t,\rho )\,\rho \,d\rho }$$

(44)

It is assumed that C_I (z,ρ) is such that c_n(0,ρ) is finite and that the origin z = 0 is chosen in such a way that m₁(0) = 0. Multiplying Eq. (41) by ζⁿ and integrating on ζ =  ± ∞ gives:

$$\psi (\rho )\frac{{\partial c_{n} }}{\partial t} = \frac{1}{{r_{t}^{2} }}\frac{\partial }{\rho \,\partial \rho }\left( {\rho \,D_{rad} (\rho )\,\frac{{\partial c_{n} }}{\partial \rho }} \right) + n(n - 1)D_{ch} (\rho )\,c_{n - 2} + n\,\psi (\rho )\,\overline{v}\,\chi (\rho )\,c_{n - 1}$$

(45)

where $\partial c_{n} /\partial \rho = 0$, ρ = 0,1. Differentiating Eq. (44) with respect to t and using Eq. (45) for $\partial c_{n} /\partial t$:

$$\frac{{dm_{n} }}{dt} = n(n - 1)\frac{2}{\varepsilon }\int_{0}^{1} {D_{{^{ch} }} (\rho )\;c_{n - 2} (t,\rho \,)\rho \,d\rho } + n\;\overline{v}\frac{2}{\varepsilon }\int_{0}^{1} {\psi (\rho )\chi (\rho )\,c_{n - 1} (t,\rho \,)\;\rho \,d\rho }$$

(46)

In particular for n = 2:

$$\tfrac{1}{2}\varepsilon \frac{{dm_{2} }}{dt} = 2\int_{0}^{1} {D_{{^{ch} }} (\rho )\;c_{0} (t,\rho \,)\rho \,d\rho } + \;2\,\overline{v}\,\int_{0}^{1} {\psi (\rho )\chi (\rho )\,c_{1} (t,\rho \,)\;\rho \,d\rho }$$

(47)

The previous formulation is now applied for the so-called axial dispersion plug flow model, which assumes transversally uniform properties, porosity ε, velocity $\overline{v}$, axial dispersion coefficient D_ch = D_ax, along with an initial distribution C_I(z) also uniform on ρ. The terms in Eqs. (40, 41, 45) involving radial derivatives vanish. Under these conditions, it becomes, in particular, c_n = m_n, $\chi (\rho )$ = 0 and c₀ = 1 (at any t). Equation (47) is then reduced to:

$$D_{ax} \, = \,\tfrac{1}{2}\varepsilon \left( {dm_{2} /dt} \right)$$

(48)

Based on this result, for the two-dimensional model Eq. (40), an effective axial dispersion coefficient is analogously defined at a given t, with dm₂/dt evaluated from Eq. (47). The moments c₀ and c₁ needed for this end are to be evaluated from Eq. (45):

$$\psi (\rho )\frac{{\partial c_{0} }}{\partial t} = \frac{1}{{r_{t}^{2} }}\frac{\partial }{\rho \,\partial \rho }\left( {\rho D_{rad} (\rho )\,\frac{{\partial c_{0} }}{\partial \rho }} \right);\quad \partial c_{0 } /\partial \rho \, = \,0\,{\text{ en }}\,\rho \, = \,0,{1}$$

(49)

$$\psi (\rho )\frac{{\partial c_{1} }}{\partial t} = \frac{1}{{r_{t}^{2} }}\frac{\partial }{\rho \,\partial \rho }\left( {\rho D_{rad} (\rho )\,\frac{{\partial c_{1} }}{\partial \rho }} \right) + \psi (\rho )\,\overline{v}\chi (\rho )\;c_{0};\quad \partial c_{{1}} /\partial \rho \, = \,0\,{\text{ en }}\,\rho \, = \,0,{1}$$

(50)

The two integrals in Eq. (47) are defined as:

$$\tfrac{1}{2}\varepsilon \,\left( {dm_{2,ch} /dt} \right) = \overline{D}_{ch,I} (t) = 2\int_{0}^{1} {D_{{^{ch} }} (\rho )\,c_{0} (t,\rho \,)\rho \,d\rho }$$

(51)

$$\tfrac{1}{2}\varepsilon \,\left( {dm_{2,v} /dt} \right)\, = \,D_{v} \left( t \right)\, = \,u\frac{2}{\varepsilon }\int_{0}^{1} {\psi (\rho )\chi (\rho )\,c_{1} (t,\rho \,)\;\rho \,d\rho }$$

(52)

where $\overline{v}$ = $u/\varepsilon$ was used. The second moment m₂ has been split as m₂ = m_2,ch + m_2,v and the overall dispersion coefficient D_ax (Eq. 48) as:

$$D_{ax} \, = \,\overline{D}_{ch,I} \, + \,D_{v}$$

(53)

$\overline{D}_{ch,I}$ is the average superficial coefficient of channel dispersion and D_v is the superficial trans-column dispersion coefficient caused by the non-uniform profile v(ρ). Both $\overline{D}_{ch,I}$ and D_v depend on the initial distribution C_I(ρ, z) and time t. However, for any C_I(ρ, z), Eq. (49) indicates that, when $t \to \infty$, c₀ becomes uniform on ρ, and from the definition in Eq. (43) it follows that the asymptotic value is c₀ = 1. Using this result in Eq. (50), it can be concluded that c₁ reaches a stationary profile (i.e. $\partial c_{1} /\partial t = 0$) as t → ∞. Denoting as $c_{1}^{\infty } (\rho )$ the stationary profile of c₁, it is obtained from Eq. (50):

$$c_{1}^{\infty } (\rho ) = c_{1}^{\infty } (0) - \overline{v}\,r_{t}^{2} \int_{0}^{\rho } {V(\rho^{\prime})} \frac{{\quad d\rho^{\prime}}}{{D_{rad} (\rho^{\prime})\rho^{\prime}}}$$

(54)

where $c_{1}^{\infty } (0)$ is the value at ρ = 0, and

$$V(\rho ) = \int_{0}^{\rho } {\psi (\rho ^{\prime})\,\chi (\rho ^{\prime})\,\rho ^{\prime}\,d\rho ^{\prime}}$$

(55)

It is noted that V(1) = 0, by virtue of Eq. (42). Equation (53) when t → ∞ is written as:

$$D_{ax,\infty } \, = \,\overline{D}_{ch} \, + \,D_{v,\infty }$$

(56)

With c₀ = 1 in Eq. (51):

$$\overline{D}_{ch} = 2\int_{0}^{1} {D_{{^{ch} }} (\rho )\rho \,d\rho }$$

(57)

Substituting $c_{1}^{\infty }$ ≡ c₁ from Eq. (54) in Eq. (52), integrating by parts, having into account definition (55) and replacing $\overline{v}$ = $u/\varepsilon$, it is obtained for D_v,∞:

$$D_{v,\infty } \, = \,\left( {\frac{{u\,r_{t} }}{\varepsilon }} \right)^{2} 2\int_{0}^{1} {\;\frac{{V^{2} (\rho )\quad d\rho }}{{D_{rad} (\rho )\rho }}}$$

(58)

where the value $c_{1}^{\infty } (0)$ (Eq. 54) has no effect on account of Eq. (42).

It is concluded that $\overline{D}_{ch}$ (Eq. 57) and D_v,∞ (Eq. 58) become independent of C_I(ρ, z).

Appendix 2: Dispersion in a bed with uniform oscillatory transversal porosity

Consider a packed bed between parallel plates with transversal porosity profile, far from the plates, given by ψ(ρ) = ε [1 + 0.5 cos(π ρ)], where ρ is the transverse coordinate in units of particle radius, d_p/2. The amplitude of the waves is similar to the maximum of the waves in Fig. 2a and its wavelength is d_p, also close to those of the actual profiles. The periodic velocity profile v(ρ) induced by ψ(ρ) will generate an asymptotic axial dispersion coefficient, D_wave, that can be evaluated from the TAMM. To this end, it is only necessary to consider an extension of half a wavelength, i.e., the range 0 < ρ < 1. Assuming the relation between v(ρ) and ψ(ρ) given by Eqs. (20–23) in the main text, applied to the Cartesian coordinate rather than the cylindrical radius, it is obtained $\psi_{c}^{2} = \int_{0}^{1} {\,\psi^{2} (\rho )\;\,d\rho }$ = $\tfrac{9}{8}$ ε², and v(ρ)/$\overline{v}$ = χ(ρ) + 1 = $\varepsilon \,\psi (\rho )/\psi_{c}^{2}$ = $\tfrac{8}{9}$[1 + 0.5 cos(π ρ)].

The equivalent expressions for V(ρ) and D_v,∞ = D_wave in the slab geometry (cfr. Eqs 8 and 4) with D_rad (ρ) = D_R (uniform) are:

$V(\rho )$ = $\int_{0}^{\rho } {\,\chi (\rho ^{\prime})\,\psi (\rho ^{\prime})\,d\rho } ^{\prime}\,$;

D_wave = $\left( {\tfrac{1}{2}u\,d_{p} /\varepsilon } \right)^{2} (1/D_{R} )\int_{0}^{1} {\;V^{2} (\rho )d\rho }$.

It is obtained D_wave $\cong \tfrac{1}{512}\,\,(u\,d_{p} )^{2} /D_{R}$, or 1/Pe_wave $\cong$$\tfrac{1}{512}$ Pe_R. Considering very large values of Pe_m, Pe_R = 9 and 1/Pe_ch = $\overline{D}_{ch}$/(u d_p) = 0.5 can be taken as a reference. The value 1/Pe_wave $\cong$ 0.0176 thus obtained is significantly lower than 1/Pe_ch = 0.5, and it can be concluded that the waves of a stationary profile v(ρ) do not contribute significantly.