INTEGRAL TRANSFORM METHOD FOR LAMINAR HEAT TRANSFER CONVECTION OF HERSCHEL-BULKLEY FLUIDS WITHIN CONCENTRIC ANNULAR DUCTS

a01v18n4

INTEGRAL TRANSFORM METHOD FOR LAMINAR HEAT TRANSFER CONVECTION OF HERSCHEL-BULKLEY FLUIDS WITHIN CONCENTRIC ANNULAR DUCTS

M.J.G.Viana¹, U.C.S.Nascimento¹, J.N.N.Quaresma¹ and E.N.Macêdo2

1Chemical Engineering Department, CT, Universidade Federal do Pará, UFPA,

Campus Universitário do Guamá. Rua Augusto Corrêa 01,

66075-900, Belém - PA, Brazil

E-mail: [e-mail]

2Mechanical Engineering Department - CT, Universidade Federal do Pará, UFPA,

Campus Universitário do Guamá. Rua Augusto Corrêa 01,

66075-900, Belém - PA, Brazil

(Received: May 16, 2001 ; Accepted: August 18, 2001)

Abstract - Related momentum and energy equations describing the heat and fluid flow of Herschel-Bulkley fluids within concentric annular ducts are analytically solved using the classical integral transform technique, which permits accurate determination of parameters of practical interest in engineering such as friction factors and Nusselt numbers for the duct length. In analyzing the problem, thermally developing flow is assumed and the duct walls are subjected to boundary conditions of first kind. Results are computed for the velocity and temperature fields as well as for the parameters cited above with different power-law indices, yield numbers and aspect ratios. Comparisons are also made with previous work available in the literature, providing direct validation of the results and showing that they are consistent.

Keywords: viscoplastic materials, concentric annular ducts, integral transform technique.

INTRODUCTION

Concentric annular ducts are important geometric configurations in the design of many fluid-flow and heat-transfer devices, the most notorious example of which is the double-pipe heat exchanger, consisting of two concentric circular tubes. In industry, several manufacturing processes such as extrusion, drawing, hot rolling, and transport of rock cuttings during the drilling of oil and gas wells involve flows of viscoplastic materials in annular geometries. The main characteristic of these materials is the existence of a yield stress, below which they do not flow. This non-Newtonian behavior is found in various materials, such as solutions or melts of polymeric materials, oils, greases, cosmetics, toothpaste, soap, detergents, paints, cement and drilling muds. The majority of these materials, especially the drilling muds, follow the Herschel-Bulkley rheological model and their rheological properties are extremely sensitive to the temperature field.

Through the years, few studies dealing with the heat and fluid flow of viscoplastic materials within concentric annular ducts have been available in the literature, despite their important industrial applications such as those listed above. Among these studies can be included the work of Laird (1957), Fredrickson and Bird (1958), and Hanks and Larsen (1979), in which the authors have determined simple equations to determine the complete characteristics of the fluid flow. On the other hand, work dealing with heat transfer problems involving viscoplastic materials, such as Herschel-Bulkley fluids, are mostly concerned with circular and parallel-plate ducts geometry. For instance, Blackwell (1985) presented a numerical solution for the Graetz problem in laminar flow of Bingham fluids within circular ducts subjected to a boundary condition of the first kind.

An experimental and theoretical heat transfer study for Herschel-Bulkley fluids was conducted by Nouar et al. (1994). In their analysis, thermally developing flow in circular ducts was assumed, and an interesting discussion regarding the influence of temperature-dependent rheological properties on the velocity profiles and Nusselt numbers can be found. Later, in another article, Nouar et al. (1995) numerically analyzed the same problem, this time assuming simultaneously developing flow, and presented some correlations for local Nusselt number and pressure gradient.

Mendes and Naccache (1995) analyzed the convective heat transfer problem in laminar flow of Herschel-Bulkley fluids within circular ducts. In this work, local Nusselt numbers were obtained as a function of the yield stress and power-law index. In a similar way, Soares et al. (1996, 1997) studied the heat transfer problem in laminar flow of Herschel-Bulkley fluids through the entrance region of circular ducts. The authors showed that thermal entry length decreases when the behavior of the rheological fluid deviates from the Newtonian case.

The thermal entry region in laminar forced convection of Herschel-Bulkley fluids was analytically solved by Quaresma and Macêdo (1998) for both circular and parallel-plate ducts, maintained at either prescribed wall temperature or prescribed wall heat flux, using the integral transform technique. Highly accurate results were obtained for local Nusselt numbers for different power-law indices and yield numbers. Recently, Nascimento et al. (2000) employed the same approach to analyze the thermal entry region in laminar flow of Bingham plastics within concentric annular ducts by adopting boundary conditions of the first kind. Nusselt numbers were calculated for both the thermal entry and fully developed regions as a function of yield numbers and aspect ratios.

Therefore, in this context, to fill in the gap in the literature for heat and fluid flow of Herschel-Bulkley fluids within concentric annular ducts, the present work aims at advancing the ideas in the classical integral transform technique and the so-called sign-count approach to determine the product of the Fanning friction factor-apparent Reynolds number and velocity profiles. It also aims at accurately computing the Nusselt numbers for both the thermal entry and fully developed regions within the range of parameters analyzed, i.e., power-law indices, yield numbers and aspect ratios. Comparisons with previous work in the literature are also made for typical situations in order to validate the numerical code developed here and to demonstrate the consistency of results produced.

ANALYSIS

The problem is geometrically defined by two circular concentric ducts between which a time independent non-Newtonian fluid flows, obeying the Herschel-Bulkley model for shear stress, according to Figure 1. The constitutive equation to describe the rheological behavior of this non-Newtonian fluid is given in the following form:

where, t_rz is the shear stress, t_o is the yield stress, is the shear rate, u is the velocity component in the axial direction z, K is the consistency index of the fluid and n is the power-law index which is less than unity for pseudoplastic fluids (shear-thinning materials) and greater than unity for dilatant fluids (shear-thickening materials). The positive sign for t₀ is used when the transport of momentum is in the positive r-direction, and in this case = -du/dr. When the momentum is being transported in the negative r-direction, t₀ is negative, while = -du/dr.

For the fully developed region of a concentric annular duct, the momentum equation for axial coordinate z is simplified to yield

subject to the following boundary conditions:

Then, substituting Equations (1) into Equation (1.a) with the appropriate signs for the yield stress and shear rate, and after the integrations are performed, the fully developed velocity profiles in the three distinct regions for a Herschel-Bulkley fluid are given by

for r_iw < r < a^*

for a^* < r < b^*

for b^* < r < r_ow

where a^* and b^* are the bounds on the plug-flow region and c^* is the value of the radial coordinate where the shear stress is zero. The velocity profile for the plug-flow region, u^o(r), is given by either Equation (3.b) or (3.c); these equations are obtained by making r = a^* in Equation (3.a) and r = b^* in Equation (3.d), respectively, and they give the same values u^o(r).

The velocity profiles given by Equations (3) for a Herschel-Bulkley fluid are split into three distinct regions. One is a^* < r < b^*, which denotes the plug-flow region where |t_rz|< t_o and the fluid behaves like a solid plug, and the other two regions are r_iw < r < a^* and b^* < r < r_ow, where |t_rz| ³ t_o and refers to those parts of the fluid which are in shear flow.

Now, the following dimensionless groups are defined:

where D_h = 2(r_ow - r_iw) is the hydraulic diameter, Y is the yield number, f is the Fanning friction factor and Re_a is the apparent Reynolds number.

Substituting the groups defined above by Equations (4) into the fully developed velocity profiles given by Equations (3) results in

The dimensionless bounds on the plug-flow region, a and b, and the dimensionless radial position where the shear stress is zero, c, are obtained by making the integration of Equation (2.a) resulting in a linear profile for the shear stress, which is evaluated at the desired points r = a^*, b^* and c^*, after that the dimensionless groups given by Equations (4) are introduced to yield the following equations:

The average flow velocity is defined as

which in dimensionless form is written as

In the analysis of the thermal problem of a non-Newtonian fluid that obeys the Herschel-Bulkley model, described by Equations (1), one assumes steady forced convection in thermally developing, hydrodynamically developed laminar flow, given by Equations (5), within a concentric annular duct, also according to Figure 1. Viscous dissipation, free convection, and axial conduction effects are neglected, and physical properties are assumed to be constant. The duct walls are subject to boundary conditions of the first kind at either the inner or outer duct wall, referred to as Case A (T_iw = T_e) and Case B (T_ow = T_e), and the fluid enters the duct at a uniform temperature, T_e.

The mathematical formulation in dimensionless form for this forced convection problem is given by

in g < R < 1, Z > 0

subject to the following inlet and boundary conditions:

where coefficient m in boundary conditions (10.c,d) identifies Case A (m = 1) or Case B(m = 0).

In Equations (10) above, the following additional dimensionless groups were used:

The problem defined by Equations (10) can be readily solved by the classical integral transform technique (Mikhailov and Özisik, 1984; Cotta, 1993). However, in order to make the boundary conditions (10.c,d) homogeneous and to obtain a convergence acceleration of the final solution in such a manner that a better computational performance in the series expansion can be obtained, the so-called splitting-up procedure (Mikhailov and Özisik, 1984; Mikhailov, 1977) is applied to this problem. Thus, a general separation into simpler problems is proposed in the following form:

In Equation (12), q_p(R) represents the separate solution for the nonhomogeneous boundary conditions (10.c,d) and q_h(R,Z) is the homogeneous solution of problem (10). After substituting Equation (12) into Equations (10), these intermediate potentials are obtained from the following formulations:

which is readily integrated to give

and the general homogeneous problem is given by

in g < R < 1, Z > 0

with the following inlet and boundary conditions:

The homogeneous problem defined by Equations (15) can also be solved by the classical integral transform technique. Then, following the procedures in the application of this technique, the appropriate auxiliary eigenvalue problem needed for its solution is

in g < R < 1

where y_i(R) and m_i are the eigenfunctions and eigenvalues, respectively. The problem defined by Equations (16) is solved by the so-called sign-count method (Mikhailov and Özisik, 1984; Mikhailov and Vulchanov, 1983), which offers safe and automatic accurate computations of as many eigenvalues and eigenfunctions as desired. The eigenvalue problem allows for the development of the following integral transform pair:

transform

where Ni, the normalization integral, is given by:

Now, Equations (15) can be integral transformed with y_i(R)dR to yield the following ordinary differential equation for the transformed potential, (Z):

with the transformed potential inlet condition given by

The solution of the transformed potential, (Z), is readily obtained as

Therefore, substituting Equation(20) into the inversion formula (17.b), the solution for q_h(R,Z) is determined in the following form:

Thus, Equation (21) for q_h(R,Z) in conjunction with Equation (14) for q_p(R) complete the solution for the potential q(R,Z), originally formulated in Equation (10), and this solution is written as:

At this point, thermal parameters of practical interest, such as average temperature and local Nusselt numbers can be obtained from Equation (22) as follows:

The average flow temperature, q_av(Z), is defined as:

Thus, substituting Equation (22) into Equation (23), q_av(Z) is readily obtained in the following form:

where

The local Nusselt numbers at the inner and outer duct walls are defined as

By taking the derivative of Equation (22) and evaluating it at R = g and R = 1 with Equation (24) for the average flow temperature, q_av(Z), the solutions for the local Nusselt numbers at the inner and outer duct walls can be readily determined from Equations (26) and (27), respectively, in the following form:

RESULTS AND DISCUSSION

Numerical results for the product of the Fanning friction factor-apparent Reynolds number and velocity profiles were calculated for different values of power-law indices, yield numbers and aspect ratios. Once these quantities had been computed, the temperature field was readily obtained from Equation (22), and consequently thermal parameters such as average temperature and local Nusselt numbers at the inner and outer duct walls were obtained from Equations (24) to (27), respectively. For this purpose, computational codes were developed in the FORTRAN programming language and implemented on a PENTIUM II 400 MHz computer platform.

The Fanning friction factor and the dimensionless bounds on the plug-flow region are functions of the power-law indices and apparent Reynolds and yield numbers. These quantities are determined by satisfying Equation (9.a) for the average flow velocity and by equalizing Equations (5.b) and (5.c), which together with Equation (6) provide a system of three equations for finding fRe_a, a and b. Once these quantities have been determined, parameter c is readily calculated from Equation (7).

Table 1 shows some results for the parameters described above within the range of the parameters analyzed (power-law index, yield number and aspect ratio). It can be noted that for different aspect ratios as the power-law index increases, the product, fRe_a, increases as well. This fact can be explained by the rheological behavior of the fluid. For n < 1 the fluid presents a shear-thinning behavior, i.e., there is a decrease in the apparent fluid viscosity as the velocity gradient increases; for n > 1 it behaves as a shear-thickening fluid and the apparent fluid viscosity increases with increasing values for the velocity gradient, so in these cases higher pressure drops are needed to make the fluid flow.

In this table, it can also be observed that for higher yield numbers the region of plug flow increases, as was expected, since they correspond to higher yield stresses. The case of n = 1 and Y = 0 corresponds to the Newtonian situation, and in this analysis there is no plug-flow region, as is verified by finding the same values for parameters a, b and c. In addition, at higher values for the aspect ratios, for example g = 0.9, fRe_a is approximately equal to that in the case of a parallel-plate channel, i.e., fRe_a» 24.

Table 2 shows a comparison between the values for parameters a, b and c calculated in the present analysis and those presented by Hanks and Larsen (1979). It can be verified that the values are in excellent agreement, providing a direct validation of the numerical code developed in the present analysis for the computation of these parameters.

Once fRe_a and parameters a, b and c were obtained, the velocity profiles could be evaluated using Equations (5), which correspond to the three distinct regions of the flow field. In Figures 2 to 5 the velocity profiles are plotted for different values for power-law indices, yield numbers and aspect ratios. The figures represent a longitudinal cutaway of the concentric annular duct. According to Figure 2, when the value of the power-law index increases, maintaining the other parameters fixed, there is an increase in the value for velocity at the centerline of the annular duct. This can be explained by a decrease in the velocity gradients far from the duct walls, and consequently, a decrease in the apparent fluid viscosity for dilatant fluids.

In Figure 3, it can be observed that for the case where Y = 0 there is no plug-flow region, as was mentioned when Table 1 was analyzed. As the yield number is increased, a correspondent enlargement of the plug-flow region occurs as well. This results in an increase in the velocity gradients close to the duct walls, mainly at the inner wall. The effect of duct geometry on the velocity profiles is evidenced in Figure 4, where it can be observed that for higher values for the aspect ratios, the velocity profile tends to become symmetric; consequently, this is expected to affect the heat transfer rates at the duct walls, resulting in symmetric distributions for the local Nusselt numbers at both duct walls as well.

In Figure 5, it can be observed that for values of the power-law index greater than unity, the increase in the values for the yield number does not considerably influence the enlargement of the plug-flow region compared with the dimension of the plug-flow region for values of n £ 1 (see Figure 3). This implies in a uniform distribution of the local Nusselt numbers for the thermal entry region.

After the parameters related to the velocity field had been determined, quantities related to the eigenvalue problem given by Equations (16) were obtained via the sign-count method (Mikhailov and Özisik, 1984), i.e., eigenvalues, eigenfunctions, and normalization integrals, which are needed to compute important thermal parameters such as average temperature q_av(Z), local Nusselt numbers and temperatures at the inner and outer duct walls Nu_iw(Z), Nu_ow(Z), q(g,.Z) º q_iw(Z) and q(1,.Z) º q_ow(Z), respectively. Here, these parameters were calculated for the case of boundary conditions of the first kind, split into two fundamental situations, represented by cases A and B.

First, the results for the local Nusselt numbers obtained in the present work were validated by comparing them with those given in the work by Blackwell (1985), who analyzed thermally developing laminar flow of Bingham plastics within circular tubes. For this particular case, the velocity profiles given by Equations (5) are represented by

where C = t_o/t_w represents the ratio between the yield stress and the wall shear stress. Equations (30) were obtained from Equations (5) by making the aspect ratio g ® 0, which represents the situation of a circular tube.

Table 3 provides this comparison, showing the excellent agreement of the results. The slight difference between the present results and those of Blackwell (1985) for axial positions near the duct inlet, can be explained by the fact that this author computed few quantities related to the eigenvalue problem (only 60 eigenquantities) which were not adequate to compute fully converged results for the Nusselt numbers at those axial positions.

Table 4 shows the results for parameters Nu_iw(Z), Nu_ow(Z) and q_av(Z), related to cases A and B computed in the present work, and comparisons with those obtained by Shah and London (1978) and Nascimento et al. (2000) for the thermal entry region in laminar flow of Newtonian fluids within concentric annular ducts. In this analysis, values of 0.1 and 0.5 were taken for the aspect ratio. The excellent agreement between the three sets of results shown in this table also validates the numerical code developed in the present work as well as indicates the consistency of the present results.

Figures 6 and 7 show another comparison of the present results for the local Nusselt numbers with those obtained by Quaresma and Macêdo (1998), when the concentric annular duct tends towards two-limit configurations, that is, circular tubes and parallel-plate channels. As can be observed in these figures, the results are also in excellent agreement, once again validating the computational code developed here.

In Figures 8 to 13 axial distributions for the local Nusselt numbers for the thermal entry region are presented for different power-law indices, yield numbers and aspect ratios. An increase in the local Nusselt number for increasing yield numbers is verified; this can be clarified by higher gradients in the velocity field near the wall when the yield number increases. This effect tends to disappear for the fully developed region where the curves are practically coincident. It can also be noticed that as the power-law index increases, the difference in the values of the Nusselt numbers, for different yield numbers, at the beginning of the thermal entry region, becomes less substantial, when compared with those obtained for the cases of power-law indices less than unity. This is due to the fact that the region of plug flow remains practically unaffected for dilatant fluids with increasing values for the yield number, as previously discussed (see Figure 5).

Finally, from these figures, it can also be observed that the local Nusselt numbers for the two cases studied here, Nu_iw(Z) and Nu_ow(Z), tend towards the same distribution when the aspect ratio increases. For example, for the case of aspect ratio g = 0.9, which practically represents a parallel-plate channel, this symmetry is clearly observed.

CONCLUSIONS

The proposed integral transform approach in conjunction with the sign-count method provided reliable and cost-effective simulations of heat and fluid flow of Herschel-Bulkley fluids within concentric annular ducts in laminar regime. Benchmark results for the product of the Fanning friction factor-apparent Reynolds number, velocity field and Nusselt numbers were systematically tabulated and graphically presented for different values for yield numbers, aspect ratios, and power-law indices, illustrating the effectiveness of the present methodology. Comparisons with previous work in the literature were also performed, demonstrating excellent agreement and furnishing direct validations of the present results as well as showing that they were consistent.

NOMENCLATURE

a^*, a lower bounds on the plug-flow region, dimensional and dimensionless, respectively b^*, b upper bounds on the plug-flow region, dimensional and dimensionless, respectively c^*, c radial positions where the shear stress is zero, dimensional and dimensionless, respectively C ratio between the yield stress and the wall shear stress c_p specific heat D_h hydraulic diameter f

Fanning friction factor

h_iw(z) local heat transfer coefficient at the inner wall h_ow(z) local heat transfer coefficient at the outer wall k thermal conductivity K consistency index of the fluid m coefficient defined in Equations (10.c,d) n power-law index N_i normalization integral Nu(Z) local Nusselt number for the case of a circular tube Nu_av average Nusselt number Nu_iw(Z) local Nusselt number at the inner wall Nu_ow(Z) local Nusselt number at the outer wall Pr_a apparent Prandtl number r, R radial coordinates, dimensional and dimensionless, respectively r_iw, r_ow inner and outer radii, respectively Re_a apparent Reynolds number s coefficient defined in Equation (4.j) T_e inlet temperature T_iw, T_ow prescribed temperatures at the inner and outer walls, respectively u_av average flow velocity u(r), U(R) velocity distributions, dimensional and dimensionless, respectively W(R) defined by Equation (11.a) Y yield number z, Z axial coordinates, dimensional and dimensionless, respectively

Greek Letters

a fluid thermal diffusivity shear rate g aspect ratio q (R,Z) dimensionless temperature distribution q_av(Z) dimensionless average temperature q_iw(Z) dimensionless temperature at the inner wall q_ow(Z) dimensionless temperature at the outer wall m_i eigenvalues of problem (16) r density t_rz shear stress t₀ yield stress y_i(R) eigenfunctions of problem (16)

• Blackwell, B.F., Numerical Solution of the Graetz Problem for a Bingham Plastic in Laminar Tube Flow with Constant Wall Temperature, J. Heat Transfer, 107, 466-468 (1985).
• Cotta, R.M., Integral Transforms in Computational Heat and Fluid Flow, CRC Press, Boca Raton (1993).
• Fredrickson, A.G. and Bird, R.B., Non-Newtonian Flow in Annuli, Ind. Eng. Chem., 50, 347-383 (1958).
• Hanks, R.W. and Larsen, K.M., The Flow of Power-Law Fluids in Concentric Annuli, Ind. Eng. Chem. Fundam., 18, 33-35 (1979).
• Laird, W.M., Slurry and Suspension Transport, Ind. Eng. Chem., 49, 138-141 (1957).
• Mendes, P.R.S. and Naccache, M.F., Heat Transfer to Herschel-Bulkley Fluids in Laminar Fully Developed Flow through Tubes, Proc. of the 13^th Brazilian Congress of Mechanical Engineering (XIII COBEM), Belo Horizonte, Brazil (1995) (on CD-ROM).
• Mikhailov, M.D., Splitting Up of Heat-Conduction Problems, Letters Heat Mass Transfer, 4, 163-166 (1977).
• Mikhailov, M.D. and Özisik, M.N., Unified Analysys and Solutions of Heat and Mass Diffusion, John Wiley, New York (1984).
• Mikhailov, M.D. and Vulchanov, N.L., Computational Procedure for Sturm-Liouville Problems, J. Comp. Phys., 50, 323-336 (1983).
• Nascimento, U.C.S., Macędo, E.N. and Quaresma, J.N.N., Solution for the Thermal Entry Region in Laminar Flow of Bingham Plastics within Annular Ducts via Integral Transformation, Hybrid Methods in Engineering, 2, 233-247 (2000).
• Nouar, C., Devienne, R. and Lebouche, M., Convection Thermique pour un Fluide de Herschel-Bulkley dans la Région d'entrée d'une Conduit, Int. J. Heat Mass Transfer, 37, 1-12 (1994).
• Nouar, C., Devienne, R., Lebouche, M. and Riou, C., Numerical Analysis of the Thermal Convection for Herschel-Bukley Fluids, Int. J. Heat Fluid Flow, 16, 223-231 (1995).
• Quaresma, J.N.N. and Macędo, E.N., Integral Transform Solution for Forced Convection of Herschel-Bulkley Fluids in Circular Tubes and Parallel-Plates Ducts, Brazilian Journal of Chemical Engineering, 15, 77-89 (1998).
• Shah, R.K. and London, A.L., Laminar Flow Forced Convection in Ducts, Supplement 1 to Advances in Heat Transfer, eds. T. F. Irvine, Jr. and J. P. Hartnett, Academic Press, New York (1978).
• Soares, M., Mendes, P.R.S. and Naccache, M.F., Heat Transfer to Herschel-Bulkley Fluids in Laminar Flow through Short Tubes, Proc. of the VI ENCIT/VI LATCYM, Florianópolis, Brazil, 1575-1580 (1996).
• Soares, M., Mendes, P.R.S. and Naccache, M.F., Heat Transfer to Viscoplastic Fluids in Laminar Flow through Isothermal Short Tubes, J. of the Braz. Soc. Mechanical Sciences, 19, 1-14 (1997).

Publication Dates

• Publication in this collection
  17 Apr 2002
• Date of issue
  Dec 2001

History

• Received
  16 May 2001
• Accepted
  18 Aug 2001