Applied Mathematics and Computation 189 (2007) 690–697 www.elsevier.com/locate/amc
Thermal stability of a reactive third grade fluid in a cylindrical pipe: An exploitation of Hermite–Pade´ approximation technique Oluwole Daniel Makinde Applied Mathematics Department, University of Limpopo, Private Bag X1106, Sovenga 0727, South Africa
Abstract This study is devoted to investigate the thermal stability of a reactive third grade fluid flowing steadily through a cylindrical pipe with isothermal wall. It is assumed that the reaction is exothermic under Arrhenius kinetics, neglecting the consumption of the material. Approximate solutions are constructed for the governing nonlinear boundary value problem using regular perturbation techniques together with a special type of Hermite–Pade´ approximants and important properties of the velocity and temperature fields including bifurcations and thermal criticality conditions are discussed. Ó 2006 Elsevier Inc. All rights reserved. Keywords: Cylindrical pipe; Third grade fluid; Arrhenius kinetics; Thermal stability; Hermite–Pade´ approximants
1. Introduction In recent years, considerable attention has been paid to flow problems of non-Newtonian fluids due to the complexity of fluids in engineering and industrial applications [2]. A large class of real fluids used in industries is chemically reactive and exhibit non-Newtonian characteristics, e.g. coal slurries, polymer solutions or melts, drilling mud, hydrocarbon oils, grease, etc. Because of the nonlinear relationship between stress and the rate of strain, the analysis of the behavior of such fluids tends to be more complicated and subtle in comparison with that of Newtonian fluids. Some of the earlier work on third grade fluids include: Ayub et al. [1], Fosdick and Rajagopal [6], Massoudi and Christie [14], Rajagopal [15], Yurusoy and Pakdemirli [18]. Meanwhile, the study of heat transfer and thermal stability of reactive non-Newtonian fluids is extremely important in order to ensure safety of life and properties during handling and processing of such fluids [5]. Mathematically speaking, the thermal boundary layer equation for a reactive third grade non-Newtonian fluids constitute a nonlinear problem and the long-term behavior of the solutions in space will provide us an insight into inherently complex physical process of thermal instability in the system. The theory of nonlinear differential equations is quite elaborate and their solution remains an extremely important problem of
E-mail address:
[email protected] 0096-3003/$ - see front matter Ó 2006 Elsevier Inc. All rights reserved. doi:10.1016/j.amc.2006.11.124
O.D. Makinde / Applied Mathematics and Computation 189 (2007) 690–697
691
practical relevance in sciences and engineering. Several numerical approaches have developed in the last few decades, e.g. finite differences, spectral method, shooting method, etc. to tackle this type of problem [4]. More recently, the ideas of classical analytical methods have experienced a revival, in connection with the proposition of novel hybrid numerical–analytical schemes for nonlinear differential equations. One such trend is related to Hermite–Pade´ approximation approach, ([8,10–13,16,17]). This approach, over the last few years, proved itself as a powerful tool and a potential alternative to traditional numerical techniques in various applications in sciences and engineering. This semi-numerical approach is also extremely useful in the validation of purely numerical scheme. Hence, the purpose of the present work is to investigate the thermal stability of a reactive third grade liquid flowing steadily through a cylindrical pipe with isothermal wall using a special type of Hermite–Pade´ approximants. The mathematical formulation of the problem is established and solved in Sections 2 and 3. In Section 4 we introduce and apply some rudiments of Hermite–Pade´ approximation technique. Both numerical and graphical results are presented and discussed quantitatively with respect to various parameters embedded in the system in Section 5. 2. Mathematical model Consider the steady flow of an incompressible reactive third grade liquid through a circular cylindrical pipe with isothermal wall. It is assumed that the fluid motion is induced by applied axial pressure gradient. We choose a cylindrical polar co-ordinate system ðr; h; zÞ where 0z lies along the axis of pipe, r is the distance measured radially such that r ¼ a is the pipe radius and h is the azimuthal angle (see Fig. 1). For hydrodynamically and thermally developed flow, both velocity and temperature fields depend on r only. Following ([1,14,18]) and neglecting the reacting viscous fluid consumption, the one dimensional governing equations for the momentum and heat balance can be written as: 3 ! l d du b3 d du dP r r ; ð1Þ þ ¼ r dr dr dr dz r dr 2 2 ! k d dT du du r l þ b3 ð2Þ þ þ QC 0 A eE=RT ¼ 0; r dr dr dr dr subject to the axial-symmetric and boundary conditions: du dT ð0Þ ¼ ð0Þ ¼ 0; dr dr
uðaÞ ¼ 0;
T ðaÞ ¼ T 0 ;
ð3Þ
where the additional Arrhenius kinetics term in energy balance equation is due to [3,7]. Here T is the absolute temperature, u the fluid velocity, U the fluid characteristic velocity, T0 the plate temperature, k the thermal conductivity of the material, Q the heat of reaction, A the rate constant, E the activation energy, R the universal gas constant, C0 the initial concentration of the reactant species, b3 the material coefficient, P the modified pressure and l is the fluid dynamic viscosity coefficient, [2,3,13]. We introduce the following dimensionless variables into Eqs. (1)–(3);
u = 0,
r
T=T0 ,
r=a
Combustible fluid
Fig. 1. Geometry of the problem.
z
692
O.D. Makinde / Applied Mathematics and Computation 189 (2007) 690–697
h¼
EðT T 0 Þ ; RT 20
r r ¼ ; a
lG2 U 2 eE=RT 0 m¼ ; QAa2 C 0
k¼
RT 0 ; e¼ E
QEAa2 C 0 eE=RT 0 ; T 20 Rk a2 dP ; G¼ lU dz
W ¼
u ; UG
b U 2 G2 c¼ 3 2 ; al
ð4Þ
and obtain the dimensionless governing equations together with their corresponding boundary conditions as (neglecting the bar symbol for clarity): 3 ! 1 d dW c d dW r r þ ¼ 1; ð5Þ r dr dr r dr dr " 2 2 !# 1 d dh dW dW ðh=ð1þehÞÞ r þm 1þc þk e ¼ 0; ð6Þ r dr dr dr dr with dW dh ð0Þ ¼ ð0Þ ¼ 0; dr dr
W ð1Þ ¼ 0;
hð1Þ ¼ 0;
ð7Þ
where k, e, c, m represent the Frank–Kamenetskii parameter, activation energy parameter, the dimensionless non-Newtonian parameter and the viscous heating parameter, respectively. In the following sections, Eqs. (5)– (7) are solved using both perturbation and multivariate series summation techniques ([11,16,17]). 3. Perturbation method Due to the nonlinear nature of the velocity and temperature fields Eqs. (5) and (6), it is convenient to form a power series expansion both in the dimensionless non-Newtonian parameter c and the Frank–Kamenetskii parameter k, i.e., 1 1 X X W ¼ W i ci ; h ¼ hi ki : ð8Þ i¼0
i¼0
Substituting the solution series (8) into Eqs. (5)–(7) and collecting the coefficients of like powers of c and k, we obtained and solved the equations for the coefficients of solution series iteratively. The solution for the velocity and temperature fields are given as 1 1 2 1 4 1 1 6 1 2 3 8 3 11 10 11 3 r r r þ W ðrÞ ¼ r þ cþ r þ c þ c þ c4 4 4 32 32 64 64 256 256 1024 1024 91 12 91 r þ ð9Þ c5 þ Oðc6 Þ 8192 8192 1 1 3 3 1 1 1 mc5 r14 þ mc4 r12 mc3 r10 þ mc2 r8 þ mcr8 mr4 r2 hðrÞ ¼ 802816 36864 12800 4096 576 64 4 1 1 3 3 1 1 1 1 mc5 þ mc4 mc3 mc2 mc þ m þ k þ mc5 r16 þ 802816 36864 12800 4096 576 64 4 205520896 1 1 3 1 1 1 1 mc4 r14 þ mc3 r12 mc2 r10 mcr8 þ mr6 þ r4 r2 7225344 614400 409600 36864 2304 64 16 3 3 1 1 1 2 1 2 r2 mc3 þ r2 mc2 r2 mc5 þ r2 mc4 þ r mc r m 51200 16384 3211264 147456 2304 256 9 1 7 9 5 1 3 2 mc5 mc4 þ mc3 c2 mc þ mþ þ k þ Oðk3 Þ ð10Þ 29360128 150528 122880 51200 12288 288 64 Using a computer symbolic algebra package (MAPLE), we obtained the first twenty terms of the above solution series (9) and (10). The physical quantities of interest in this problem are the skin-friction parameter ðC f Þ and the Nusselt number (Nu) which are defined by
O.D. Makinde / Applied Mathematics and Computation 189 (2007) 690–697
aiw dW ð1Þ; ¼ dr lUG 2aEqw dh Nu ¼ ¼ 2 ð1Þ; 2 dr kRT 0
693
ð11aÞ
Cf ¼
ð11bÞ
where iw ¼ l du=dr and qw ¼ k dT =dr are the shear stress and the heat flux evaluated at the wall (i.e. r ¼ a), respectively. We are aware that the power series solutions (9), (10) are valid for very small parameter values. However, using Hermite–Pade´ approximation technique, we have extended the usability of the solution series beyond small parameter values as illustrated in the following section. 4. Thermal criticality and bifurcation study The thermal instability properties of the materials under consideration can easily be determined by the appearance of thermal runaway phenomenon in the system or through non-existence of steady-state solution for certain parameter values. In order to determine the appearance of thermal instability in the system together with the evolution of temperature field (i.e. k > 0), we employ a special type of Hermite–Pade´ approximation technique. Suppose that the partial sum U N 1 ðkÞ ¼
N 1 X
ai ki ¼ U ðkÞ þ OðkN Þ
as k ! 0;
ð12Þ
i¼0
is given. We are concerned with the bifurcation study by analytic continuation as well as the dominant behavior of the solution by using partial sum (12). We expect that the accuracy of the critical parameters will ensure the accuracy of the solution. It is well known that the dominant behaviour of a solution of a differential equation can often be written as Guttamann [8], ( b H ðkc kÞ for b 6¼ 0; 1; 2; . . . U ðkÞ ð13Þ ask ! kc ; b H ðkc kÞ ln jkc kj for b ¼ 0; 1; 2; . . . where H is some constant and kc is the critical point with the exponent a. However, we shall make the simplest hypothesis in the context of nonlinear problems by assuming that U ðkÞ is the local representation of an algebraic function of k. Therefore, we seek an expression of the form F d ðk; U N 1 Þ ¼ A0N ðkÞ þ Ad1N ðkÞU ð1Þ þ Ad2N ðkÞU ð2Þ þ Ad3N ðkÞU ð3Þ ;
ð14Þ
such that AdiN ðkÞ ¼
A0N ðkÞ ¼ 1;
dþi X
bij kj1 ;
ð15Þ
j¼1
and F d ðk; U Þ ¼ OðkN þ1 Þ as k ! 0;
ð16Þ
where d P 1; i ¼ 1; 2; 3. The condition (15) normalizes the Fd and ensures that the order of series AdiN increases as i and d increase in value. There are 3ð2 þ dÞ undetermined coefficients bij in the expression (15). The requirement (16) reduces the problem to a system of N linear equations for the unknown coefficients of Fd. The entries of the underlying matrix depend only on the N given coefficients ai. Henceforth, we shall take N ¼ 3ð2 þ dÞ;
ð17Þ
so that the number of equations equals the number of unknowns. Eq. (16) is a new special type of Hermite– Pade´ approximants. Both the algebraic and differential approximants forms of Eq. (16) are considered. For instance, we let U ð1Þ ¼ U ;
U ð2Þ ¼ U 2 ;
U ð3Þ ¼ U 3 ;
ð18Þ
694
O.D. Makinde / Applied Mathematics and Computation 189 (2007) 690–697
and obtain a cubic Pade´ approximant. This enables us to obtain solution branches of the underlying problem in addition to the one represented by the original series. In the same manner, we let U ð1Þ ¼ U ;
U ð2Þ ¼ DU ;
U ð3Þ ¼ D2 U ;
ð19Þ
in Eq. (15), where D is the differential operator given by D ¼ d=dk. This leads to a second order differential approximants. It is an extension of the integral approximants idea by Hunter and Baker [9] and enables us to obtain the dominant singularity in the flow field i.e. by equating the coefficient Ad3N ðkÞ in the Eq. (16) to zero. Meanwhile, it is very important to know that the rationale for chosen the degrees of AdiN in Eq. (15) in this particular application is based on the simple technique of singularity determination in second order linear ordinary differential equation with polynomial coefficients as well as the possibility of multiple solution branches for the nonlinear problem, [4]. In practice, one usually finds that the dominant singularities are located at zeroes of the leading polynomial Ad3N coefficients of the second order linear ordinary differential equation. Hence, some of the zeroes of Ad3N may provide approximations of the singularities of the series U and we expect that the accuracy of the singularities will ensure the accuracy of the approximants. The critical exponent bN can easily be found by using Newton’s polygon algorithm. However, it is well known that, in the case of algebraic equations, the only singularities that are structurally stable are simple turning points. Hence, in practice, one almost invariably obtains bN ¼ 1=2. If we assume a singularity of algebraic type as in Eq. (13), then the exponent may be approximated by bN ¼ 1
A2N ðkCN Þ : DA3N ðkCN Þ
ð20Þ
5. Results and discussion The bifurcation procedure in section (4) above was applied to the first 20 terms of the solution series and we obtained the results as shown in Tables 1 and 2 below: Table 1, demonstrates the rapid convergence of our procedure with respect to the non-Newtonian parameter ðcc Þ dominant singularity and its critical exponent ðbc Þ for a third grade fluid flowing steadily in a cylindrical pipe. A bifurcation point (i.e. a turning point) occurs in the flow field at ðcc ; C f Þ ¼ ð16=27; 3=4Þ as shown in Fig. 4. It is noteworthy that the convergence of our procedure improves with gradual increase in the number of series coefficients utilized in the approximants. Table 2, shows the thermal criticality conditions ðkc Þ for a reactive third grade liquid with respect to pipe flow. Generally, the magnitude of thermal criticality increases while the rate of heat transfer across the wall decreases with an increase in value non-Newtonian parameter. This implies that increasing values of the non-Newtonian parameter enhances thermal stability of a reactive third grade liquid. The fluid velocity profile with respect to increasing value of c is shown in Table 1 Computations showing the procedure rapid convergence and bifurcation point in the velocity field d
N
Cf
cc
bcN
1 2 3 4
9 12 15 18
0.750095 0.750000 0.750000 0.750000
0.592617 0.592592 0.592592 0.592592
0.499987 0.500000 0.500000 0.500000
Table 2 Computations showing thermal ignition criticality for different parameter values m
c
Nu ðe ¼ 0Þ
kc ðe ¼ 0Þ
Nu ðe ¼ 0:1Þ
kc ðe ¼ 0:1Þ
bcN
1.0 1.0 1.0 1.0
0.0 0.1 0.2 0.3
4.221370958 4.217362228 4.212964486 4.208225870
1.94543584443 1.94607136100 1.94675509035 1.94748243186
5.374510971 5.369541143 5.364098368 5.358240483
2.20683839 2.20747552 2.20816105 2.20889033
0.5 0.5 0.5 0.5
O.D. Makinde / Applied Mathematics and Computation 189 (2007) 690–697
695
Fig. 2. A Poiseuille parabolic velocity profile is observed, however, a gradual decrease in the magnitude of fluid velocity is noticed with an increase in value of non-Newtonian parameter. Fig. 3 shows a transverse increase in the fluid temperature with maximum temperature along the centreline. Furthermore, it is noteworthy that the fluid temperature generally increases with an increase in the value of Frank–Kamenetskii
Fig. 2. Velocity profile for (___) c ¼ 0; () c ¼ 0:4; (+++) c ¼ 0:8.
Fig. 3. Temperature profile for m ¼ 1; e ¼ 0; c ¼ 0:1; (___) k ¼ 0:1; () k ¼ 0:2; (+++) k ¼ 0:3.
Fig. 4. A bifurcation diagram in the ðc; C f Þ plane for a third grade liquid.
696
O.D. Makinde / Applied Mathematics and Computation 189 (2007) 690–697
Nu 8.4 II 4.2 λ c =1.946071361 I 0
1.0
2.0
3.0
λ
Fig. 5. A slice of approximate bifurcation diagram in the ðk; Nu ðm ¼ 1; c ¼ 0:1; e ¼ 0ÞÞ plane.
parameter ðkÞ due to Arrhenius kinetics. A slice of the bifurcation diagram for c > 0 in the ðk; NuÞ plane is shown in Fig. 5. It represents the variation of wall heat flux (Nu) with the Frank-Kamenetskii parameter ðkÞ. In particular, for every 0 6 e 6 0:1 there is a critical value kc (a turning point) such that, for 0 6 k < kc there are two solution branches (labelled I and II). The upper and lower solution branches occur due to Arrhenius kinetics in the governing thermal boundary layer equation (Eq. (3)). When kc < k the system has no real solution and displays a classical form indicating thermal runaway. The magnitude of kc increases with a decrease in the fluid activation energy ðe ¼ 0:1Þ, hence, preventing the early development of thermal runaway and enhancing thermal stability. 6. Conclusion In this paper, a novel hybrid numerical–analytical scheme based on a special type of Hermite–Pade´ approximants is employed to investigate the thermal stability of a reactive third grade fluid flowing steadily through a cylindrical pipe with isothermal wall. The procedure reveals accurately the thermal stability conditions and various solution branches. It is observed that an increase in non-Newtonian parameter enhances the thermal stability of a reactive third grade fluid. Acknowledgements The author gratefully appreciates the financial support received from the National Research Foundation (NRF) of South Africa Thuthuka programme. References [1] M. Ayub, A. Rasheed, T. Hayat, Exact flow of a third grade fluid past a porous plate using homotopy analysis method, Int. J. Eng. Sci. 41 (2003) 2091. [2] G.K. Batchelor, An Introduction to Fluid Dynamics, Cambridge University Press, 1967. [3] J. Bebernes, D. Eberly, Mathematical Problems from Combustion Theory, Springer-Verlag, New York, 1989. [4] C. Bender, S.A. Orszag, Advanced Mathematical Methods for Scientists and Engineers, McGraw-Hill, 1978. [5] P.C. Bowes, Self-heating: Evaluating and Controlling the Hazard, Elsevier, Amsterdam, 1984. [6] R.L. Fosdick, K.R. Rajagopal, Thermodynamics and stability of fluids of third grade, Proc. Roy. Soc. London A 339 (1980) 351. [7] D.A. Frank Kamenetskii, Diffusion and Heat Transfer in Chemical Kinetics, Plenum Press, New York, 1969. [8] A.J. Guttamann, Asymptotic analysis of power – series expansions, in: C. Domb, J.K. Lebowitz (Eds.), Phase Transitions and Critical Phenomena, Academic Press, New York, 1989, pp. 1–234. [9] D.L. Hunter, G.A. Baker, Methods of series analysis III: integral approximant methods, Phys. Rev. B 19 (1979) 3808–3821. [10] O.D. Makinde, Exothermic explosions in a slab: a case study of series summation technique, Int. Commun. Heat Mass Transfer 31 (8) (2004) 1227–1231. [11] O.D. Makinde, Strong exothermic explosions in a cylindrical pipe: a case study of series summation technique, Mech. Res. Commun. 32 (2005) 191–195.
O.D. Makinde / Applied Mathematics and Computation 189 (2007) 690–697
697
[12] O.D. Makinde, Laminar falling liquid film with variable viscosity along an inclined heated plate, Appl. Math. Comput. 175 (2006) 80– 88. [13] O.D. Makinde, Thermal ignition in a reactive viscous flow through a channel filled with a porous medium, J. Heat Transfer 128 (2006) 601–604. [14] M. Massoudi, I. Christe, Effects of variable viscosity and viscous dissipation on the flow of a third grade fluid in a pipe, Int. J. NonLinear Mech. 30 (1995) 687. [15] K.R. Rajagopal, On boundary conditions for fluids of the differential type, in: A. Sequira (Ed.), Navier–Stokes Equations and Related Non-linear Problems, Plenum Press, New York, 1995, p. 273. [16] A.V. Sergeyev, D.Z. Goodson, Summation of asymptotic expansions of multiple valued functions using algebraic approximationsapplication to anharmonic oscillators, J. Phys. A31 (1998) 4301–4317. [17] Y. Tourigny, P.G. Drazin, The asymptotic behaviour of algebraic approximants, Proc. Roy. Soc. London A456 (2000) 1117–1137. [18] M. Yurusoy, M. Pakdemirli, Approximate analytical solutions for the flow of a third grade fluid in a pipe, Int. J. Non-Linear Mech. 37 (2002) 187.