Appl. Math. Mech. -Engl. Ed., 33(8), 1015–1034 (2012) DOI 10.1007/s10483-012-1602-8 c Shanghai University and Springer-Verlag Berlin Heidelberg 2012

Applied Mathematics and Mechanics (English Edition)

Mixed convection flow in vertical channel with boundary conditions of third kind in presence of heat source/sink∗ J. C. UMAVATHI,

J. PRATHAP KUMAR,

JAWERIYA SULTANA

(Department of Mathematics, Gulbarga University, Gulbarga 585106, Karnataka, India)

Abstract The effects of viscous dissipation and heat source/sink on fully developed mixed convection for the laminar flow in a parallel-plate vertical channel are investigated. The plate exchanges heat with an external fluid. Both conditions of equal and different reference temperatures of the external fluid are considered. First, the simple cases of the negligible Brinkman number or the negligible Grashof number are solved analytically. Then, the combined effects of buoyancy forces and viscous dissipation in the presence of heat source/sink are analyzed by a perturbation series method valid for small values of the perturbation parameter. To relax the conditions on the perturbation parameter, the velocity and temperature fields are solved by using the Runge-Kutta fourth-order method with the shooting technique. The velocity, temperature, skin friction, and Nusselt numbers at the plates are discussed numerically and presented through graphs. Key words mixed convection, viscous fluid, perturbation method, Runge-Kutta shooting method, heat source/sink Chinese Library Classification O357.1 2010 Mathematics Subject Classification

76D50

Nomenclature Bi1 , Bi2 , Br, cp , g, Gr, h1 , h2 , k, L, Nu1 , Nu2 , p, P, P r, Q, Re,

Biot numbers; Brinkman number; specific heat at constant pressure; acceleration due to gravity; Grashof number; external heat transfer coefficients; thermal conductivity; channel width; Nusselt numbers; pressure; difference between the pressure and the hydrostatic pressure; Prandtl number; rate of internal heat absorption/ generation; Reynolds number;

RT , T, T1 , T 2 , T0 , u, u ¯, U, U0 , x, X, y, Y,

temperature difference ratio; temperature; reference temperatures of the external fluid; reference temperature; dimensionless velocity in the X-direction; mean value of u; velocity component in the X-direction; reference velocity; dimensionless streamwise coordinate; streamwise coordinate; dimensionless transverse coordinate; transverse coordinate.

∗ Received Mar. 16, 2011 / Revised Mar. 29, 2012 Corresponding author J. C. UMAVATHI, Professor, Ph. D., E-mail: jc [email protected]

1016

J. C. UMAVATHI, J. PRATHAP KUMAR, and JAWERIYA SULTANA

Greek letters α, β, ε, φ, λ,

1

thermal diffusivity; thermal expansion coefficient; dimensionless parameter as defined by (54); dimensionless parameter of heat absorption/generation coefficient; dimensionless parameter as defined by (19);

μ, ν, θ, θb , ρ, ρ0 ,

viscosity; kinematic viscosity; dimensionless temperature; dimensionless bulk temperature; mass density; value of the mass density when T = T0 .

Introduction

Buoyancy induced flows in ducts deserve wide attention mainly for their engineering applications in several thermal control devices ranging from electronics to nuclear plants. Indeed, in passive or semi-passive thermal control systems, either purely free convection flows or mixed convection flows are involved. Several studies on mixed convection problems for a Newtonian fluid in a vertical channel have already been presented in the literature. In particular, some analytical solutions for the fully developed flow have been performed. The boundary conditions of uniform wall temperatures have been analyzed by Aung and Worku[1] . The boundary conditions of uniform wall temperatures on a wall and uniform wall heat fluxes have been studied by Hamadah and Wirtz[2] and by Cheng et al.[3] . The effects of viscous dissipation on the velocity and on the temperature profiles have been analyzed by Barletta[4] for the boundary conditions of uniform wall temperatures and by Zanchini[5] for boundary conditions of the third kind. The effect of viscous dissipation may become very important in several flow configurations occurring in the engineering practice. In fact, viscous dissipation affects strongly the heat transfer process whenever the operating fluid has a low thermal conductivity, a high viscosity, and flows in ducts with a small cross-section and a small wall heat flux. All these features may occur, for instance, in the micro-channel flows. As is well known, the effect of viscous heating increases with the square of the mass flow rate and, as a consequence, becomes specially important under conditions of forced convection. In the past, the laminar forced convection heat transfer in the hydrodynamic entrance region of a flat rectangular channel has been investigated either for the temperature boundary condition of the first kind, characterized by prescribed wall temperature[6–8] , or for the boundary condition of the second kind, expressed by the prescribed wall temperature heat flux[8–9] . A more realistic condition in many applications, however, is the temperature boundary condition of the third kind: the local wall heat flux is a linear function of the local wall temperature. The natural convection process in the presence of heat source/sink is presented in various physical phenomena such as fire engineering, combustion modeling, nuclear energy, heat exchangers, petroleum reservoir, etc. The liquid metal having low Prandtl number (because of very large thermal conductivity) are generally used as coolants and have applications in manufacturing processes such as the cooling of the metallic plate, nuclear reactor, etc. The liquid metal has ability to transport heat even if small temperature difference exists between the surface and fluid. For this reason, the liquid metal is used as a coolant in the nuclear reactor to transfer waste heat from the core region. Heat transfer in the laminar region of a flat channel for the temperature boundary conditions of the third kind was explored by Javeri[10]. Javeri[11] investigated the influences of the temperature boundary condition of the third kind on the laminar heat transfer in the thermal entrance region of a rectangular channel. Later in 1998, Zanchini[5] analyzed the effect of viscous dissipation on mixed convection in a vertical channel with boundary conditions of third kind. Kumari and Nath[12] analyzed the effect of localized cooling/heating and injection/suction on mixed convection flow on a thin vertical cylinder. Mixed convection flow in a

Mixed convection flow in vertical channel with boundary conditions

1017

vertical channel embedded in a porous medium with asymmetric wall heating conditions in the presence of source/sink was analyzed by Umavathi et al.[13] . Mahanthi and Gaur[14] investigated the effect of the viscosity and thermal conduction on the steady free-convective flow of a viscous incompressible fluid along an isothermal vertical plate in the presence of heat sink. Umavathi and Prathap[15] found the exact solutions for the mixed convection flow of a micropolar fluid in a vertical channel with symmetric and asymmetric conditions in the presence of source/sink. The aim of this paper is to extend the analysis performed by Zanchini[5] in the presence of heat source/sink and to solve it numerically. Both equal and different reference temperature of the external fluid, as well as both different and equal Biot numbers are considered. In the absence of source/sink, the solutions obtained in this paper coincide with Zanchini[5] for small values of perturbation parameter. However, the solutions for large values of perturbation parameters are obtained using the Runge-Kutta fourth-order method with the shooting technique.

2

Mathematical formulation

Consider the steady and laminar flow of a Newtonian fluid in the fully developed region of a parallel-plate vertical channel as shown in Fig. 1. The X-axis lies on the axial plane of the channel, and its direction is opposite to the gravitational field. The Y -axis is orthogonal to the walls. The channel occupies the region of space −L/2  Y  L/2. The thermal conductivity, the thermal diffusivity, the dynamic viscosity, and the thermal expansion coefficient of the fluid are assumed to be constants.

Fig. 1

Physical model and coordinate system

As customary, the Boussineq approximation and the equation of state are adopted as ρ = ρ0 (1 − β(T − T0 )).

(1)

Moreover, it is assumed that the only nonzero component of the velocity field U is the X-component of U . Thus, since ∇ · U = 0, one has ∂U =0 ∂X

(2)

such that U depends only on Y . The momentum balance equations along the X- and Y -axes yield[16] d2 U 1 ∂P +ν βg(T − T0 ) − = 0, (3) ρ0 ∂X dY 2

1018

J. C. UMAVATHI, J. PRATHAP KUMAR, and JAWERIYA SULTANA

∂P = 0, ∂Y

(4)

where P = p+ρ0 gX. Since on account of (4), P depends only on X. Hence, (3) can be rewritten as 1 dP ν d2 U T − T0 = . (5) − βgρ0 dX βg dY 2 From (5), one obtains ∂T 1 d2 P = , ∂X βgρ0 dX 2

(6)

∂T ν d3 U =− , ∂Y βg dY 3

(7)

∂2T ν d4 U = − . ∂Y 2 βg dY 4

(8)

Both the walls of the channel are assumed to have a negligible thickness and to exchange heat by convection with an external fluid. In particular, at Y = −L/2, the external convection coefficient is considered as uniform with the value h1 , and the fluid in the region Y < −L/2 is assumed to have a uniform reference temperature T . At Y = L/2, the external convection coefficient is considered as uniform with the value h2 , and the fluid in the region Y > L/2 is supposed to have a uniform reference temperature T2  T1 . Therefore, the boundary conditions on the temperature field can be expressed as −k

 ∂T  = h1 (T1 − T (X, −L/2)), ∂Y Y =−L/2

−k

 ∂T  = h2 (T (X, L/2) − T2 ). ∂Y Y =L/2

On account of (7), (9) and (10) can be rewritten as  d3 U  βgh1 (T1 − T (X, −L/2)), =  3 dY Y =−L/2 kν  d3 U  βgh2 (T (X, −L/2) − T2 ). = dY 3 Y =L/2 kν

(9)

(10)

(11)

(12)

∂T is zero both at Y = −L/2 and at It is easily verified that (11) and (12) imply that ∂X ∂T Y = L/2. Since (6) ensures that ∂X does not depend on Y , which leads to the conclusion that ∂T ∂X is zero everywhere. Therefore, the temperature T depends only on Y , i.e., T = T (Y ). Thus, on account of (6), there exists a constant A such that

dP = A. dX

(13)

Mixed convection flow in vertical channel with boundary conditions

1019

For the problem under exam, the energy balance equation in the presence of viscous dissipation can be written as[16] ∂ 2T ν d2 U  ν  dU 2 Q  1 dP − . (14) =− ∓ 2 ∂Y αcp dY k βgρ0 dX βg dY 2 Equations (8) and (14) yield a differential equation for U , i.e., d4 U gβ  dU 2 Q  d2 U A  . = ∓ − 4 2 dY αcp dY k dY νρ0

(15)

The boundary conditions on U are U (−L/2) = U (L/2) = 0.

(16)

Together with (11) and (12), and on account of (5), (11) and (12) can be rewritten as   d3 U  h1 d2 U  − dY 3 Y =−L/2 k dY 2 Y =−L/2 βgh1 Ah1 − (T0 − T1 ), kμ kν   d3 U  h2 d2 U  + dY 3 Y =L/2 k dY 2 Y =L/2

=−

=

βgh2 Ah2 − (T2 − T0 ). kμ kν

(17)

(18)

Equations (15)–(18) determine the velocity distribution. They can be written in a dimensionless form by means of the following dimensionless parameters: ⎧ U Y T − T0 ⎪ ⎪ u= , y= , , θ= ⎪ ⎪ U0 ΔT D ⎪ ⎪ ⎪ ⎪ ⎪ 3 ⎪ ⎪ gβΔT D μU02 U0 D ⎪ ⎪ Gr = , Br = , , Re = ⎪ ⎪ ν2 ν kΔT ⎪ ⎪ ⎪ ⎪ ⎨ Gr μ cp T2 − T1 , λ= , RT = , Pr = ⎪ k Re ΔT ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ Bi = h1 D , Bi = h2 D , ⎪ 1 2 ⎪ ⎪ k k ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ QD2 Bi1 Bi2 ⎪ ⎩S = . , φ= Bi1 Bi2 + 2Bi1 + 2Bi2 k

(19)

In (19), D = 2L is the hydraulic diameter, while the reference velocity U0 and the reference temperature T0 are given by ⎧ AD2 ⎪ ⎪ , ⎪ U0 = − ⎨ 48μ (20) ⎪  1  ⎪ + T T 1 ⎪ 1 2 ⎩ T0 = +S (T2 − T1 ). − 2 Bi1 Bi2

1020

J. C. UMAVATHI, J. PRATHAP KUMAR, and JAWERIYA SULTANA

The reference temperature ΔT is given either by ΔT = T2 − T1 or by ΔT =

ν2 cp D 2

if

if

T1 < T2

(21)

T1 = T2 .

(22)

Therefore, as in [4], the value of the dimensionless parameter RT can be either 0 or 1. More precisely, RT equals 1 for asymmetric fluid temperatures (T1 < T2 ), and equals 0 for symmetric fluid temperatures (T1 = T2 ). The dimensionless mean velocity u¯ and the dimensionless bulk temperature θb are given by

1/4

u ¯=2

udy,

(23)

−1/4

2 θb = u ¯



1/4

uθdy.

(24)

−1/4

On account of (13), for upward flow A < 0, U0 , Re, and λ are positive. For downward flow A > 0, U0 , Re, and λ are negative. By employing the dimensionless quantities defined in (19), (15)–(18) can be rewritten as  du 2 d4 u d2 u = λBr ∓ φ ∓ 48φ, dy 4 dy dy 2

(25)

u(−1/4) = u(1/4) = 0,

(26)

  1 d3 u  d2 u  − dy 2 y=−1/4 Bi1 dy 3 y=−1/4 4  RT λS  1+ , 2 Bi1   d2 u  1 d3 u  + dy 2 y=1/4 Bi2 dy 3 y=1/4

= −48 +

= −48 −

4  RT λS  1+ . 2 Bi2

(27)

(28)

Similarly, (14) and (19) yield  du 2 φ  d2 θ d2 u  48 + , = −Br ∓ dy 2 dy λ dy 2

(29)

while from (5) and (19), one obtains θ=−

1 d2 u  48 + 2 . λ dy

(30)

Equations (25)–(30) show that the dimensionless velocity profile and the dimensionless temperature profile depend on five parameters: the ratio λ = Gr/Re, the Brinkman number Br,

Mixed convection flow in vertical channel with boundary conditions

1021

the temperature difference ratio RT , the Biot numbers Bi1 and Bi2 , and the heat absorption/generation coefficient φ. The Nusselt numbers can be defined at each boundary as follows:  ⎧ dT  1 ⎪ ⎪ Nu1 = , ⎪ ⎪ RT (T (L/2) − T (−L/2)) + (1 − RT )ΔT dY Y =−L/2 ⎨  ⎪ ⎪ dT  1 ⎪ ⎪ Nu = . 2 ⎩ RT (T (L/2) − T (−L/2)) + (1 − RT )ΔT dY Y =L/2

(31)

The Nusselt numbers Nu1 and Nu2 can be employed to evaluate the heat fluxes at the walls. In fact, the heat flux per unit area is given by q1 = −k(dT /dY )|Y =− L 2

at the left wall and by q2 = −k(dT /dY )|Y = L 2

at the right wall. Let us first assume RT = 1. Then, from (31), one obtains ⎧ kNu1 ⎪ ⎪ ⎨ q1 = − D (T (L/2) − T (−L/2)),

(32)

⎪ ⎪ ⎩ q = − kNu2 (T (L/2) − T (−L/2)). 2 D The heat flux densities q1 and q2 can be also expressed as

q1 = −h1 (T (−L/2) − T1 ),

(33)

q2 = −h2 (T2 − T (L/2)). Equations (32) and (33) yield T (L/2) − T (−L/2) =

T2 − T1  k Nu1 1+ D h1 +

Nu2 h2

.

(34)

Let us now assume RT = 0. Equation (31) yields ⎧ kNu1 ⎪ ⎪ ⎨ q1 = − D ΔT, ⎪ ⎪ ⎩ q = − kNu2 ΔT, 2 D

(35)

where ΔT is the reference temperature difference defined by (22). By employing (19), (31) can be written as  ⎧ dθ  1 ⎪ ⎪ Nu1 = , ⎪ ⎨ RT (θ(1/4) − θ(−1/4) + (1 − RT )) dy y=−1/4  ⎪ dθ  1 ⎪ ⎪ . ⎩ Nu2 = RT (θ(1/4) − θ(−1/4) + (1 − RT )) dy y=1/4

(36)

1022

3

J. C. UMAVATHI, J. PRATHAP KUMAR, and JAWERIYA SULTANA

Solutions

3.1 Separated effects of buoyancy forces and viscous dissipation In this section, the simpler cases of either negligible viscous dissipation or negligible buoyancy forces will be solved analytically. The case of negligible viscous dissipation can be obtained by setting Br = 0 in the dimensionless energy equation (29). Equations (25)–(28) can be easily solved and yield   u = C1 + C2 y + C3 sinh( φy) + C4 cosh( φy) − 24y 2 (37) for the case of heat absorption and   u = C1 + C2 y + C3 sin( φy) + C4 cos( φy) − 24y 2

(38)

for the case of heat generation, respectively. With Bi1 = Bi2 = Bi, (37) and (38) can be rewritten as  3 u = C2 y + C3 sinh( φy) + − 24y 2 2

(39)

for the case of heat absorption and  3 u = C2 y + C3 sin( φy) + − 24y 2 2

(40)

for the case of heat generation, respectively. In the limit Bi → +∞, one obtains the special case in which the boundaries of the channel are kept at the temperatures T1 and T2 , respectively. In this limit, (19) yields S → 1 such that (37) and (38) reduce to (i.e., for S = 1)  3 u = C2 y + C3 sinh( φy) + − 24y 2 2

(41)

for the case of heat absorption and  3 u = C2 y + C3 sin( φy) + − 24y 2 2

(42)

for the case of heat generation, respectively. That is the velocity profile determined by Barletta[4]. By substituting (37) and (38) into (23) and (24), one obtains ⎧  √φ   √φ   4 ⎪ ⎪ u ¯ = C4 √ sinh − cosh + 1, ⎪ ⎪ 4 4 φ ⎨ (43) b1 + b2 + b3 + b4 + b5 ⎪ ⎪        θ = ⎪ √ √ b ⎪ ⎩ +1 λ C4 √4φ sinh 4φ − cosh 4φ for the case of heat absorption and ⎧  √φ   √φ   4 ⎪ ⎪ √ u ¯ = C − cos + 1, sin 4 ⎪ ⎪ 4 4 φ ⎨ b1 + b2 + b3 + b4 + b5 ⎪ ⎪  √   √  θ =   ⎪ ⎪ ⎩ b 4 √ +1 λ C4 φ sin 4φ − cos 4φ for the case of heat generation, respectively.

(44)

Mixed convection flow in vertical channel with boundary conditions

1023

Moreover, (30) and (36)–(38) yield   φ θ = − (C3 sinh( φy) + C4 cosh( φy)) λ

(45)

for the case of heat absorption and θ=

  φ (C3 sin( φy) + C4 cos( φy)) λ

for the case of heat generation, respectively. In the absence of heat absorption/generation equations, (30), (36), and (37) reduce to

θ = 2SRT y, Nu1 = Nu2 = 2RT ,

(46)

(47)

which are the same solutions obtained by Zanchini[5] . Plots of u versus y evaluated through (37) are reported in Fig. 2 for λ = 0, λ = 1 000, and Bi1 = Bi2 = 10 for different values of heat absorption coefficient φ.

Fig. 2

Plots of u versus y in case of RT = 1 for Br = 0 and Bi1 = Bi2 = 10

Let us now consider the case of negligible buoyancy forces with a relevant viscous dissipation, which corresponds to λ = Gr/Re = 0. Since a purely forced convection occurs in this case, the following Hagen-Poiseuille velocity profile is presented within the channel: u=

3 − 24y 2. 2

(48)

Indeed, for both symmetric and asymmetric fluid temperatures, (48) is the solution to (25)–(30) when λ = 0. Equations (9), (10), and (19) yield the boundary conditions on θ as follows:  ⎧  SR   −1  dθ  4  T ⎪ ⎪ 1 + , + θ = Bi 1 ⎪ ⎨ dy y=−1/4 2 Bi1 4 (49)   SR    1  ⎪  dθ 4 ⎪ T ⎪  1+ . −θ = Bi2 ⎩ dy  2 Bi2 4 y=1/4

On account of (29), (48), and (49), the temperature field can be expressed as    Br  2 2  θ = C5 sinh( φy) + C6 cosh( φy) + 482 y + φ φ

(50)

1024

J. C. UMAVATHI, J. PRATHAP KUMAR, and JAWERIYA SULTANA

for the case of heat absorption and    Br  2 2  y − θ = C5 sin( φy) + C6 cos( φy) − 482 φ φ

(51)

for the case of heat generation, respectively. Equations (31), (50), and (51) yield ⎧ √

√ √

√ φC5 cosh φ/4 − φC6 sinh φ/4 − 482 Br/(2φ) ⎪ ⎪

⎪ √ Nu1 = , ⎪ ⎨ 2RT C5 sinh φ/4 + (1 − RT ) √

√ √

√ ⎪ φC5 cosh φ/4 + φC6 sinh φ/4 + 482 Br/(2φ) ⎪ ⎪ ⎪ Nu2 = √

⎩ 2RT C5 sinh φ/4 + (1 − RT )

(52)

for the case of heat absorption and ⎧ √

√ √

√ φC5 cos φ/4 + φC6 sin φ/4 − 482 Br/(2φ) ⎪ ⎪ ⎪ √

Nu1 = , ⎪ ⎨ 2RT C5 sin φ/4 + (1 − RT ) √

√ √

√ ⎪ φC5 cos φ/4 + φC6 sin φ/4 + 482 Br/(2φ) ⎪ ⎪ ⎪ √

= Nu 2 ⎩ 2RT C5 sin φ/4 + (1 − RT )

(53)

for the case of heat generation, respectively. Equation (53) ensures that, for RT = 0, Nu1 > 0 and Nu2 < 0 for every value of Br. On the other hand, for RT = 1, Nu1 and Nu2 may become singular, and their sign depends on the values of Br, Bi1 , and Bi2 . Plots of θ versus y for RT = 1 (T1 < T2 ) evaluated through (51) are reported in Figs. 3 and 4 for some values of Br. Figure 3 refers to Bi1 = Bi2 = 10, while Fig. 4 refers to Bi1 = 1 and Bi2 = 10. Indeed, in Fig. 4, the plots for Br = 0 and Br = 1 are such that T (−L/2) > T (L/2).

Fig. 3

Plots of θ versus y in case of RT = 1 for λ = 0 and Bi1 = Bi2 = 10

Fig. 4

Plots of θ versus y in case of RT = 1 for λ = 0, Bi1 = 1, and Bi2 = 10

3.2 Combined effects of buoyancy forces and viscous dissipation In this section, both buoyancy forces and viscous dissipation are considered as non-negligible. First, (25)–(28) are solved by a perturbation series method. Then, the dimensionless temperature field is determined by means of (30). As in [4], let us consider the dimensionless

Mixed convection flow in vertical channel with boundary conditions

1025

parameter ε = λBr = Re P r

βgD , cp

(54)

which is independent of the reference temperature difference ΔT . For fixed values of RT , λ, Bi1 , Bi2 , and φ, the solution to (25)–(28) can be expressed by the perturbation expansion as follows: u (y) = u0 (y) + εu1 (y) + ε2 u2 (y) + · · · .

(55)

To obtain the solution to (25)–(28), with the form of (55), one first substitutes (55) into (25)–(28) and collects terms having power of ε. Then, one equates the coefficient of ε to zero[17] . Thus, one obtains a sequence of boundary value problems which can be solved in succession and yields the unknown functions un (y). The boundary value problems for n = 0 and n = 1 are d4 u0 d2 u0 = ∓φ ∓ 48φ, dy 4 dy 2 u0 (−1/4) = u0 (1/4) = 0,   d2 u0  1 d3 u0  − dy 2 y=−1/4 Bi1 dy 3 y=−1/4  R λS  4  T 1+ , 2 Bi1   d2 u0  1 d3 u0  + dy 2 y=1/4 Bi2 dy 3 y=1/4

= −48 +

 R λS  4  T 1+ , 2 Bi2 d4 u1 d2 u1  du0 2 = ∓φ + , dy 4 dy 2 dy

= −48 −

u1 (−1/4) = u1 (1/4) = 0,   1 d3 u1  d2 u1  − dy 2 y=−1/4 Bi1 dy 3 y=−1/4  R λS  4  T 1+ , 2 Bi1   d2 u1  1 d3 u1  + dy 2 y=1/4 Bi2 dy 3 y=1/4

= −48 +

= −48 −

 R λS  4  T 1+ , 2 Bi2

(56) (57)

(58)

(59) (60) (61)

(62)

(63)

where plus sign relates to heat absorption and minus sign relates to heat generation. The solutions to (56)–(63) are given by   u0 = C1 + C2 y + C3 sinh( φy) + C4 cosh( φy) − 24y 2 ,

(64)

1026

J. C. UMAVATHI, J. PRATHAP KUMAR, and JAWERIYA SULTANA

 u1 = C5 + C6 y + (C7 + l10 y + l11 y 2 )sinh( φy)  + (C8 + l9 y 2 + l12 y)cosh( φy)   + l7 cosh(2 φy)+l8 sinh(2 φy) + l13 y 4 + l14 y 3 + l15 y 2

(65)

for the case of heat absorption and   φy + C4 cos( φy) − 24y 2 ,

(66)

 u1 = C5 + C6 y + (C7 + l4 y 2 + l6 y)sin( φy)  + (C8 + l3 y 2 + l5 y)cos( φy)   + l1 cos(2 φy)+l2 sin(2 φy) + l7 y 4 + l8 y 3 + l9 y 2

(67)

u0 = C1 + C2 y + C3 sin



for the case of heat generation, respectively. The right-hand sides of (64) and (66) coincide with those of (37) and (38), respectively, and give the dimensionless velocity field in the case of Br = 0. The dimensionless temperature θ can be written in the following form:  1 ((C3 φ + ε(C7 φ + l17 y + l19 + l11 φy 2 ))sinh( φy) λ  + (C4 φ + ε(C8 φ + l18 y + l16 + l9 φy 2 ))cosh( φy)   + ε(4l7 φcosh(2 φy) + 4l8 φsinh(2 φy)+12l13 y 2 + 6l14 y + 2l15 ))

θ =−

(68)

for the case of heat absorption and  1 θ = ((C3 φ + ε(C7 φ + l4 φy 2 + l11 y − l13 )) sin( φy) λ  + (C4 φ + ε(C8 φ + l3 φy 2 − l10 y − l12 )) cos( φy)   + ε(4l1 φ cos(2 φy) + 4l2 φ sin(2 φy) − 12l7 y 2 − 6l8 y−2l9))

(69)

for the case of heat generation, respectively. Using the solutions given in (68) and (69) yields the following expressions for Nu1 and Nu2 : √

√

√ √ (C3 + εC7 )φ φcosh φ/4 − (C4 + εC8 )φ φsinh φ/4 + εn1





√ √ , Nu1 = RT 2φC3 sinh φ/4 + ε 2C7 φsinh φ/4 + a1 + λ(1 − RT )

(70)

√

√

√ √ (C3 + εC7 )φ φcosh φ/4 + (C4 + εC8 )φ φsinh φ/4 + εn2

√

√

Nu2 = RT 2φC3 sinh φ/4 + ε 2C7 φsinh φ/4 + a1 − λ(1 − RT )

(71)

for the case of heat absorption and √

√

√ √ (C3 + εC7 )φ φsin φ/4 + (C4 + εC8 )φ φcos φ/4 + εn1

√

√

, Nu1 = RT 2φC3 sin φ/4 + ε 2φC7 sin φ/4 + a1 + λ(1 − RT )

(72)

Mixed convection flow in vertical channel with boundary conditions

√

√

√ √ (C3 + εC7 )φ φsin φ/4 − (C4 + εC8 )φ φcos φ /4 + εn2 √

√

Nu2 = RT 2φC3 sin φ/4 + ε(2φC7 sin φ/4 + a1 ) + λ(1 − RT )

1027

(73)

for the case of heat generation, respectively. Equations (23), (55), (64), and (65) yield the expression of the mean dimensionless velocity as follows: u ¯ = (C4 + εC8 )((4/

   φ)sinh( φ/4) − cosh( φ/4))

+ 1 + ε(−(d1 + d2 )/2 + ub1 )

(74)

for the case of heat absorption. Equations (23), (55), (66), and (67) yield the expression of the mean dimensionless velocity as follows:    u ¯ = (C4 + εC8 )((4/ φ)sin( φ/4) − cos( φ/4)) + 1 + ε(−(d1 + d2 )/2 + ub )

(75)

for the case of heat generation. The constants n1 , n2 , a1 , a2 , d1 , d2 , b1 , and ub appeared in all the above solutions are not presented due to brevity. 3.3 Numerical solutions The analytical solutions obtained in the preceding section include only two terms of the series, which are not applicable for large values of λ, i.e., for large values of buoyancy force. In many practical problems, the values of λ are usually large. Therefore, a numerical scheme is used to solve the non-linear boundary value problem by using the Runge-Kutta shooting method. The validity of the numerical scheme is justified by comparing the analytical solutions with the numerical solutions, and it is found that this agreement is good enough.

4

Results and discussion

In this paper, we analyze the problem of mixed convection flow and heat transfer in a vertical channel with boundary conditions of the third kind. The analytical solutions are obtained by using the regular perturbation method with the product of λ(Gr/Re) and Br as the perturbation parameter valid for small values of perturbation parameters. The restriction on the perturbation parameters to be small is relaxed by finding the numerical solutions of the basic equations by using the Runge-Kutta shooting method. Two cases are considered depending on the thermal characteristics of the problem such as heat absorption and heat generation. The results are depicted graphically in Figs. 2–11 for the heat absorption case and in Figs. 12–16 for the heat generation case. Equation (37) for the velocity field u is evaluated for different values of λ and φ as shown in Fig. 2. When Br = 0 for equal values of Biot numbers, it is seen that the case for λ = 1 000 shows a flow reversal near the cold wall at y = −1/4 with U0 > 0, i.e., for upward flow. A flow reversal induced by the buoyancy forces occurs at the cold wall. By performing a reflection of the Y -axis of Fig. 2, plots of u for λ = −1 000 can be obtained. Further, it is also observed from Fig. 2 that as the heat absorption coefficient φ increases, the velocity increases in the cold wall and the velocity decreases near the hot wall. Equation (30) for the non-dimensional temperature θ is evaluated for different values of Brinkman number Br and heat absorption coefficient φ with λ = 0. We notice that the temperature field in linear indicats that the heat transfer is purely by conduction in the absence of dissipation (Br = 0). Also, the temperature field increases with the increase in the value of Br. In the absence of viscous dissipation, the temperature increases at the cold wall and decreases

1028

J. C. UMAVATHI, J. PRATHAP KUMAR, and JAWERIYA SULTANA

at the hot wall. The heat absorption coefficient decreases subsequently both at the cold and hot walls in the presence of viscous dissipation for equal Biot numbers. It is also observed from Fig. 3 that the temperature is almost linear in the middle of the channel for Br = 0 which indicates that the convection dominates in the boundary layer region. Figure 4 for unequal Biot numbers shows the similar effect of Brinkman number and heat absorption coefficient φ on the temperature field as that for equal Biot numbers except that the temperature field in the presence of viscous dissipation is effective near the cold wall compared to the hot wall. Figure 5 shows the variation of the non-dimensional velocity and temperature profiles u and θ for λ = ±500 with φ = 16 for different values of ε and λ with equal Biot numbers. It is seen that, for upward flow, the velocity and temperature at each position are increasing functions of ε. Moreover, the effect of ε on u is stronger for higher values of λ, while that on θ is weaker. It is also observed from Fig. 5 that, for downward flow at each position, u is a decreasing function for ε < 0, while θ is an increasing function for ε < 0. Further, it is also shown that for small values of λ, there is no flow reversal, while as for large values of λ, flow reversal occurs at both walls. This is because that enhancement of viscous dissipation results in higher values of temperature which in turn enhances the values of the buoyancy force. The increase of the buoyancy force yields an increase of the fluid velocity for λ > 0 and decreases the fluid velocity for λ < 0. Figures 2–5 compare the analytical solutions with the numerical solutions. We notice that the agreement between the numerical and analytical solution is good enough for small values of perturbation and the difference increases as ε increases. Further, it is also observed that the results obtained in Figs. 2–5 are in good agreement with the result reported by Zanchini[5] in the absence of heat absorption coefficient φ and for small values of ε as they use series solution method.

Fig. 5

Plots of u and θ versus y in case of RT = 1 for φ = 16 and Bi1 = Bi2 = 10

The effect of heat absorption coefficient φ on the flow field is shown in Fig. 6. An increase in the values of φ decreases the velocity field near the cold wall and increases the velocity field near the hot wall for downward flow (λ < 0), while the velocity increases the velocity field near the cold wall and decreases near the hot wall for upward flow (λ > 0). Similar effect of heat absorption is observed on temperature. However, the effect of φ on temperature is less operative when compared to the velocity. The effects of |ε| on Nu1 and Nu2 are the same as the results obtained by Zanchini[5] and hence not shown graphically. Figure 7 is drawn for the Nusselt number at both the walls for different values of φ and λ for equal Biot numbers. The Nusselt number is an increasing function of heat absorption coefficient φ at both the walls. However, as |λ| increases, the Nusselt number decreases at the cold wall, while as |λ| increases, the Nusselt number increases at the hot wall.

Mixed convection flow in vertical channel with boundary conditions

Fig. 6

Fig. 7

1029

Plots of u and θ versus y in case of RT = 1 for Bi1 = Bi2 = 10

Plots of Nu1 and Nu2 versus φ for heat absorption in case of RT = 1 for |ε| = 5 and Bi1 = Bi2 = 10

The effect of u¯ is the same as the result obtained by Zanchini[5] and hence not shown graphically. Figure 8 displays the effect of λ for various values of φ on the average velocity u ¯, which shows that the larger value of heat absorption coefficient φ reduces the average velocity for upward flow, while it increases for the downward flow.

Fig. 8

Plots of u ¯ versus φ in case of RT = 1 for |ε| = 5 and Bi1 = Bi2 = 10

Figure 9 shows the plots for velocity and temperature for unequal Biot numbers for different values of ε. It shows that as ε increases, both velocity and temperature increase, but the effect

1030

J. C. UMAVATHI, J. PRATHAP KUMAR, and JAWERIYA SULTANA

of ε on temperature is more operative at the cold wall. Figure 9 tells that fixing ε and λ, the temperature variations are significant for small values of Bi which is again the similar result obtained by Zanchini[5] .

Fig. 9

Plots of u and θ versus y in case of RT = 1 for λ = 300, φ = 9, Bi1 = 0.1, and Bi2 = 10

Figures 10 and 11 show the velocity and λθ for symmetric wall heating conditions for equal and unequal Biot numbers. The effect of ε on the flow field shows that both u and λθ are increasing functions of ε for equal and unequal Biot numbers for symmetric wall heating conditions. Also, the velocity and temperature profiles are more significant in the case of upward flow (ε > 0) than in the case of downward flow (ε < 0) and become stronger if either Bi1 or Bi2 becomes smaller (similar results obtained by Zanchini[5] ). It is also observed in Figs. 10 and 11 that for ε > 0, the flow field is reduced for large values of heat absorption coefficient φ, while as for ε < 0, the flow field increases for large values of φ for both equal and unequal Biot numbers.

Fig. 10

Plots of u and λθ versus y in case of RT = 0 for Bi1 = Bi2 = 10

It is found that the effect of various non-dimensional parameters on velocity and temperature in the case of heat generation is qualitatively the same as that in the case of heat absorption. However, we discuss the results for heat generation where there are deviations from heat absorption. The effect of heat generation coefficient φ on the flow field is to increase the velocity and temperature field as φ increases for upward flow at both the walls and to decrease the flow field for downward flow at both the walls as seen in Fig. 12. However, the flow reversal is observed at the left wall for upward flow and at the right wall for downward flow. Figure 13 shows the effect of heat generation coefficient φ on the Nusselt numbers at both the walls. The Nusselt number at the cold wall is a decreasing function of heat generation

Mixed convection flow in vertical channel with boundary conditions

1031

coefficient φ for downward flow, while it is an increasing function of φ for upward flow. The Nusselt number at the hot wall is a decreasing function of heat generation coefficient φ at both cold and hot walls. However, as |λ| decreases, the Nusselt number at the cold wall decreases, while as |λ| increases, the Nusselt number increases at the hot wall.

Fig. 11

Fig. 12

Fig. 13

Plots of u and λθ versus y in case of RT = 0 for Bi1 = 0.1 and Bi2 = 10

Plots of u and θ versus y in case of RT = 1 for some values of ε and λ, and Bi1 = Bi2 = 10

Plots of Nu1 and Nu2 versus φ for heat generation in case of RT = 1 for |ε| = 5 and Bi1 = Bi2 = 10

1032

J. C. UMAVATHI, J. PRATHAP KUMAR, and JAWERIYA SULTANA

Figure 14 shows that the average velocity is an increasing function for upward flow by increasing values of φ and a decreasing function for downward flow by increasing values of φ. It is also observed that the average velocity increases as λ increases for upward flow, while it decreases as λ increases for downward flow. Figure 15 displays the effect of ε on the flow field for unequal Biot numbers for upward flow. It shows that as ε increases, both the velocity and temperature fields decrease. However, the effect of ε is more significant at the cold wall. Figure 16 shows the flow field for symmetric wall heat condition for unequal Biot number for variation of heat generation coefficient φ. The velocity and temperature profiles are more significant in the case of downward flow (ε < 0) than in the case of upward flow (ε > 0) and become stronger if either Bi1 or Bi2 becomes stronger.

Fig. 14

Fig. 15

5

Plots of u ¯ versus φ in case of RT = 1 for Bi1 = Bi2 = 10

Plots of u and θ versus y in case of RT = 1 for Bi1 = 0.1 and Bi2 = 10

Conclusions

The effect of viscous dissipation and heat source/sink on a fully developed mixed convection for the laminar flow in a parallel-plate vertical channel has been analyzed. The boundary condition of the convective heat exchange with an external fluid at each boundary plane is considered. The simple cases of negligible viscous dissipation or negligible buoyancy forces are solved analytically and numerically. Then, the combined effects of viscous dissipation and buoyancy forces in the presence of source/sink are analyzed. Both the cases of asymmetric fluid temperature (RT = 1) with equal and different Biot numbers and the case of symmetric fluid temperature (RT = 0) with equal and different Biot numbers are considered. The velocity and temperature distributions are discussed numerically. For upward flow, the velocity and temperature at each

Mixed convection flow in vertical channel with boundary conditions

Fig. 16

1033

Plots of u and λθ versus y in case of RT = 1 for Bi1 = 0.1 and Bi2 = 10 when φ = 20 or φ = 25

position are increasing functions of |ε|. For downward flow at each position, the velocity is a decreasing function of |ε|, and the temperature is an increasing function of |ε|. The effect of heat absorption on the flow field is to decrease the velocity and temperature as φ increases. The effect of heat generation on the flow field is to increase the velocity and temperature as φ increases. The effect of u and θ increases when at least one of the Biot number decreases. Acknowledgements The authors thank University Grant Commission in New Delhi for the financial support under UGC-Major Research Project and Maulana Azad National Fellowship for Minority Students.

References [1] Aung, W. and Worku, G. Theory of fully developed, combined convection including flow reversal. ASME Journal of Heat Transfer, 108, 485–588 (1986) [2] Hamadah, T. T. and Wirtz, R. A. Analysis of laminar fully developed mixed convection in a vertical channel with opposing buoyancy. ASME Journal of Heat Transfer, 113, 507–510 (1991) [3] Cheng, C. H., Kou, H. S., and Huang, W. H. Flow reversal and heat transfer of fully developed mixed convection in vertical channels. Journal of Thermophysics and Heat Transfer, 4, 375–383 (1990) [4] Barletta, A. Laminar mixed convection with viscous dissipation in a vertical channel. International Journal of Heat and Mass Transfer, 41, 3501–3513 (1998) [5] Zanchini, E. Effect of viscous dissipation on mixed convection in a vertical channel with boundary conditions of the third kind. International Journal of Heat and Mass Transfer, 41, 3949–3959 (1998) [6] Sparrow, E. M. Analysis of Laminar Forced-Convection Heat Transfer in the Entrance Region of Flat Rectangular Ducts, Washington Technical Note, National Advisory Committee for Aeronautics, Vol. 3331 (1955) [7] Stephan, K. Warmenbergang a Druckabfall Bei nicht ausgebildeter laminarstromung in Rohren und in ebenen spalten. Chemie Ingeniear Technik, 31, 773–778 (1959) [8] Hwang, C. L. and Fan, L. T. Finit difference analysis of forced convection heat transfer in entrance region of a flat rectangular duct. Applied Scientific Research, 13, 401–422 (1964) [9] Siegal, R. and Sparrow, E. M. Simultaneous development of velocity and temperature distributions in a flat duct with uniform wall heating. American Institute of Chemical Engineers Journal, 5, 73–75 (1959) [10] Javeri, V. Heat transfer in laminar entrance region of a flat channel for the temperature boundary condition of the third kind. W¨ arme-und Stoff¨ ubertragung, 10, 137–144 (1977) [11] Javeri, V. Laminar heat transfer in a rectangular channel for the temperature boundary condition of the third kind. International Journal of Heat and Mass Transfer, 10, 1029–1034 (1978)

1034

J. C. UMAVATHI, J. PRATHAP KUMAR, and JAWERIYA SULTANA

[12] Kumari, M. and Nath, G. Mixed convection boundary layer flow over a thin vertical cylinder with localized injection/suction and cooling/heating. International Journal of Heat and Mass Transfer, 47, 969–976 (2004) [13] Umavathi, J. C., Mallikarjun, P. B., and Pop, I. On laminar mixed convection flow in a vertical porous stratum with symmetric wall heating conditions. International Journal of Transport Phenomena, 8, 127–140 (2006) [14] Mahanthi, N. C. and Gaur, P. Effect of varying and thermal conductivity on steady free convective flow and heat transfer along an isothermal vertical plate in the presence of heat sink. Journal of Applied Fluid Mechanics, 2(1), 23–28 (2009) [15] Umavathi, J. C. and Prathap Kumar, J. Mixed convection flow of micropolar fluid in a vertical channel with symmetric and asymmetric wall heating conditions. International Journal of Applied Mechanics and Engineering, 16(1), 141–159 (2010) [16] Arpaci, V. S. and Larsen, P. S. Convection Heat Transfer, Pentice-Hall/Englewood Cliffs, New Jersey, 51–54 (1984) [17] Aziz, A. and Na, T. Y. Perturbation Methods in Heat Transfer, Hemisphere Publishing Corporation, New York (1984)