# Natural Convection Flow along an Isothermal Vertical Flat Plate with Temperature Dependent Viscosity and Heat Generation

Authors: Md. Mamun Molla; Anita Biswas; Abdullah Al-Mamun; Md. Anwar Hossain
Journal: Journal of Computational Engineering
 (2014)

---

## Abstract

The purpose of this study is to investigate the natural convection laminar flow along an isothermal vertical flat plate immersed in a fluid with viscosity which is the exponential function of fluid temperature in presence of internal heat generation. The governing boundary layer equations are transformed into a nondimensional form and the resulting nonlinear system of partial differential equations is reduced to a convenient form which are solved numerically using an efficient marching order implicit finite difference method with double sweep technique. Numerical results are presented in terms of the velocity and temperature distribution of the fluid as well as the heat transfer characteristics, namely, the wall shear stress and the local and average rate of heat transfer in terms of the local skin-friction coefficient, the local and average Nusselt number for a wide range of the viscosity-variation parameter, heat generation parameter, and the Rayleigh number. Increasing viscosity variation parameter and Rayleigh number lead to increasing the local and average Nusselt number and decreasing the wall shear stress. Wall shear stress and the rate of heat transfer decreased due to the increase of heat generation.


---

## Body

## 1. Introduction

A large number of physical phenomena involve natural convection driven by heat generation. The study of heat generation in moving fluids is important in view of several physical problems such as those dealing with chemical reactions and those concerned with dissociating fluids. Possible heat generation effects may alter the temperature distribution and, therefore, the particle deposition rate. This may occur in such applications related to nuclear reactor cores, fire and combustion modeling, electronic chips, and semiconductor wafers. In fact, the literature is replete with examples dealing with the heat transfer in laminar flow of viscous fluids. Vajravelu and Hadjinicolaou [1] studied the heat transfer characteristics in the laminar boundary layer of a viscous fluid over a linearly stretching continuous surface with viscous dissipation or frictional heating and internal heat generation. In this study, Vajravelu and Hadjinicolaou [1] considered that the volumetric rate of heat generation, q′′ [W/m3], should be
(1)q′′={Q0(T--T∞),forT-≥T∞,0,forT-<T∞,
where Q0 is the heat generation constant. The above relation, explained by Vajravelu and Hadjinicolaou [1], is valid as an approximation of the state of some exothermic process and having T∞ as the onset temperature. Following Vajravelu and Hadjinicolaou [1], Molla et al. [2–5] investigated the natural convection with heat generation along a vertical wavy surface, horizontal circular cylinder, and sphere, respectively.In all of the above studies, the authors assumed that the viscosities of the fluids are constant throughout the flow regime. But this physical property may change significantly with temperature. For example, the viscosity of air is 0.6924 × 10−5 kg/m·s, 1.3289 kg/m·s, 2.286 kg/m·s, and 3.625 kg/m·s at 100 K, 200 K, 400 K, and 800 K temperature, respectively. (See Cebeci and Bradshaw [6].) In order to predict accurately the flow behaviour, it is necessary to take into account the temperature dependence of viscosity. Gary et al. [7] and Mehta and Sood [8] found that the flow characteristics change substantially when the effect of temperatur-dependent viscosity was considered. Hady et al. [9] have also investigated the mixed convection boundary layer flow on a continuous flat plate with variable viscosity. Kafoussias and Williams [10] and Kafoussias et al. [11] have studied the effects of variable viscosity on the free and mixed convection flow from a vertical flat plate in the region near the leading edge. Hossain et al. [12, 13] and Molla et al. [14] considered the natural convection with variable viscosity, where the viscosity is inversely proportional to a linear function of temperature, a model that was proposed by Ling and Dybbs [15]. Besides these, it has been found that for liquid, such as organ gas, the viscosity varies with temperature in an approximately linear manner. On assuming that the viscosity of the fluid is linear functions of temperature, a semiempirical formula was proposed by Charraudeau [16]. Following him Molla et al. [17, 18] have investigated the effects of chemical reaction, heat, and mass diffusion on natural convection flow from an isothermal sphere and circular cylinder with temperature-dependent viscosity.Natural convection flows were analyzed applying a new method in order to obtain the approximate equations by Gray and Giorgini [19]. The established Boussinesq approximation was verified under this method in the case of a given Newtonian liquid or gas. Considering the constant property (CP) model and the variable property (VP) model the laminar natural convection flow along a horizontal rectangular duct taking into account the significant buoyancy effects was studied by Wang et al. [20]. The Boussinesq approximations were introduced for the CP model and the properties were used at four different reference temperatures. On the other hand, the density and the transport properties were worked out with the help of the state equation of an ideal gas and power law correlations.Saha et al. [21] investigated the natural convection flow along an isothermal vertical plate surrounded in a stratified medium with uniform heat source. The numerical solutions of the nonsimilar equations were carried out by the implicit finite difference method and it is mentioned that the same method will be used in the present investigation. The nonsimilar laminar free convection adjacent to a vertical wall applying exceptional prescribed boundary conditions was analyzed by Kao [22]. The locally nonsimilar method for the second order truncation was employed in the study in order to investigate the nonsimilar free convection problems. The numerical simulations of the buoyancy forced flow over an isothermal vertical surface having a leading edge were carried out by Wright and Gebhart [23]. In this study the coordinate transformations were taken under consideration in order to allow efficient calculation of both the boundary layer flow and the extensive ambient entrainment in-flow.But thermal convection in a fluid whose viscosity is strongly temperature dependent is of interest in geophysical problems as well as in various engineering applications. Numerous investigations have been focused on theoretical and experimental aspects of this problem, such as those of Torrance and Turcotte [24], Booker [25], Stengel et al. [26], Richter et al. [27], and Bottaro et al. [28]. Torrance and Turcotte [24] have investigated the influence of large variations of viscosity on convection in a layer of fluid heated from below. Solutions for the flow and temperature fields were obtained numerically assuming infinite Prandtl number, free-surface boundary conditions, and two-dimensional motion of fixed horizontal wavelength. The effects of a temperature-dependent and a depth-dependent viscosity were each studied. Thermal convection in a fluid with viscosity which is strong function of temperature may be directly applicable to the earth’s mantle. So, it is known that the viscosity of the mantle can be formulated as the exponential function of temperature.The effect of the heat generation on natural convection flow along a vertical flat plate with the exponentially varying temperature-dependent kinematic viscosity has not been studied yet and the present work demonstrates the issue. For the experimentalists, this investigating numerical result helps a lot to predict the flow behaviour. In this investigation the focus is on the boundary layer regime promoted by the combined events along a flat plate with temperature-dependent viscosityΠ(T) (as in Torrance and Turcotte [24]) in presence of heat generation while the surface is at a uniform temperature. The basic equations are transformed to the convenient form of the boundary layer equations, which are solved numerically using a very efficient marching order implicit finite-difference method. Consideration is given to the situation where the buoyancy forces assist the natural convection flow for various combinations of the viscosity-variation parameter ε, heat generation parameter λ, and the Rayleigh number Ra. The results allow us to predict the different behavior that can be observed when the relevant parameters are varied.


## 2. Formulation of Problem

Consider a steady two-dimensional laminar natural convective flow from a uniformly heated semi-infinite vertical flat plate, which is immersed in a viscous and incompressible fluid having temperature-dependent kinematic viscosity. Here kinematic viscosity is the exponential function of the fluid temperature. It is assumed that the surface temperature of the flat plate isTw, where Tw>T∞. Here T∞ is the ambient temperature of the fluid; the configuration considered is as shown in Figure 1. In Figure 1, δM and δT represent the momentum and thermal boundary layer thickness, respectively.Figure 1
Physical model and coordinate system.The equations governing the flow are(2)∂u-∂x-+∂v-∂y-=0,u-∂u-∂x-+v-∂u-∂y-=∂∂y-(ν∂u-∂y-)+gβ(T--T∞),u-∂T-∂x-+v-∂T-∂y-=α∂2T-∂y-2+Q0ρCp(T--T∞).The boundary conditions of (2) are(3a)u-=v-=0,T-=Tw,aty-=0,(3b)u-⟶0,T-⟶T∞asy-⟶∞,where(u-,v-) are velocity components along the (x-,y-) axes, g is the acceleration due to gravity, ρ is the density, ν is the dimensional kinematic viscosity of the fluid depending on the fluid temperature T-, β is the coefficient of thermal expansion, α is the thermal diffusivity of the fluid, and Q0 being a constant, which may take either positive or negative. The source term represents the heat generation when Q0>0 and the heat absorption when Q0<0 and Cp is the specific heat at constant pressure.We now introduce the following nondimensional variables:(4)x=x-l,y=y-l,u=lαu-,v=lαv-,θ=T--T∞Tw-T∞,Π=νν0Ra=gβ(Tw-T∞)l3αν0=GrPr,Pr=ν0α,Gr=gβ(Tw-T∞)l3ν02,
where ν0 is the reference kinematic viscosity at the ambient, l is the characteristic length of the vertical flat plate, Gr is the Grashof number, Pr is the Prandtl number, Ra is the Rayleigh number, and T is the nondimensional temperature.Out of the many forms of kinematic viscosity variation, which are available in the literature, we will consider only the following form proposed by Torrance and Turcotte [24]:
(5)ν=ν0exp⁡[ε(12-T--T∞Tw-T∞)],
where ε is the viscosity-variation parameter. When ε=0.0, the kinematic viscosity of the fluid is constant. If ε<0.0, the kinematic viscosity is larger than the constant kinematic viscosity and if ε>0.0, the kinematic viscosity is smaller than the case of ε=0.0.Substituting variables (4) into (2) leads to the following nondimensional equations:
(6)∂u∂x+∂v∂y=0,u∂u∂x+v∂u∂y=PrΠ[∂2u∂y2-ε∂u∂y∂T∂y]+RaPrθ,u∂θ∂x+v∂θ∂y=∂2θ∂y2+λθ,
where Π is the form of variable kinematic viscosity and λ is the heat generation parameter defined, respectively, as
(7)Π=exp⁡⁡[ε(12-θ)],λ=Q0l2ραCp.The boundary conditions (3a) and (3b) become(8a)u=v=0,θ=1aty=0,(8b)u⟶0,θ⟶0asy⟶∞.Here we introduce new transformations for the numerical scheme(9)X=x,Y=yx1/4,U=ux1/2V=x1/4v.Using (9) into (6), we get
(10)X∂U∂X+12U-14Y∂U∂Y+∂V∂Y=0,(11)XU∂U∂X+(V-14YU)∂U∂Y+12U2=PrΠ[∂2U∂Y2-ε∂U∂Y∂θ∂Y]+RaPrθ,(12)XU∂θ∂X+(V-14YU)∂θ∂Y=∂2θ∂Y2+X1/2λθ.The corresponding boundary conditions are(13a)U=V=0,θ=1atY=0,(13b)U⟶0,θ⟶0asY⟶∞.Now (10)–(12) subject to the boundary conditions (13a) and (13b) are discretised for numerical scheme using central-difference for diffusion terms and the forward-difference for the convection terms; finally we get a system of tridiagonal algebraic equations. The algebraic equations have been solved by Gaussian elimination technique. In computation we directly solve the continuity equation for the normal velocity V by the following discretisation:
(14)Vi,j=Vi,j-1+Yj14(Ui,j-Ui,j-1)-ΔY14(Ui,j-1+Ui,j)-XiΔYΔX(Ui,j-Ui-1,j).The computation is started atX=0.0 and then marches up to the point (X=4.0). Here Δx=Δy=0.01 are used for the X- and Y-grids, respectively. The physical quantities, namely, the local skin-friction coefficient and the local Nusselt number which are important from the application point of view, are calculated from the following dimensionless relations:
(15)CfGr1/4=exp⁡(-ε2)X1/4(∂U∂Y)Y=0,(16)NuGr-1/4=-X-1/4(∂θ∂Y)Y=0.For calculating the stream function Ψ, the fluid velocity over the whole boundary layer has been integrated, which can be written as(17)Ψ=∫0YX1/2UdY.The average Nusselt number is calculated from the following relation:(18)NumGr-1/4=-1l∫0lX-1/4(∂θ∂Y)Y=0dx.


## 3. Results and Discussion

Equations (10)–(12) subject to the boundary conditions (13a) and (13b) are solved numerically using marching order implicit finite-difference method. The numerical solutions are started at X=0.0 and proceeded to the downstream region. Here it should be mentioned that the solutions of the parabolic equations (10)–(12) can be terminated anywhere at X>0. In this investigation we have taken solutions up to X=4. Solutions are obtained for the Rayleigh number Ra (50.0, 100.0, and 200.0), heat generation parameter λ (=0.0, 2.5, 5.0), and for a wide range of values of the variable viscosity-variation parameter ε (=−2.0, −1.0, 0.0, 1.0, 2.0, 3.0, and 5.0). Numerical values of the shearing stress in terms of the local skin-friction coefficient CfGr1/4 from (15) and the heat transfer rate in terms of the local Nusselt number Nu Gr−1/4 from (16) are calculated. Numerical values of CfGr1/4 and Nu Gr−1/4 are depicted in Figures 2–4, 8, and 10. For the whole computations Pr is fixed at 1.Grid independence test in terms of (a) the local skin-friction coefficient and (b) Nusselt number while Pr = 1.0, Ra = 200,ε=5, and λ=2.5.
(a)
(b)(a) Skin-friction coefficient. (b) Rate of heat transfer for different values ofε while Ra = 100.0 and λ=0.0.
(a)
(b)(a) Skin-friction coefficient. (b) Rate of heat transfer for different values ofε while Ra = 100.0 and λ=2.5.
(a)
(b)In order to validate our numerical results the comparison with the published results is not straightforward, because the present boundary layer formulation is a little bit different from the available literature. For comparison, following Ede [29], it has been assumed Ra = 1 in (11) and (1/Pr) accompanied with the diffusion term (∂2θ/∂Y2) of (12). For ε=λ=0.0, the numerical results of the modified Nusselt number Nu(Gr/X)−1/4 are presented in Table 1 for the three Prandtl numbers along with the results of Ede [29] and Huang et al. [30]. The comparison shows that the present results obtained using marching order implicit finite difference method are in an excellent agreement with the solutions of Ede [29] and Huang et al. [30]. A grid independence test has been conducted for three different grid arrangements as Case 1: 200 × 1000 (x×y), Case 2: 400 × 2000, and Case 3: 800 × 4000 which is shown in Figure 2. The results are shown in terms of the local skin-friction coefficient as well as the local Nusselt number and the agreement between the results is quite excellent. The whole simulation has been done using the grid arrangement as in Case 2.Table 1
Comparison of the Nusselt number Nu(Gr/X)−1/4 with Ede [29] and Huang et al. [30].
Pr
Ede [29]
Huang et al. [30]
Present results
10
0.8269
0.8268
0.82683
100
1.5506
1.5486
1.54946
1000
2.8047
2.8084
2.80356Figure3 shows the effect of the viscosity-variation parameter ε (=−2.0, −1.0, 0.0, 1.0, and 2.0) on the local skin-friction coefficient CfGr1/4 and the local rate of heat transfer Nu Gr−1/4 without having any heat generation. For ε<0.0, the local skin-friction coefficient increases and the local Nusselt number decreases. Because, for ε<0.0, the kinematic viscosity of the fluid increases this increases wall shear stress. On the other hand, for ε>0.0, the local skin-friction coefficient decreases and the local Nusselt number increases due to small kinematic viscosities of the fluid. For ε≫0.0, the kinematic viscosities are very small which produce zero wall shearing stress. This is an impractical situation to get an inviscid boundary layer. In the present investigation for ε>5, the values of the kinematic viscosity are approximately zero. For ε≪0.0, the kinematic viscosities are so high which are appropriate for the highly non-Newtonian shear-thickening fluids.The effect of the viscosity-variation parameter,ε (=0.0, 1.0, 3.0, and 5.0) having heat generation effects λ (=2.5) on the local skin-friction coefficient CfGr1/4 and the local rate of heat transfer Nu Gr−1/4, is illustrated in Figures 4(a)-4(b), respectively, while Ra = 100.0. From these figures it is seen that the skin-friction coefficient decreases and the rate of heat transfer increases while the viscosity-variation parameter increases. It is happening due to the fact that neat surface viscosity of the fluid increases and, for high viscous fluid, the velocity gradient decreases and hence the skin-friction coefficient decreases. For high viscosity of the fluid, the temperature distribution decreases within the boundary layer (see Figure 7). Owing to the temperature difference between the plate and fluid, the rate of heat transfer increases significantly. For example at X=2.0, the local Nusselt number Nu Gr−1/4 increase 183.12% as increasing values of ε from 0 to 5.Figures5–7 show the isolines viscosity, streamlines, and isotherms, respectively, for different values of the viscosity-variation parameter while Ra = 100.0 and λ=2.5. From Figures 5(a)–5(c), it was clearly seen that the viscosity of the fluid increased at the vicinity of the surface which indicates that the viscosity of the fluid is strongly dependent on temperature. Due to high viscosity of the fluid, the velocity distribution decreases within the boundary layer and hence the values of the stream function Ψ decrease which are shown by the dashed lines in Figure 6. Also the temperature distribution reduces slightly for large values of ε. Finally, it can be concluded that the momentum and thermal boundary layer become thin for high viscous fluid and which is expected.Isolines viscosity for (a)ε=1.0, (b) ε=3.0, and (c) ε=5.0 while Ra = 100.0 and λ=2.5.
(a)
(b)
(c)Streamlines (a) solid forε=0.0 and dashed for ε=1.0 and (b) solid for ε=3.0 and dashed for ε=5.0 while Ra = 100.0 and λ=2.5.
(a)
(b)Isotherms (a) solid forε=0.0 and dashed for ε=1.0 and (b) solid for ε=3.0 and dashed for ε=5.0 while Ra = 100.0 and λ=2.5.
(a)
(b)(a) Skin-friction coefficient. (b) Rate of heat transfer for different values ofλ while Ra = 100.0 and ε=3.0.
(a)
(b)Figures8(a)-8(b) show the local skin-friction coefficient CfGr1/4 and the local rate of heat transfer Nu Gr−1/4 for the different values of the heat generation parameter λ (=0.0, 2.5, and 5.0) while Ra = 100.0 and ε=3.0. From these figures, it is seen that the increase of heat generation parameter λ leads to decreasing the local skin-friction coefficient CfGr1/4 and the local Nusselt number Nu Gr−1/4. These are expected, since the heat generation mechanism creates a layer of hot fluid near the surface and finally the resultant temperature of the fluid exceeds the surface temperature (see isotherms in Figure 9 for λ=5.0). For this reason the rate of heat transfer from the surface decreases. Owing to enhanced temperature, viscosity of the fluid increases and the corresponding local skin-friction coefficient decreases slightly, which is also expected.Figure 9
Isolines viscosity (a), streamlines (b), and isotherms (c) for different values ofλ while Ra = 100.0 and ε=3.0.(a) Skin-friction coefficient. (b) Rate of heat transfer for different values of Ra whileλ=2.5 and ε=3.0.
(a)
(b)Figure9 illustrates the effect of the heat generation parameter λ on the development of isolines viscosity (a), streamlines (b), and isotherms (c), which are plotted for Ra = 100.0 and ε=3.0. For large heat generation parameter the viscosity of the fluid increases which are shown in the figures of the isolines of viscosity. From Figure 9(b), it is seen that without effect of heat generation (i.e., λ=0.0) the nondimensional value of Ψmax⁡ within the computational domain is about 0.8 when the boundary layer thickness is the highest, but Ψmax⁡ increases with the increment of λ and it attains about 0.87 and 0.95 for λ=2.5 and 5.0, respectively. This phenomenon fully coincides with the early discussion made on Figure 8(a); the fluid speeds up as λ increases and the thickness of the velocity boundary layer also grows. The isotherm patterns for corresponding values of λ are shown in Figure 9(c). From all these fames, it is seen that the growth of thermal boundary layer over the surface of the flat plate is significant. As X increases from the leading edge (X≈0.0), the hot fluid raises and hence the thickness of the thermal boundary layer increases. But this phenomenon is not very straightforward as can be seen in 2nd and 3rd frames of Figure 9(c) for λ=2.5 and for λ=5.0 that the levels of isotherm are noticeably higher than the surface level, and in 3rd frame the fluid temperature exceeds the surface level and the surface heat transfer rate became negative which was also noticed in Figure 8(b).The effect of different values of the Rayleigh Ra on skin friction coefficientCfGr1/4 and the rate of heat transfer Nu Gr−1/4 with the heat generation parameter λ=2.5 and for ε=3.0 are illustrated in Figures 10(a)-10(b), respectively. It can easily be seen that the effect of the Rayleigh number Ra leads to increasing the skin friction coefficient CfGr1/4 and the rate of heat transfer Nu Gr−1/4. This phenomenon can easily be understood from the fact that with the increasing values of the Rayleigh number Ra, the buoyancy force increase, which speeds up the fluid velocity (see Figure 11(b)), that means, increases the velocity gradient and hence the local skin-friction coefficient CfGr1/4 enhances. Owing to the increase in the values of Ra leads to decrease the fluid temperature (see in Figure 11(c)) within the boundary layer and the associated thermal boundary layer becomes thinner. For decreasing fluid temperature, the temperature difference between fluid and surface increases and the corresponding rate of heat transfer Nu Gr−1/4 increases.(a) Viscosity, (b) velocity, and (c) temperature profiles for different values of Ra whileλ=2.5 and ε=3.0 at X=1.0.
(a)
(b)
(c)Figures12(a)–12(c) depict the average rate of heat transfer in terms of the average Nusselt number with the effect of the viscosity variation, Rayleigh number, and the heat generation, respectively. Average Nusselt number increases with the increase of viscosity-variation parameter and the Rayleigh number. On the other hand, with increase of heat generation the average Nusselt number decrease like local Nusselt number.Average Nusselt number for (a) viscosity-variation parameter while Ra = 200 andλ=2.5 (b) heat generation parameter while Ra = 100 and ε=3.0 (c) Rayleigh number while λ=2.5 and ε=3.0.
(a)
(b)
(c)


## 4. Conclusion

Numerical solutions for the steady laminar free convection boundary layer flow along a vertical flat plate subjected to a uniform surface temperature in presence of heat generation effect with fluid having viscosity which is the exponential function of temperature have been investigated theoretically. The governing boundary layer equations of motion are transformed into a nondimensional form and the resulting nonlinear systems of partial differential equations are reduced to convenient boundary layer equations, which are solved numerically by using marching order implicit finite difference method. From the present investigation the following conclusions may be drawn.(i)
With effect of viscosity-variation parameterε and the Rayleigh number Ra, the skin-friction coefficient decreases and the local and average rate of heat transfer increase.(ii)
An increase in values ofλ leads to slightly decreasing the skin-friction coefficient CfGr1/4 but decrease of the local and average rate of heat transfer is more significant.(iii)
For increase values of viscosity-variation parameterε, the momentum and thermal boundary layer become thinner.(iv)
With the effect of the Rayleigh number Ra both the viscosity and velocity distribution increase and the temperature distributions decrease significantly and the thickness of the momentum boundary layer is enhanced.



---

*Source: 712147-2014-05-27.xml*