6200
Ind. Eng. Chem. Res. 2005, 44, 6200-6211
Design Criteria for Soil Cleaning Operations in Electrokinetic Remediation: Hydrodynamic Aspects in an Annular Geometry Mario A. Oyanader,† Pedro E. Arce,*,‡ and Andrew Dzurik§ Department of Chemical Engineering, Universidad Cato´ lica del Norte, Avda. Angamos 0610, Antofagasta, Chile, Department of Chemical Engineering, Prescott Hall 214, Tennessee Technological University, Cookeville, Tennessee 38505, and Civil and Environmental Engineering Department, Florida State University, 2525 Pottsdamer Street, Tallahassee, Florida 32310-6046
The electrokinetic process of soil remediation is affected by different transport driving forces that are responsible for the motion of the bulk fluid and solute species. In particular, the electromechanisms, e.g., electroosmosis, electrophoresis, and electromigration, may compete with buoyancy and advection, promoting distinct flow regimes. The earlier applications of electrokinetic phenomenase.g., electroassisted drug delivery, electrophoretic separations, and material processing, just to name a few mainly in the area of electrophoresissneglected this competition, and therefore, the hydrodynamics of electrokinetic systems was considered simpler. Field test results demonstrate that this is not the case with soil cleanup processes. The unique characteristics of soil porous media call for a different approach and are in need of further analysis. In this contribution, the basic aspects of the behavior of such a system are captured by using an annular capillary model. Under the proposed geometry, a differential model is formulated using simplifying assumptions to maintain the mathematical aspects at a minimum level and a solution is presented for the different fields, e.g., the temperature and the velocity. Several numerical examples are shown to portray the flow situations found in the system for a selection of values of the parameter space. From the analysis of these graphical representations, a qualitative and semiquantitative description of the different flow regimes inside the annular channel is obtained. Particularly interesting in this study is the inclusion of a resistive heating effect in the core of the annular capillary channel as a force term. Temperature developments are explained and analyzed under different scenarios. This information is useful to identify further aspects for the investigation and to delineate a systematic approach for a more rigorous description. Implications for the design of devices and cleaning strategies are also included. The results obtained in this study are useful to promote a deeper understanding of the behavior of the system, to have a better idea about the experimental effort needed for validation of the different trends, and to lead to important guidelines for improving the separation or cleaning efficiency in a given application. 1. Introduction and Motivation The use of electrical fields in soil decontamination has generated different versions of the intended application, e.g., cation selective membrane, ceramic casting, Lasagna, electrochemical ion exchange, electrokinetic bioremediation, electrochemical geooxidation, and electrosorb, to name some.1,2 All these versions of technology constitute the so-called “electrokinetic remediation” in which electrical fields are applied to the impacted area to induce controlled contaminant movement. The inability of most technologies to succeed in heavy metals soil remediation has invoked other approaches for soil cleanup involving electrical fields. The test results of such approaches indicate that this promising technique is not only effective with heavy metals but with radioactive species, DNAPLs, petroleum hydrocarbons, and several organic compounds as well2-4 and that it can be applied in unsaturated and saturated contaminated soils.2 Nevertheless, and despite these promising re* To whom correspondence should be addressed. Tel.: (931) 372-3189. Fax: (931) 372-6352. E-mail:
[email protected]. † Universidad Cato ´ lica del Norte. ‡ Tennessee Technological University. § Florida State University.
sults, uncertainty, lack of protocols, inability to extrapolate results, and scaling problems are some of the technology’s drawbacks in in situ implementations.2 From a practical point of view, it seams that the missing link to success is the understanding of the fundamentals applied to electrokinetic processes in the soil matrix. Understanding the principles of electrokinetics is useful for design and operation purposes in the processes involved. In the area of electrophoresis, these principles are fairly well-understood for most of the earlier applications of electrokinetic phenomena. For instance, it is very well-known that a combination of different transport driving forces are responsible for the motion of solute through a capillary channel. Although the movement of the solute species is the most important result, it can be caused by different mechanisms and not necessarily by electrophoresis only, which is characteristic of separation processes. On the contrary, environmental applications in porous media are mostly characterized by electroosmosis and electromigration as the driving mechanisms. For simplicity, the way in which these mechanisms collaborate with other driving forces, e.g., buoyancy and hydrodynamics, has been either ignored or neglected, as is the usual rule of thumb in chemical processes. The characteristics of soil, e.g.,
10.1021/ie049045p CCC: $30.25 © 2005 American Chemical Society Published on Web 07/08/2005
Ind. Eng. Chem. Res., Vol. 44, No. 16, 2005 6201
heterogeneous, anisotropic, low-permeability, etc., call for a different approach, where the competition among buoyancy, hydrodynamics, and electroosmosis/electromigration cannot be over looked because it may promote distinct flow regimes and, therefore, affect solute transport. The authors in previous works have demonstrated and reported these hydrodynamics aspects using rectangular and cylindrical geometries.5-7 Although the basic equations for the transport of solutes, under an applied electrical field, are wellknown,3,8,9 the systematic analysis of the different parameter effects and value ranges has not been completed so as to know how to avoid the pitfalls of electrokinetic field testing. To the best of our knowledge, the two contributions made by Oyanader et al.5,6 are the first and latest efforts to provide design criteria for the different electrokinetic applications involving the three main driving forces, i.e., pressure, buoyancy, and electrical. Initially, these authors conducted their systematic analysis using a rectangular geometry, and recently, they concentrated their study in a cylindrical capillary system. Although the objectives of these works have been rather modest, considering the limitation of rectangular and cylindrical systems, the use of both geometries has revealed interesting hydrodynamics aspects to be considered in the design and operation of electrokinetic processes. These findings justify the continuation of a similar study in an annular representation of a porous media. The authors are aware that no single geometry can completely described field conditions; however, there are further aspects not yet considered that can be studied by adopting this proposed annular geometry. In fact, this approximation is more realistic for a porous media than a rectangular or cylindrical geometry in that it provides additional space for the analysis of the effects of porosity and resistive heating of the core on the hydrodynamics of the system. Consequently, this work concentrates on developing the respective annular model and seeks the understanding of the different flow regimes that may be possible in an annular channel. The work presented in this contribution follows the spirit of Professor Dudukovic’s style in that important technological applications can be further improved by the use of fundamental techniques and approaches; we are grateful for the opportunity to contribute to this issue honoring his many contributions, as we have benefited tremendously from his work during the years. 2. Model Formulation Consider the system being analyzed as an annular channel of length L, inner radius R1, outer radius R2, and an inclination of an angle R with respect to the horizontal line (Figure 1). The annular channel is exposed to a constant electrical field E. The inner and outer walls of the channel have net but constant and uniform charges, Ψ*1 and Ψ*2, respectively. The notation ψ(ξ), in Figure 1, is used to indicate the electrostatic potential across the radial domain of the channel. The outer wall surface presents a uniform temperature as defined by its interaction with the temperature of the surroundings T∞. The axes (r and x) have been placed coincidently with the lower end of the capillary channel, x, and the origin of the r-axis at the center of channel. This choice of the coordinate axis is the normal convention in cylindrical coordinates. The system being modeled requires one to consider at least three main aspects of transport phenomena for
Figure 1. Geometrical sketch of the annular capillary channel and coordinate system used in the analysis.
the proposed analysis: heat transfer, electrostatics, and hydrodynamics. From the respective transport equations, profiles of temperature, electrostatic potential, and axial velocity will be obtained for the study of the system behavior. In the sections below, the description of these different aspects is included. 2.1. The Heat Transfer Model. The annular channel (see Figure 1) is assumed to interchange heat with its surroundings within a range of efficiencies, i.e., the Nusselt number may take values within a considerable range of values. In this study, the most illustrative and representative values have been considered (Nu ) 5 and Nu ) 10). The situation being described may also lead to the scenario in which the temperature of the outer wall surface of the annular channel at the position r ) R2 (ξ ) 1) reaches the value of the temperature of its environment T∞; this value is considered to be constant. In addition, the ratio R2/L is assumed to be small enough to neglect any end effects on the temperature profile inside the channel. Furthermore, the conductiondominated regime 10 is assumed to be valid for the present analysis. An applied electrical field E in the axial direction (e.g., the x-axis) of the channel is present, and because of the fact that the fluid inside the annular channel and the core show a nonzero resistance to the electrical current, heat generation takes place. For this study, the generated heat is assumed to be constant with time and uniform across the annular channel. This situation yields two cases of heat transfer that need to be analyzed separately; these are the heat generation in the core and in the fluid. The analysis of the core presents similarities with the electric heating of a wire for which Bird et al.11 conducted a preliminary modeling. In this case, the wire is the core of the annular channel. By applying Fourier’s law, the following equation is obtained:
-
∂θ ) φe2ξ ∂ξ
(1)
The following definitions have been used for the nondimensional expressions:
(T - T∞) SeR22 r 2 ; φe ≡ ξ≡ ; θ≡ R2 T∞ 2KCT∞
(2)
where φe2 is the resistive heating number, Se is the heat generation rate per volume unit, KC is the thermal conductivity of the core of the annular channel, and T∞ is the ambient temperature, beyond r ) R2 (ξ ) 1).
6202
Ind. Eng. Chem. Res., Vol. 44, No. 16, 2005
Equation 1 can be easily solved by assuming that KC is constant and using the boundary condition
θ ) θ(b) @ ξ ) b
(3)
Under the assumptions described above, the energy equation11 reduces to
-
1 ∂ ∂θ ξ ≡ φ2 ξ ∂ξ ∂ξ
( )
(4)
In the right-hand side of eq 4, the Joule heating number has been identified as
φ2 ≡
QR22 KFT∞
(5)
where Q is the Joule heating generation and KF is the thermal conductivity of the fluid of the annular channel. The conservation of energy equation (eq 4) needs boundary conditions at both the inner wall and the outer wall surface of the annular channel; although different types of boundary conditions are possible12,13 in this analysis, we are interested in a simple case. Therefore, the flux or Robin boundary conditions12 are selected and they are given by the following equations:
-
∂θ ) φe2ξ @ ξ ) b ∂ξ
(6a)
-
∂θ ) Nuθ @ ξ ) 1 ∂ξ
(6b)
In the previous expression, eq 6b, the Nusselt number, Nu, is defined as
Nu ≡
hR2 KF
(7)
where h is the convective coefficient and all the other parameters have been defined above. 2.2. Electrostatic Model. As was described above, the annular channel is assumed to have charged walls and the electrostatic potential across them is described by the Poisson-Boltzmann expression.3 This equation assumes that the applied electrical field will produce an electrical work that is small compared to the thermal energy. A common practice is solving the PoissonBoltzmann equation using the Debye-Hu¨ckel approximation,3 which leads to the following form of the simplified model:
1 ∂ ∂ψ ξ - λ2ψ ) 0 ξ ∂ξ ∂ξ
( )
(8)
The notation λ corresponds to the dimensionless inverse of the Debye length, k,3 and is computed by multiplying k by the channel outer radius, R2, e.g., λ ) kR2. The inverse of the Debye length can be viewed as a measure of the power of the action of the electrostatic forces. Following the description on Figure 1, the boundary conditions for eq 8 are
ψ ) ψ/1 @ ξ ) b
(9a)
ψ ) ψ/2 @ ξ ) 1
(9b)
Now, the solution to the Debye-Hu¨ckel approximation, eq 8, is useful to obtain the electroosmotic force term in the Navier-Stokes equation of motion. This will be the subject matter of the next section. 2.3. Hydrodynamic Model. The fluid in the annular channel, as described above, is assumed to be Newtonian, incompressible for the mass conservation aspects, and under steady-state conditions. This fluid is also assumed to have constant properties everywhere except for the density in the buoyancy force term. This is, in fact, the assumption suggested by Boussinesq.14 All the assumptions described in Section 2 are assumed to be valid for the hydrodynamic flow problem as well. In particular, the “no end effects” and the conductiondominated regime (e.g., small magnitude of velocity field) are invoked here. Moreover, a pressure gradient is assumed to be present, but its magnitude must be relatively small to comply with the assumption of a small velocity field. Under these assumptions, the axial or x-component of the Navier-Stokes equation11 is given by
(
)
∂Vx 1 ∂ ∂p - F(T)gx - Fe(ψ)Ex rµ ) r ∂r ∂r ∂x
(10)
where the applied electrical field in the axial direction, Ex, is assumed to be constant at a porosity of 100% (Ex0 @ b ) 0) but related to the area of the porous, annular area of the channel, according to15
Ex )
Ex0 1 - b2
(11a)
and the function F(T) is computed by a first-order Taylor approximation around a mean temperature Tm of the system11
F(T) ) F(Tm) - βmFm(T - Tm)
(11b)
and where βm is the volumetric compressibility of the fluid at a mean temperature Tm. The parameter Tm is determined by the total mass conservation condition that may be stated as
2π
∫RR
2
1
F(Tm)Vx(y,Tm)r dr ) 0
(12a)
or, as a dimensionless equation,
2π
∫b1 FmVx+(ξ,θm)ξ dξ ) 0
(12b)
The condition given by eq 12 requires the computation of the hydrodynamic velocity profile prior to its solution. Finally, the function Fe, the electrostatic density, is defined by the following expression
Fe ) -
k2 ψ(ξ) 4π
(13)
where is the media permittivity and k is the inverse of the Debye length. Equation 13 features the electrostatic potential, ψ(ξ), obtained by solving the complete Poisson-Boltzmann equation using the method proposed by Oyanader and Arce,16 the ×c4AO correction function. To have a convenient way of analyzing the different aspects related to the velocity profile Vx(r), the following dimensionless variables and numbers are proposed.
Ind. Eng. Chem. Res., Vol. 44, No. 16, 2005 6203
VE )
Ex0ψ2 4πµ
(14a)
Vx VE
(14b)
Vx+ )
Gr )
βmFm2R23T∞g µ2 Gr βmT∞
(14d)
VER2Fm µ
(14e)
∂hP Gr* Gr* + sin(R) ∂x Re Re
(14f)
Gr* ) Re ) Pm )
(14c)
The Grashof numbers, Gr and Gr*, represent the buoyancy to viscous forces due to changes in temperature and density, respectively, while in eq 14e, the Reynolds number, Re, represents the inertia to viscous forces. A convenient combination of the Grashof and Reynolds numbers has been mathematically applied to dimensionally reduce the total hydraulic head gradient, which yields the dimensionless number Pm, eq 14f. By using these dimensionless numbers and variables in the Navier-Stokes component, eq 10, the following dimensionless differential equation is obtained
( )
+ Gr λ2 ψ(ξ) 1 ∂ ∂Vx ξ ) Pm sin(R) (θ - θm) + ξ ∂ξ ∂ξ Re 1 - b2 ψ2 (15)
This equation shows several terms that account for the different forces present in the system. The left side of the equation is the viscous term; on the right side, the first term represents the pressure-driven force, the second term is related to the buoyancy effects, and the last one is associated with the electroosmosis term. To complete the problem for the electro-hydrodynamic velocity profile, the nonslip-boundary conditions, at the inner and outer wall surfaces, will be assumed for the annular channel
Vx+ ) 0 @ ξ ) b
(16a)
Vx+ ) 0 @ ξ ) 1
(16b)
Equation 15 can be trivially solved by simple substitution of the temperature and electrostatic potential profile equations followed by mathematical integration and constants’ evaluation using the boundary conditions of eqs 16a and 16b. The mathematical result is described in the next section. 3. Solution To Model Equations This section focuses on the development of the solution to the differential models for the heat transfer, the electrostatic field, and the hydrodynamics velocity profile. These models, presented in Section 2, require a simple but specific procedure to be solved since they are sequentially coupled. First of all, the heat transfer and electrostatic models must be solved to obtain the tem-
perature and electrostatic potential profile equations. Second, the velocity profile must be computed using the hydrodynamic model. After an expression for the velocity profile has been determined, the result is introduced in the mass conservation equation in order to isolate the characteristic mean temperature of the system, Tm. This last step constitutes the closing loop of the solution procedure for the entire system. The detailed results of the described strategy are discussed in the next subsections. 3.1. Heat Transfer Model Solution. As described in Section 2.1, in terms of heat transfer, there are two differential models to solve. The first one is related to the core of the annular channel represented by eq 1 with boundary condition eq 3. The integrated solution is computed as
θ(ξ) ) θ(b) +
{ ( )}
φe2b2 ξ2 1; 0eξeb 2 b
(17)
The second solution, related to the fluid and described by eq 4 with boundary conditions eqs 6a and 6b, is readily computed as
θ(ξ) )
{
}
φ2 2 φ2 φ2 2 - ξ + b - φe2b2 ln|ξ| + 4 4 2 1 φ2 φ2 2 - b + φe2b2 ; b e ξ e 1 (18) Nu 2 2
{
}
Equation 18 is an analytical function of the position of the cylindrical channel, across the radial direction, and it is very useful in the computation of the hydrodynamic velocity profile to be described in a later section. However, some interesting limiting cases can be derived from eq 15. For example, the situation of a high convective cooling system with high Nusselt number values leads to
θ(ξ) )
{
}
φ2 2 φ2 φ2 2 - ξ + b - φe2b2 ln|ξ| (19a) 4 4 2
and for any Nusselt number value, the minimum temperature is located at the wall surface, location ξ ) 1, which is exposed to ambient temperature T∞. This is
θmin )
{
1 φ2 φ2 2 - b + φe2b2 Nu 2 2
}
(19b)
This situation produces the lowest temperature in the fluid system for any values of the Joule heating parameter, φ2, and resistive heating number, φe2. Also, the temperature difference between any value and the lowest value is readily given by
∆θ(ξ) )
{
}
φ2 2 φ2 φ2 2 - ξ + b - φe2b2 ln|ξ| (19c) 4 4 2
This equation becomes useful to predict the temperature difference at any location of the annular channel and the surface of any such domain. It is particularly interesting that eq 19c is equivalent to eq 19a. The significance of this finding is that the limiting behavior of a high convective cooling system (eq 19a) shows the same numerical value as the temperature difference (eq 19c). 3.2. Electrostatic Model Solution. The simplified version of the Poisson-Boltzmann equation, e.g., eq 8,
6204
Ind. Eng. Chem. Res., Vol. 44, No. 16, 2005
under boundary conditions given by eq 9, has been studied, and its solution has been reported elsewhere.15-17 The analytical solution for this system is given by ψ(ξ) ) I0(λξ){ψ1K0(λ) - ψ2K0(λb)} + K0(λξ){ψ2I0(λb) - ψ1I0(λ)} I0(λb)K0(λ) - I0(λ)K0(λb)
(20a) where I0 and K0 are the modified Bessel function of the first and second kind, respectively. Equation 20a is only valid in the range -1 e ψ e +1 as a direct consequence of the Debye-Hu¨ckel approximation; however, using the strategy proposed by Oyanader and Arce,16 this range can be extended to, at least, -4 e ψ e +4, where more practical electrostatic potential values could be found.9,17 The associated and corrected explicit solution to eq 8 is given by
ψ(ξ) )
+ I0(λfAOb)K0(λfAO) - I0(λfAO)K0(λfAOb) K0(λfAOξ){ψ2I0(λfAOb) - ψ1I0(λfAO)}
(20b)
where λ is modified by the correction function ×c4AO16 not only to extend the valid range but also to obtain more accurate ψ(ξ) values. The electrostatic force term in the Navier-Stokes equation (see the next section) can now be computed using eq 20b, which only invokes Bessel functions. 3.3. Hydrodynamic Model Solution. The dimensionless Navier-Stokes equation, eq 15, can be integrated after the temperature function θ(ξ) has been replaced by eq 18, which is valid for the fluid part, and the electrostatic potential function, ψ(ξ), has been replaced by eq 20b. The solution, after application of the boundary condition eq 16, yields
A0 2 A1 A2 (ξ - 1) + (ξ4 - 1) - (2ξ2 ln|ξ| 4 16 8 ψ12 - 1 ln|ξ| 1 ψ*(ξ) + ξ2 + 1) + A3 ln|ξ| + 2 1-b 1 - b2 ln|b| (21a)
{
} {
}
where the following parameters and expressions have been identified in the function above:
[
]
2
φe φ2 Gr A0 ) Pm sin(R) (1 - b2) + b + θmin - θm Re 4 2 (21b) A1 )
A2 )
A3 )
(
Gr φ2 sin(R) Re 4
(
(21c)
)
φ2 2 Gr sin(R) b - φe2b2 Re 2
(21d)
A1 A0 (1 - b2) + (1 - b4) + 4 16 A2 2 1 (2b ln|b| - b2 + 1) (21e) 8 ln|b|
)
(21f)
)
φ2 2 b - φe2b2 - C2Pm + C3 - 16C4 2 + θm ) Gr sin(R) C2 Re 2 2 φe φ 1 φ2 φ2 2 (1 - b2) + b+ - b + φe2b2 (22) 4 2 Nu 2 2 C0φ2 + C1
I0(λfAOb)K0(λfAO) - I0(λfAO)K0(λfAOb)
ψ1 ψ(ξ) ; ψ12 ) ψ2 ψ2
The parameters A0 and A2 are closely related to the buoyancy-driven term, affected by the Joule and resistive heating effect, and are porosity dependent. Only A0 is related to the pressure-driven term. On the other hand, A1 is affected by resistive heating and buoyancy forces. Finally, A3 is a combination of the previous parameters and, consequently, inherits their characteristics driven by porosity. Some of the qualitative and semiqualitative information about the flow given by the hydrodynamic velocity profile is analyzed in the next section. 3.4. Solution of the Mass Conservation Condition. The main objective in the solution of eq 12 is to obtain an expression for the system mean temperature, Tm, or its equivalent dimensionless form, θm. The expression in this case corresponds to
(
I0(λfAOξ){ψ1K0(λfAO) - ψ2K0(λfAOb)}
Vx+(ξ) )
ψ*(ξ) )
[
}]
{
where the following parameters have been isolated in the function above:
C0 ) -
C1 )
1 Gr sin(R) (4b6 ln|b| + 3b2 - 3b6 - 3 + 48 Re 1 3b4 - 4 ln|b|) (23a) ln|b|
1 Gr sin(R) (4b4(ln|b|)2 - 5b4 ln|b| + 4b2 ln|b| + 4 Re 1 2b4 - 4b2 + 2 + ln|b|) (23b) ln|b|
C2 ) {b4 ln|b| - b4 + 2b2 - 1 - 2 ln|b|}
{
C3 ) 2b2 - 2 +
1 ln|b|
(23c)
}
ψ12 - 1 4 (2b2 ln|b| - b2 + 1) ln|b| 1 - b2 (23d)
C4 ) {[I0(λb) - ψ12I0(λ)][bK1(λb) - K1(λ)] + [ψ12K0(λ) - K0(λb)][bI1(λb) - I1(λ)]}/[(I0(λb)K0(λ) I0(λ)K0(λb))(1 - b2)λ] (23e) In eqs 22 and 23, it is observed that θm is a linear function of the parameter Pm, the Reynolds number Re, the Joule heating generation f 2, and the resistive heating generation φe2. In contrast, θm is a hyperbolic function of both the Nusselt, Nu, and the Grashof, Gr, numbers, as well as of the inclination angle, R. The inverse dimensionless Debye length, λ, and the core radius, b, show nonlinear relationships. 4. Development of New Design Criteria To compute meaningful numerical values for the temperature and velocity profiles, a number of new
Ind. Eng. Chem. Res., Vol. 44, No. 16, 2005 6205
criteria must be developed. This is first accomplished by restricting the range of the feasible Joule heating generation number, φ2, and of the resistive heating number, φe2, to those values that, in combination, do not imply a change in the fluid phase. Under ambient conditions, a change of phase will occur at temperature values of approximately θ g 0.28, yielding a range of the Joule heating generation number of 0 e φ 2 e 0.9 and of the resistive heating number of 0 e φe2 e 1.0. Another physical constraint, taken into consideration, is that the mean temperature cannot exceed the maximum temperature in the system. By using this physical constraint, a range of feasible values for the dimensionless hydraulic head gradient has been derived. In particular, eqs 24a and 24b are used to identify the location of the maximum temperature and its value, respectively. The described criteria are implicit in all the calculations herein presented.
ξmax ) b θmax )
{
φ2 (1 - b2) + 4
}
(24a)
{
}
φ2 2 1 φ2 φ2 2 b - φe2b2 ln|b| + - b + φe2b2 2 Nu 2 2
(24b)
In addition, the minimum temperature is located at the external surface, i.e., ξmin ) 1, and therefore, eq 19b gives the θmin. By using the expressions for the maximum and minimum temperatures in the system, eqs 24b and 25b, respectively, in combination with the expression for the mean temperature, eq 22, the following dimensionless hydraulic head gradient criteria are established.
Pm,min )
(
)
φ2 2 b - φe2b2 + C3 - 16C4 2 + C2
C0φ2 + C1
{
(
) | |}
φe2 2 Gr φ2 2 sin(R) b b - φe2b2 ln b Re 2 2
Pm,max )
(
C0φ2 + C1
)
(25)
φ2 2 b - φe2b2 + C3 - 16C4 2 + C2
{
}
φe2 2 φ2 Gr sin(R) b + (1 - b2) Re 2 4
(26)
The use of these design equations is very useful to identify parameter values for the numerical illustration that is presented in the next section. Once again, the criteria represented by eqs 27 and 28 are implicit in all the calculations herein presented. 5. Ilustrative Results and Discussion This section particularly focuses on the use of the developed criteria (see Section 4) in combination with the model solutions (see Section 3) to obtain the temperature and the hydrodynamics velocity profiles. Several case scenarios of the main system variables have been portrayed for analysis and discussion. Figure 2 shows the temperature profile for different combination of two cases of the resistive heating number (φe2), two cases of the core radius number (b), and two cases of the Nusselt number (Nu) with the Joule heating
number (φ2) as a parameter ranging between the values of 0.05 and 0.90. A general view of Figure 2 suggests that, because of the symmetrical characteristics of the system, the locations of the maximum temperature values are found at the center region of the annular capillary. In addition, an increase in the Joule heating parameter yields a more pronounced parabolic type curve, flattened at the top when no resistive heating occurs. For example at center locations ξ ) 0, an incremental increase of φ 2 from 0.25 to 0.70 produces an increase in temperature of approximately 180%. At the same location, an incremental increase of φ 2 from 0.25 to 0.50 produces an increase in temperature of approximately 100%, and an incremental increase of φ2 from 0.05 to 0.25 increases the temperature by 400%. This observed effect indicates that a small variation on the Joule heating parameter may cause an important impact on the system temperature, especially in the lower range of φ 2. When comparing parts a and b of Figure 2 (or parts e and f of Figure 2), one can clearly observe the important role of porosity in temperature development. When the core presents no resistive heating, a lower porosity reduces the effect of joule heating and, therefore, a flattened-at-the-top temperature profile is developed. On the contrary, when the core does present resistive heating as shown in Figure 2c,d, a steeper temperature profile is observed, dominated by resistive heating. A less-sharp profile is observed at the center of the annular core for low φ 2 values where only resistive heating dominates the process. In all these cases, the effect of a lower porosity amplifies the resistive heating effect just described. When comparing parts a and e of Figure 2 (or parts b and f of Figure 2), one can see that the role of the heat exchange forces is clearly working; for example, in the case of the smaller Nusselt values, there is a “lifting” effect on all the temperature values, equally distributed along the radial position ξ. To illustrate the roles of Joule heating, resistive heating, porosity, buoyancy, electroosmosis, and pressure-driven forces on the velocity profile, six sets of figures are presented and analyzed. On all these plots, the inverse dimensionless Debye length is used as a parameter between λ ) 1 and λ ) 20. This particular range covers the most typical values of inverse dimensionless Debye lengths that have effects on the electrostatic potential and, therefore, on velocity. Parts a, b, and c of Figure 3 show velocity profiles for three cases of the Joule heating parameter, respectively, with inverse dimensionless Debye length as a parameter. On this set of plots, the Grashof number, Gr, has been held at 1 and the electrical potential ratio, Ψ12, has been held at 0.8. These particular values are the most illustrative and representative of the phenomenon under study. In Figure 3a, it is clearly observed that, at low values of φ 2 (0.05), electroosmosis dominates, controlling the flow for all λ values. In particular, as λ values increase from 1 to 20, the location of the maximum punctual velocity moves from a position near the inner wall to a position near the outer wall. In the next two plots, parts b and c of Figure 3, where the Joule heating parameter is gradually increased, buoyancy forces begin to compete, changing the direction of the flow at the midway point of the annular region. This phenomenon causes clearly observed flow-reversal regimes at higher magnitudes of the inverse dimensionless Debye length; however, this trend tends to vanish at
Figure 2. Dimensionless temperature profiles (inside the annular capillary channel) for various values of the heat generation parameter. Results include cases of Nusselt number Nu ) 10 (cooling system), with no resistive heating φe2 ) 0, (a) b ) 0.2 and (b) b ) 0.4; with resistive heating φe2 ) 0.7, (c) b ) 0.2 and (d) b ) 0.4. Results also include cases of Nusselt number Nu ) 5 (less cooling system), with no resistive heating φe2 ) 0, (e) b ) 0.2 and (f) b ) 0.4.
6206 Ind. Eng. Chem. Res., Vol. 44, No. 16, 2005
Ind. Eng. Chem. Res., Vol. 44, No. 16, 2005 6207
Figure 3. Typical hydrodynamic behavior of the physical system for Grashof number Gr ) 1, electrical field in the axis direction, electrical potential ratio Ψ12 ) 0.8, and no resistive heating effect φe2 ) 0. Dimensionless axial velocity profiles (inside the annular capillary channel) are shown for different values of the Joule heating number, (a) φ 2 ) 0.05; (b) φ 2 ) 0.5; and (c) φ 2 ) 0.9. In all three situations, the dimensionless inverse Debye-Hu¨ckel length, λ, is used as a parameter with the values indicated in the figure.
higher Joule heating parameters, causing buoyancydriven flows.
Figure 4. Typical hydrodynamic behavior of the physical system for Grashof number Gr ) 1, electrical field in the axis direction, electrical potential ratio Ψ12 ) 0.8, and resistive heating effect φe2 ) 1.0. Dimensionless axial velocity profiles (inside the annular capillary channel) are shown for different values of the Joule heating number, (a) φ 2 ) 0.05; (b) φ 2 ) 0.5; and (c) φ 2 ) 0.9. In all three situations, the dimensionless inverse Debye-Hu¨ckel length, λ, is used as a parameter with the values indicated in the figure.
A variation of Figure 3 is Figure 4, which describes the same variables; however, resistive heating has been
6208
Ind. Eng. Chem. Res., Vol. 44, No. 16, 2005
introduced and held at 1. In general, the set of figures allows the verification of the influence of the heat generation due to the electrical resistivity of the core. The overall trend is very similar to that of Figure 3; however, the main difference is that all the locations of maximum punctual velocities have moved to a position near the inner wall. This clearly indicates that buoyancy is the dominating force, because higher temperature values are developed in the inner wall region. Figure 5 also illustrates the same variables as in Figure 3; however, resistive heating has been introduced and held at 1, as in Figure 4. Additionally, the electrical potential ratio, Ψ12, has been held high (1.8). By stressing this last parameter, in general terms, the moving effect of the maximum punctual velocity location to the inner wall is amplified. In this case, buoyancy and electroosmotic forces collaborate to produce a “lifting” effect on all the positive velocity values, starting at the inner wall and gradually vanishing toward the outer wall. On the contrary, the negative velocity values are uniformly amplified, yielding a pressure-driven type of flow. Overall, the occurrence of flow-reversal regimes is more evident in Figure 5 than in the two previous figures. For example, Figure 5a presents a flow-reversal regime for λ ) 1, which was not observed before. Additionally, Figure 5b,c shows amplified versions of the flow reversal identified in Figure 3b,c. In another view, parts a, b, and c of Figure 6 illustrate velocity profiles for three cases of Joule heating parameters, respectively, with inverse dimensionless Debye length as a parameter. On this set of plots, the Grashof number, Gr, has been held at 1 and the electrical potential ratio, Ψ12, has been held high (1.8). No resistive heating effect is considered, and the porosity has decreased (b ) 0.4). Comparing parts a, b, and c of Figure 6 for low λ values (
