ARTICLES
A solution for a heat transfer model in a moving bed through the self-adjoint operator method.
H. F. Meier, D. Noriler and S. L. Bertoloi
Department of Chemical Engineering, Regional University of Blumenau, São Paulo Street 3250, Zip Code 89030-000, Blumenau, Brazil
meier@furb.br
Abstract - Usually, heat and/or mass transfer models with time dependence, in a fixed, moving or cross-flow beds, are solved analytically by the use of the Laplace transform method. When the determination of the character of the poles is not on easy problem, this method presents the transform inversion using the residue theorem as the major application difficulty. In this work, an alternative method is discussed which casts the system of equations into a matrix problem of the Sturm-Liouville type. As an example, the solution of a heat transfer model in a moving bed is presented. The advantage this approach is a direct solution of the temperature profiles in the particle and in the bulk fluid near the solid-fluid interface by using a spectral expansion in terms of the self-adjoint matrix operator involved, with guaranteed convergency, and it can be used easily as an interpolation scheme to solve numerically advection/diffusion problems.
Keywords - Heat Transfer; Moving Bed; Multiphase Reactors; Analytical Solution; Selfadjoint.
I. GENERAL CONSIDERATIONS
As for the heat and/or mass transfer between a particulated solid phase and one or more fluid phases (gaseous or liquid), a great variety of technological applications is observed, mainly due to the high eciency of those transfer processes. We can, for example, emphasize the pyrolysis of oil shale fine particles in a moving bed reactor (Bertoli, 2000), the Fischer Tropsch synthesis in well-mixed slurry reactors (Ahn et al., 2005), the cross-flow continuous fluidized bed dryer (Izadifar and Mowla, 2003), the catalytic cracking of petroleum in a circulating fluidized bed reactor (Michalopoulos et al., 2001; Becerril et al., 2004), and others.
The mathematical modeling of the above exemplified transfer processes follows basically models of two or more phases that are related by interphase transfer processes; in these models, the spatial variation of the variables are considered in the distributed parameter model or neglected by the lumped parameter model. Real processes are distributed, but some of them can be approximated as lumped ones (in one or more phases). A situation in which the temperature gradient in the solid phase may be neglected will be commented in the following.
Bertoli (1989, 2000) presents a comparison between a lumped parameter model and a distributed parameter one, for the heat transfer process in a moving bed, based on experimental data obtained by Lisbôa (1987). In that analysis, the analytical solutions were obtained by the Laplace transform method with inversion by the residue theorem, and the results for both models proved a good concordance between the experimental data and those predicted by the model. We intend to demonstrate the application of an alternative method for the solution of partial differential equation systems from self-adjoint operators, that permits a maximum preservation of the physical characteristics of the problem without, however, needing excessive numerical calculations. This method was mathematically detailed by Arce and Ramkrishna (1986, 1988) for the heat transfer in a fixed bed reactor.
The motivation to apply analytical methods for this kind of advection/diffusion problem resumes, in fact, in two aspects:
the first one in a mathematical point of view, the analytical solution can be easily extended for numerical methods as interpolation scheme in the CFD techniques to increase convergence rate and stability of the numerical solutions;
in a physical point of view, this heat transfer model represents the interface between gas-solid phases and it can be applied for several physical situations.
II. MATHEMATICAL MODELING
Figure 1 illustrates a differential volume element of the reactor on which the energy conservation equation will be applied.
Figure 1: Volume element of the moving bed reactor.
The simplifying hypotheses adopted in the model are the following
the velocities and physical properties of the fluid and the particles are constant and uniform;
the temperature of the internal reactor wall remains constant along the axial axis, characterizing an infinite heat reservoir;
the dragging gas can be considered as transparent to radiation;
the heat flux term, due to radiation mechanism, from the reactor wall to the particles, defined constitutively by the Stefan-Boltzmann equation, is linearized for an expression of the type:
where,
the particles are geometrically characterized as uniform spheres, by adopting a particle size distribution model and consequently by determining the Sauter average diameter;
the heat dispersion mechanism in both, radial and axial directions, may be neglected;
the form factor for radiation, , equals 1 (one);
the heat due the chemical reaction is negligible in comparison with the overall demand of the process.
From the energy balance to the solid phase, considering the Fourier equation as constitutive for the heat diffusion, one obtains the following expression for the spheric particles:
![]() |
(1) |
with
For the fluid, having in mind the advective, accumulative, and interface heat transfer terms, the energy equation may be expressed as follows:
![]() |
(2) |
with
The boundary conditions for the Eqs. 1 and 2 are:
![]() |
(3) |
And the initial conditions:
![]() |
(4) |
Introducing the residence time concept and the following non-dimensional parameters, as a manner of reducing the involved parameter number:
![]() |
(5) |
Besides, considering the reactor operating in the steady state, the Eqs. 1 and 2 and their respective boundary conditions will be
![]() |
(6) |
where S is the spheric Laplacian operator:
and,
![]() |
(7) |
With respect to the physical interpretation of the non-dimensional model parameters occurs that: αp and αw quantify the effects of the particles on the fluid and of the wall on the fluid, respectively, with α being the combination of both effects; Bp and Br relate the convective heat transfer between the particle surface and the fluid, and the radiative heat transfer between the particle surface and the wall to the diffusive heat transfer in the particle interior.
It can be observed that the higher the values of αp and αw, the higher are the interactions. For diluted systems with their porosity approaching the unit (which is the present case), the interaction represented by αp is small but very significant for the solid phase.
Analyzing the limit cases with respect to parameter B, it is concluded: when B → ∞ the diffusion is the limiting transfer process. In other words, the interphase transfer is much faster than the intraphase transfer and necessarily the model must be considered for distributed parameters relative to the solid phase; when B → 0, the transfer limitation comes from the interphase effects, enabling the utilization, without great information losses, of a lumped parameter model.
III. ANALYTICAL SOLUTION
The partial differential equation system above can be written in a compact manner in the matrix form by,
![]() |
(8) |
where is a function vector with two elements:
![]() |
(9) |
L is a differential matrix expression:
![]() |
(10) |
and, a( ) is a functional defined by:
![]() |
(11) |
Note that the matrix L contains only spatial operations.
We now introduce formally a Hilbert space, adequate for the vectors
![]() |
with, E constant |
and domains for S and L.
Note that the first element is a real number, e.g. 1, because it corresponds to the functions θf , that do not depend on any spatial variable. The second one is a function of r, v(r), and corresponds to the functions θp(r) that depend on r.
Firstly, the Hilbert space is defined for v(r),
with the inner product:
Now the domain,
and in this way the operator
is attainded.
Following,
and the operator
It must be observed that S' is not self-adjoint but has an extension S that is self-adjoint in Hr, i.e.,
for any vi,vj ∈ D(S), and corresponds to the isolated particle.
Following, we define the differential matrix
where a is such a functional that:
and,
L' operates on vectors
in a Hilbert space of the direct sum type,
that has as its domain
where the last condition incorporates the interaction on the particle surface.
Finally, the operator of our system is
We now define an internal product in D(L '),
with q being a positive number that we try to determine in such a way that L becomes self-adjoint D(L'), i.e.:
![]() |
(12) |
This is possible and the explicit calculation of the Eq. 12 results in:
![]() |
(13) |
L being an self-adjoint operator in the previously defined Hilbert space, it has the following fundamental properties (Ramkrishna and Amundson, 1985):
1. A denumerable infinite number, of real and distinct eigenvalues;
2. An infinite set of orthogonal eigenvectors, one for each eigenvalue;
3. The set of eigenvectors is complete, i.e., any vector
) may be expanded as an infinite series in the eigenvectors
), the so called spectral expansion, with guaranteed convergence.
The next step will be the determination of the eigenvalues and eigenvectors of L, solving the equation,
![]() |
(14) |
since with this equation we can write the Eq. (9) like:
![]() |
(15) |
Thus, the solution of the temporal part becomes immediate,
![]() |
(16) |
An arbitrary vector can now be expanded like:
![]() |
(17) |
The Ai coefficients are easily determined using the initial condition:
![]() |
(18) |
Taking the scalar product with follows:
![]() |
(19) |
since the eigenvectors are orthogonal.
Substituting 19 in 17, the solution of the problem is:
![]() |
(20) |
Therefore, it is necessary to determine the eigenvalues λi and eigenvectors ) through the Eq. 14 and later to calculate the coefficients
.
From the Eqs. 14 and 6 results:
![]() |
(21) |
The solution of 21 is standard and can be written like:
![]() |
(22) |
Using the boundary condition in r = 1, the constant C can be determined:
![]() |
(23) |
Introducing and substituting 23 in 21 results:
![]() |
(24) |
which is the transcendental equation to determine the eigenvalue . This equation has to be solved numerically.
Also normalizing the eigenvectors )
and requiring that
we obtain:
Now,
![]() |
(25) |
Resuming the solution:
![]() |
(26) |
or,
![]() |
(27) |
with
IV. RESULTS AND ANALYSIS
To concretize the solution it is necessary to determine a set of eigenvalues with enough terms to guarantee sufficient precision for the series convergence. These eigenvalues correspond to the zeros of the Eq. 24.
Figure 2 illustrates qualitatively the behavior of the eigenvalue equation.
Figure 2: Behavior of the eigenvalue equation.
Note that the eigenvalues approach the singularities of , initially from the left, later passing to an approach from the right. It is foreseeable that for very large values of σ, the eigenvalues approach the singularities of f(σ), given by the zeros of the equation:
In face of this behavior, it was opted for a combination of methods for the numerical eigenvalues calculation. In other words: the Muller method is applied to find the singularities of f(σ) contained in the interval ; and the bisection method for the eigenvalues contained between two singularities, beginning the approach from the left singularity.
In order to verify the accuracy of the model prediction, as model validation activity, the analytical solution was compared with the experimental data in pilot scale, obtained by Lisbôa (1987) for heat transfer to a pneumatically conveyed mixture of a oil shale fine particles and air. A brief description of the experimental unit utilized by Lisbôa (1987) and more details in the experimental data may also be found in Bertoli (2000).
As the experimental information resumes to average conditions at the reactor inlet and outlet, it becomes necessary to evaluate the average volumetric temperature of the solid particles. This evaluation can be made analytically through the expression:
![]() |
(28) |
The conditions of the experiments from Lisbôa (1987) are presented in Table 1. The respective results and predicted by the model, are presented in Table 2.
Table 1: Conditions of the experiments
Table 2: Comparison between experimental values and the model prediction values.
As criterion, the relative error defined by;
was adopted, and a good concordance between the experimental values and the values predicted by the model is observed. In the calculations with the model solution of Table 2, the fluid-particle heat transfer coefficient, hp, was obtained from Kato's correlation (Kato et al., 1993) and the wall-mixture heat transfer, by the following,
![]() |
(29) |
obtained in the experimental unit (Lisbôa, 1987).
The numerical verification of the model has been done through the dependence of the series convergence with number of series terms. Evaluating the series convergence, one reaches easily by inspection and approximation of the nth term of the series, to:
While for r = 1 and r = 0 this last approach becomes,
Hence,
Consequently, it can be concluded that even under the initial condition at the particle center, where we would have the most critical situation, the convergence is guaranteed. This can be seen in Figs. 3, 4 and 5.
Figure 3: Dependence of the series convergence with the number of series terms: average particle and fluid temperature (run 6).
Figure 4: Dependence of the series convergence with number of series terms: particle profile temperature at τ = 0 (run 6).
Figure 5: Series convergence with 260 terms: particle
the number of series terms: average particle and fluid profile temperature at τ = 0.
temperature (run 6).
The Fig. 3 shows the dependence of the series convergence with the number of series terms for run 6 along the reactor height. It is observed that the series for fluid temperature and for particle average temperature converges with approximately 3 terms, even in the critical region, where τ = 0.
In Fig. 4, it can be seen the dependence of the particle temperature series convergence with number of series terms for run 6 at the inlet reactor (τ = 0). In comparison with Fig. 3, it is possible to note that the series convergence for the particle temperature is more difficult than particle average and fluid temperature. For this case, a satisfactory convergence was obtained with about 200 series terms. The Fig. 5 presents the particle temperature profile at the inlet reactor evaluated with a 260 series terms. It can be noted that in the particle surface and at the center of the particle, the convergence series is slower. Figure 6 shows the particle temperature profile for run 6, and for other parameter sets obtained through the thermal conductivity variation and/or the wall temperature, so that the hydrodynamic conditions of run 6 could be maintained.
Figure 6: Radial profile of the particle temperature for various radiation Biot numbers at the reactor outlet
Rearranged in the form of the Biot number for radiation may physically be interpreted as the ratio of the internal conduction resistance and external radiation resistance. In Fig. 6, one notes that, as Br increases, the other parameters remaining constant (including τ ≡ Fourier number = Fo), the temperature gradient in the particle interior increases with Br until reaching a maximum value in
, decreasing afterward.
This noteworthy behavior is explained based on the fact that in a moving bed heating system with a dragging fluid transparent to thermal radiation, under certain conditions, (specially when the radiation is dominant), the particles can become hotter than the fluid and an increase in Br (for example, through the Tw increase) contributes to a temperature gradient diminution.
One comes to the important conclusion that in this system, with the other parameters kept constant, (also Fo), different numbers of compounded Biot, B = Bp + Br, of the same order of magnitude, can not be taken separately as a form for the temperature gradient classification in the particle interior. An approximate criterion to determine this temperature inversion is developed in Appendix I for B < < 1, at the limit of long time scales.
Applying the criterion for αp = 0.0668, αw = 0.0087, and Bp = 0.0037, one obtains:
= 0.0049,
Corresponding approximately to what is observed in Fig. 6, despite that in this case (i. e., Br = 0.0050) τ = 0.742, being under the value required for the criterion application, which, in this system is.
The axial distribution of the fluid temperature and the average particle temperature can be visualized through Fig. 7. In this case, the average particle temperature variation is approximately linear with the reactor height, what is explained by the fact that the total heat rate transferred to the particles, (qT )PS, was approximately constant along the course. At the beginning, (qC)PS was greater due to the quick heating of the fluid with relation to the particles, but approximately from 0.3 meters, the radiant heat rate (qR)PS increases, compensating the decrease in (qC)PS (an analysis of the relative radiant rate variation with the reactor length is presented in Bertoli (2000).
Figure 7: Axial profile of the gas and average particle temperature for run 6.
V. CONCLUSIONS
While examining Table 2 it is clearly noted that the experimental data confirm the mathematical model adopted in this work for the heat transfer prediction in a moving bed reactor/exchanger. Similar results were also obtained by Bertoli (1989, 2000) through different analytical techniques. Particularly, the reader is referred to a comparison with the lumped capacity model solution of Bertoli (2000).
Among the results of Section IV, it was possible to conclude that in a moving bed heating system with the dragging fluid transparent to thermal radiation, different numbers of compounded Biot, B, of the same order of magnitude, cannot be taken separately as a way for the classification of the temperature gradient magnitude in the particle interior for a same Fo, while keeping constant the rest of the parameters.
The method of self-adjoint operators, which enables an alternative solution for the problem, has shown to be highly applicable in view of the few analytical and numerical calculations necessary for a solution consolidation, a characteristic that is not observable when the Laplace transform method is used.
The use of the same method for the analytical solution of non homogeneous problems, compound of three phases (dragging fluid, particles and fluid of external heating), in conjunction with the generalized integral transform method, is seen as possible.
These solutions contribute in an appreciable way to the understanding of heat transfer phenomena by the fact that the operators preserve in their proper formulation aspects related to the physical reality.
As a final comment, we note that recently Sazhin et al. (2004) have obtained an analytical solution for a model that has allowance for a variable fluid-particle heat transfer coefficient (i. e. hp = hp(t)). However, the attained solution does not present an equation for the temperature evolution of the coupled fluid phase, but it must be determined by other means such as, for example, by computational fluid dynamics techniques.
1The authors express their gratitude to the Department of Chemical Engineering of Regional University of Blumenau, which enabled that this work could be carried out and also are very indebted to Prof. Gerardus L. M. Dassen (deceased) who inspired the present work.
APPENDIX
In this appendix, criteria are developed to evaluate the Br value at which the particle surface temperature, θp(1, τ), becomes equal to the fluid temperature, θf(τ), for systems with B < < 1 and long time scales. With this purpose, approximations for the roots (i. e., for the characteristic values, σi) and for the singularities of Eq. 24, respectively, are obtained in sections A and B.
A. Approximation for the roots of Eq. 24 when B < < 1
Rewriting Eq. 24 in the way,
![]() |
(30) |
and approximating B - 1 ≈ -1; Br - 1 ≈ -1 for Equation 30) results,
![]() |
(31) |
The roots of Eq. 31 are the approximations of the characteristics values:
Since σ = 0 is not the root of Eq. 24, a root σi = ε is sought, so that ε be a small number and, therefore,
![]() |
(32) |
Replacing these approximations in Eq. 24 abolishing the term of the 4th order, one has,
![]() |
(33) |
B. Approximation for the singularities of Eq. 24 when B < < 1
The singularities of Eq. 24 are the singularities of f(σ) and are in the roots of
![]() |
(34) |
Approximating sin(σ) ≈ σ and , one obtains for the first singularity,
![]() |
(35) |
and the remaining singularities can be approximated by roots of
and, therefore,
Comparing with the results of Section A.1, one notes a trend of the characteristic values to approximate the singularities as can be seen in Fig. 2; in other hand, once the characteristic values are among the singularities, one concludes that:
![]() |
(36) |
In this way, if B < < α < 1 is possible to approximate , the following ordination is obtained:
![]() |
(37) |
C. Criteria for the equality between θp(1,τ) and θf(τ) for systems with B < < 1 and long time scales
Of Eqs. 24, 27,
![]() |
(38) |
Then, admitting that the time scale is sufficiently long as to make acceptable the approximation of Eq. series (37) by the first term (i. e., and
< <
, a situation that, according to Section A.2, occurs, for example, when α > > B and
then the conditions for the temperature equality are:
(a)
![]() |
(39) |
(b)
![]() |
(40) |
To obtain these results, the approximation for σ1 (i. e., σ1 ≈ ε) was used, as well as sin(σ1); , besides, one admitted I = 1 (i. e., Tfi = Tpi.
If αw = αp, both criteria transformed into
![]() |
(41) |
REFERENCES
1. Ahn, V.R., E.F. Costa Jr., J.E.P. Monteagudo, C.E. Fontes, E.C. Biscaia Jr. and P.L.C. Lage , "A comprehensive mathematical model for the Fischer-Tropsch synthesis in well-mixed slurry reactors," Chemical Engineering Science, 60, 677-694 (2005).
2. Arce, P. and D. Ramkrishna, "Self-adjoint operators of transport in interacting solid-fluid sistems: I," Chemical Engineering Science, 41, 1539-1547 (1986).
3. Arce, P. and D. Ramkrishna, "Self-adjoint operators of transport in interacting solid-fluid sistems: II," Chemical Engineering Science, 43, 933-944 (1988).
4. Becerril, E.L., R.M. Yescas and D.S. Sotelo, "Effect of modelling pressure gradient in the simulation of industrial FCC risers," Chemical Engineering Journal, 100, 181-186 (2004).
5. Bertoli, S.L., Transferência de calor convectiva e radiante em leito de arrasto, M.Sc. thesis, COPPE, Federal University of Rio de Janeiro (1989).
6. Bertoli, S.L., "Radiant and convective heat transfer on pneumatic transport of particles: an analytical study," Int. Journal of Heat and Mass Transfer, 43, 2345-2363 (2000).
7. Izadifar, M. and D. Mowla,"Simulation of a cross-flow continuous fluidized bed dryer," Journal of Food Engineering, 58, 325-329 (2003).
8. Lisbôa, A.C.L., Transferência de calor em leito de arrasto de xisto, M.Sc. thesis, COPPE, Federal University of Rio de Janeiro (1987).
9. Michalopoulos, J., S. Papadokonstadakis, G. Arampatzis and A. Lygeros, "Modeling of an industrial fluid catalitic cracking unit using neural networks," Trans IChemE, 79, 137-142 (2001).
10. Ramkrishna, D. and N.R. Amundson, Linear operator methods in chemical engineering with applications to transport and chemical reaction sistems, Prentice-Hall, Englewood Cliffs, (1985).
11. Kato, K., L. Onozawa and Y. Noguchi, "Gas particle heat transfer in a dispersed bed," Journal of Chemical Engineering of Japan, 16, 178-182 (1983).
12. Sazhin, S.S., P.A. Krutitskii, W.A. Abdelghaffar, E.M. Sazhina, S.V. Mikhalovsky, S.T. Meikle and M.R. Heikal, "Transient heating of diesel fuel droplet," Int. Journal of Heat and Mass Transfer, 47, 3327- 3340 (2004).
Received: February 27, 2008
Accepted: November 4, 2008
Recommended by Subject Editor: Walter Ambrosini