Conjugate mixed convection in the entrance region of a symmetrically heated vertical channel with thick walls

Conjugate mixed convection in the entrance region of a symmetrically heated vertical channel with thick walls

International Journal of Thermal Sciences 98 (2015) 245e254 Contents lists available at ScienceDirect International Journal of Thermal Sciences jour...

2MB Sizes 0 Downloads 31 Views

International Journal of Thermal Sciences 98 (2015) 245e254

Contents lists available at ScienceDirect

International Journal of Thermal Sciences journal homepage: www.elsevier.com/locate/ijts

Conjugate mixed convection in the entrance region of a symmetrically heated vertical channel with thick walls G. Yang, J.Y. Wu* Institute of Refrigeration and Cryogenics, Key Laboratory for Power Machinery and Engineering of M.O.E, Shanghai Jiao Tong University, Shanghai, 200240, China

a r t i c l e i n f o

a b s t r a c t

Article history: Received 29 June 2014 Received in revised form 17 July 2015 Accepted 17 July 2015 Available online 10 August 2015

Conjugate mixed convection with buoyancy assisted laminar flow in the entrance region of a vertical channel is considered numerically. The problem is solved by a finite volume method for a thick walled, two-regional channel which has constant and uniform outside wall temperatures. The effects of wall thermal conduction as well as assisted buoyancy force on the flow and heat transfer are discussed in detail. Results are presented for a Prandtl number of 0.7, solid-to-fluid thermal conductivity ratios of 1  k* < ∞, wall thickness-to-channel length ratios of 0  l*  5, Reynolds numbers of 200  Re  1000, and for various Grashof numbers. The critical buoyancy parameter (Gr/Re), above which the flow reversal occurs, increases linearly with the increasing l*/k*, while it is independent on the Reynolds number. © 2015 Elsevier Masson SAS. All rights reserved.

Keywords: Mixed convection Conjugate heat transfer Flow reversal Vertical channel

1. Introduction Mixed convection heat transfer in vertical pipe and channel flows has been extensively studied in the past few decades due to its importance in industrial and engineering applications, such as heat exchanger systems, nuclear reactors, electronic cooling, fluid transport, building works and so on. It has been recognized that, when the flow velocity is low and the temperature difference between the channel wall and the fluid is large, the direction and magnitude of the buoyancy force may significantly affect the flow structure and heat transfer characteristics in the channel. Most of the previous literatures concerned with the mixed convection flow and heat transfer in a vertical channel with imposed heat flux or temperature at the channel wall, neglected wall thermal conduction [1,2]. The flow reversal, which occurs when the buoyancy parameter exceed a threshold value, was one of the most extensively investigated subjects, as it determined the flow structure in the channel and, consequently, heat transfer, pressure drop, fluid friction and entropy generation, etc. The regime of such buoyancy induced reversed flow has been presented comprehensively by researchers for fully developed or developing flow, accounting for cases of symmetric or asymmetric heated

* Corresponding author. Tel./fax: þ86 21 34206776. E-mail addresses: [email protected] (G. Yang), [email protected] (J.Y. Wu). http://dx.doi.org/10.1016/j.ijthermalsci.2015.07.023 1290-0729/© 2015 Elsevier Masson SAS. All rights reserved.

boundary conditions [3e10]. In a recent study, Desrayaud and Lauriat [11] numerically investigated the flow reversal phenomenon for laminar mixed convection of air in a vertical parallel-plate channel. The effects of the assisted buoyancy on the flow pattern and temperature profiles were discussed in detail, and the regime clet number of reversed flow was identified for high values of the Pe in a Pe-Gr/Re map. Their study was thereafter extended to three dimensional mixed convection flow by the present authors [12,13], who investigated the flow reversal and heat transfer in a three dimensional symmetrically heated rectangular duct. However, in most practical situation, such as for hot/cold fluid transport and heat exchangers, the boundary conditions of the fluid zone are not known initially but depend on the coupling between convection and wall conduction at the fluidesolid interface, and the effect of wall conduction is even pronounced in the thermal entrance region [14]. The earlier studies of conjugate heat transfer problems were mainly concerned with the coupling of wall conduction and pure forced convection flow or the coupling of wall conduction with natural convection. For laminar convection flows in parallel plates or in circular pipes, the effect of axial wall conduction was examined by many authors such as Davis and Gill [15], Mori et al. [16,17], Faghri and Sparrow [18], etc. The wall conduction was considered as one dimensional in these studies. Bilir [19] numerically analyzed the conjugate heat transfer problem within thermal developing laminar pipe flow, involving two dimensional (axial and radial) wall and axial fluid conduction. Adelaja et al. [20]

246

G. Yang, J.Y. Wu / International Journal of Thermal Sciences 98 (2015) 245e254

Nomenclature cp f g Gr H k l L N Nu p P Pe Pr Re Ri T U, V X, Y

specific heat of the fluid, J/(kg K) friction factor gravitational acceleration, m/s2 Grashof number height of the duct, m thermal conductivity, W/(m K) thickness of the solid zone, m length of the fluid zone, m grid number Nusselt number pressure, Pa non-dimensional pressure clet number Pe Prantdl number Reynolds number Richardson number temperature, K dimensionless velocity dimensionless coordinate system

analytically studied the effect of two dimensional wall and fluid clet number laminar flow heat transfer with conduction for low Pe convective boundary conditions of the third kind for the thermal entrance region problem. They found that an increase in the wall thickness resulted in reduced heat flux, while increases in the Biot number and the ratio of the wall-to-fluid thermal conductivity resulted in increased heat flux. For studies considering the coupling of natural convection, forced convection and wall conduction (i.e. conjugate mixed convection), Chou and Lien [21] numerically studied the effect of wall conduction on laminar mixed convection in the thermal entrance region of horizontal rectangular channels. Results showed that the flow and heat transfer characteristics are affected significantly by the wall thermal conduction and the buoyancy-induced secondary flow. Bernier and Baliga [22] investigated the conjugate conduction and laminar mixed convection in vertical pipes for upward flow and uniform heat flux. Their results were presented for the conditions of Pr ¼ 5, Gr ¼ 5000, Re ¼ 1 and 10, and different values of solid-tofluid thermal conductivity ratios and wall thickness-to-pipe diameter ratios. Laplante and Bernier [23] presented a numerical study aimed at quantifying the effects of wall conduction on laminar mixed convection in vertical pipes for a downward flow and a uniform wall heat flux boundary condition. Results indicated that a significant amount of heat is redistributed upstream of the heated section when the solid-to-fluid thermal conductivity and/or the wall thickness-to-pipe diameter ratios are high. Omara et al. [24,25] numerically studied transient heat transfer for downward and upward mixed convection in a vertical thick pipe partially exposed to a constant heat flux, focusing on the transient distribution of heat flux and friction coefficient. However, most of these studies considered the boundary conditions of the solid wall surface as constant heat flux, and the buoyancy induced flow reversal regime in the vertical channel with thick walls of constant temperature was not clearly revealed. The aim of the present study is to assess the combined effects of the wall conduction and assisted buoyancy on the flow and heat transfer characteristics in symmetrically heated vertical channels with isothermal thick walls, for the solid-to-fluid thermal conductivity ratios of 1  k* < ∞, wall thickness-to-channel length ratios of 0  l*  5, Reynolds numbers of 200  Re  1000, a Prandtl

Greek symbols b coefficient of volumetric expansion, 1/K q dimensionless temperature m dynamic viscosity of the fluid, Pa s n kinematic viscosity of fluid, m2/s r density of the fluid, kg/m3 Subscripts avg average value b bulk value cr critical value ∞ free-stream or inlet condition f fluid s solid w wall x, y Cartesian coordinates Superscript 0 dimensional variable * relative value

number of 0.7 and for various Grashof numbers. The aim is achieved by numerically solving the governing equations in the fluid and solid regions. Another important objective of this study is to highlight the flow reversal mechanism in the vertical channel with conjugate mixed convection heat transfer.

2. Physical description and mathematical formulation The physical model of mixed convection flow in a vertical channel with thermal conducting walls, along with the relevant dimensions considered in the present study is illustrated in Fig. 1. 0 and Newtonian fluids enter the channel with uniform velocity V∞ temperature T∞ from the bottom. The channel walls are of the thickness l (0  l* ¼ l/L  5) with the outer surface maintained at a

Fig. 1. Schematics of the upward flow in a parallel-plate channel with heated thick walls.

G. Yang, J.Y. Wu / International Journal of Thermal Sciences 98 (2015) 245e254

constant and uniform temperature Tw, which is higher than the temperature of the inlet fluid. The fluid zone in the channel is of length L and height H ¼ 50L, which has been proved to be large enough to have no influence on the occurrence of the reversal flow [11]. The physical properties of the fluid are taken to be constant except for the density in the body forces, for which the Boussinesq approximation has been adopted. The flow is assumed to be two dimensional, laminar and steady, and the effect of viscous dissipation is neglected. In these conditions, the continuity, momentum and energy equations for fluid can be written as follows:

vU vV þ ¼0 vX vY U

(1)

vU vU vP 1 þV ¼  þ V2 U vX vY vX Re

(2)

U

vV vV vP 1 Gr þV ¼  þ V2 V þ 2 q vX vY vX Re Re

(3)

U

vq vq 1 2 þV ¼ V q vX vY RePr

(4)

247

U ¼ V ¼ 0;

(11)

Solid-fluid interface (X ¼ 1/2, 0  Y  H/L):

 qs ¼ qf and

vq vX



¼ k* f



vq vX

 (12) s

where k* ¼ ks/kf is the thermal conductivity ratio of solid to fluid. In some numerical cases, we also used the simplified one dimensional wall conduction model for the solid walls, in order to obtain the results without considering the wall condition in the axial direction (Y). The integrated one dimensional form of Eq. (5) and the boundary condition of Eq. (11) then can be combined to the non-dimensional equation of:



vq vX

 ¼ f

k* ð1  qÞ l*

(13)

The heat transfer between the wall and fluid is calculated by the Nusselt number on the wall surface. The local bulk Nusselt number and ambient Nusselt number, which are based on the wall-to-bulk

and for solid is

V2 q ¼ 0

(5)

where the non-dimensional variables are respectively defined as:

U¼ ¼

U0 ; 0 V∞



V0 ; 0 V∞



X0 ; L



Y0 ; L



p 02 rV∞

and

q

T  T∞ Tw  T∞

The non-dimensional parameters applied in the above equations are respectively defined as

Re ¼

0 L rV∞ ; m

Pr ¼

mcp k

and

Gr ¼

gbðTw  T∞ ÞL3 v2

Due to the axial symmetry, only a half of the physical model is considered as the calculation domain (0  x  L/2 þ l and 0  y  H in Fig. 1). The dimensionless physical boundary conditions of interest in this study maybe written as follows: Inlet boundary (0  X  1/2, Y ¼ 0):

U ¼ 0; V ¼ 1 and q ¼ 0

(6)

Right boundary (X ¼ 1/2 þ l*, 0  y  H/L):

q¼1

(7)

Top/bottom solid boundary (1/2  X  1/2 þ l*, Y ¼ 0 or Y ¼ H/L):

vq=vY ¼ 0

(8)

Symmetry line (X ¼ 0, 0  Y  H/L):

U ¼ 0; vV=vX ¼ 0; and vq=vX ¼ 0

(9)

Outlet boundary (0  X  1/2, Y ¼ H/L):

vU=vY ¼ vV=vY ¼ vq=vY ¼ 0 Solid regions (1/2  X  1/2 þ l*, 0  Y  H/L):

(10) Fig. 2. Comparison of present (a) centerline axial velocity and (b) local ambient Nusselt number with benchmark values.

248

G. Yang, J.Y. Wu / International Journal of Thermal Sciences 98 (2015) 245e254

Table 1 Grid independence study for l* ¼ 0.5, k* ¼ 1, Re ¼ 300, Gr ¼ 3.22  105 and Pr ¼ 0.7. Grid size (Nx  Ny) 20 35 55 75

   

200 300 400 500

q

Vy

Nuavg

f

0.011745 0.011505 0.011221 0.011097

0.782863 0.783175 0.783898 0.786963

1.083776 1.081918 1.081027 1.081033

0.074126 0.075288 0.076019 0.076101

3.1. Solution methodology

and wall-to-ambient temperature differences, respectively, are defined as:

Nub ¼

 1 vq  1  qb vnw

 vq  Nu∞ ¼   vn w

(14)

(15)

where the qb is the dimensionless bulk temperature, which is calculated by

Z

1

qb ðYÞ ¼ ZX¼0 1

rVqdX

X¼0

3. Numerical details

(16) rVdX

Physically, as Nub is the ratio of the convective heat flux to the wall-to-bulk temperature difference that drives the heat transfer, it reflects the convective heat transfer coefficient, or the heat transfer efficiency under a certain temperature difference. On the other hand, Nu∞ is the ratio of the convective heat flux to the wall-toambient temperature difference, and it directly reflects the actual heat flux.

The governing equations along with the corresponding boundary conditions for the incompressible fluid flow is solved by using a finite volume based method according to SIMPLEC algorithm [26]. A third order accurate QUICK scheme [27] has been used to discretize the convective terms, and the diffusive terms are discretized using the central difference scheme for both the momentum and energy equations. The domain decomposition method is used to treat the fluid and solid zones separately with patched grids [28], and the equations are solved in the fluid and the solid regions simultaneously with an interface-value updating algorithm [29,30]. Since the grids in solid and fluid domains are non-overlapping patched and share the common collocation points at the interface, the continuity of the temperature is naturally satisfied. On the other hand, the heat flux continuity at the interface is fulfilled by solving tri-diagonal systems of equations, whose coefficient matrix is derived based on Eq. (12), with the Thomas method [31]. The updating of the interface value is carried out at each iterative step, to maintain Eq. (12) [32]. The grid distribution employed is non-uniform in Y direction, providing a finer mesh at the channel inlet. The size ratio between the largest cell and the smallest cell is 10 in the vertical direction, with a hyperbolic tangent function stretching between them [33]. In X direction, the grid distribution is uniform in the solid zone. While in the fluid zones, the grid is denser near the solid-fluid interface, due to the sensitivity of the result to the interface condition; and is also denser near the centerline of the channel, in order to capture the flow reversal phenomenon which occurs at the centerline. The solution is considered to converge when the relative variation of any dependent variable is less than 106.

Fig. 3. Contour lines of axial velocity (left) in the fluid zone and isotherms (right) in the channel with thick walls for (a) k* ¼ 2, (b) k* ¼ 10, (c) k* ¼ 50 and (d) k* ¼ 100 at l* ¼ 1, Re ¼ 300, Gr ¼ 1.98  105 (Gr/Re ¼ 660 or Ri ¼ 2.2) and Pr ¼ 0.7.

G. Yang, J.Y. Wu / International Journal of Thermal Sciences 98 (2015) 245e254

3.2. Code validation To validate the accuracy of the present numerical model, two direct comparisons are made between the present data with those available in the literature. Firstly, the streamwise variation of the centerline axial velocity in a symmetrically heated channel with constant wall temperature, along with the results by Ingham et al. [5] and Desrayaud and Lauriat [11] are presented in Fig. 2(a). Comparisons are carried out at Re ¼ 250 and H ¼ 50L for the buoyancy parameters of Gr/Re ¼ 1200 and 636 (all the parameters in the references are transformed to be characterized by the width of the channel, L). It can be seen that the present results agree well with those reported data, and the overall deviations are less than 6% and 4% on average, with respect to that by Ingham et al. [5] and Desrayaud and Lauriat [11], respectively. Furthermore, it is clear that the flow reversal occurs at the centerline of the channel with a large buoyancy force, and the minimum axial velocity decreases with the increase of the buoyancy force. Secondly, to verify the present scheme and code solving heat transfer problem with coupled solid thermal conduction and fluid convection, the natural convection flow and heat transfer in a square enclosure with a centered internal conducting square, which is reported by Das and Reddy [34], is recalculated here. The comparison of local Nusselt number at the hot wall for Ra ¼ 105 and k* ¼ 5 with the body size of L/W ¼ 0.5 and inclination angel of 90 is presented in Fig. 2(b). A good agreement can be seen with an average deviation of about 1.5%. This further validates the present numerical solution procedure.

249

with thick walls for various thermal conductivity ratios of k* ¼ 2, 10, 50 and 100 at fixed values of l* ¼ 1, Re ¼ 300, Gr ¼ 1.98  105 (Gr/ Re ¼ 660 or Ri ¼ 2.2) and Pr ¼ 0.7. Such flow parameters are chosen for the reason that the effect of the assisted buoyancy is large enough to induce the flow reversal (negative axial velocity) phenomenon if the wall thickness is neglected [11]. However, considering the thermal conduction of the channel walls, Fig. 3 clearly shows that both the accelerated and the decelerated zones shrink with the decrease of k*, and it is indicated from the contour lines of axial velocity in Fig. 3 that the assisted buoyancy induced flow reversal occurs at k*  50, and disappears at k*  10 in the channel of l* ¼ 1 for Re ¼ 300 and Gr ¼ 1.98  105. This phenomenon is due to the increase of the wall thermal resistance with the decrease of k*, which decreases the temperature at the fluidesolid interface, as can be inferred from the isothermal lines in Fig. 3. In order to further obtain the relations between the thermal conductivity ratio and the assisted buoyancy induced reversed flow regime, the critical buoyancy for the occurrence of flow reversal under various thermal conductivity ratios of k*  1 is investigated. The regime of reversed flow occurrence as a function of k* for the wall thickness of l* ¼ 1 is illustrated in Fig. 4(a). It is seen that with increasing k* from 1 to 100, the critical value of Gr/Re decreases

3.3. Grid independence study A grid independence test was conducted to decide the appropriate mesh resolution. The sensitivity studies for the influence of grid parameters were carried taking the case of l* ¼ 1/2 as an example. The variation of local temperature and axial velocity at the point of (X ¼ 0.5, Y ¼ 6), and the average Nusselt number and average friction factor at the fluidesolid interface with grid sizes for k* ¼ 1, Re ¼ 300, Gr ¼ 3.22  105 and Pr ¼ 0.7 are presented in Table 1. It shows that all the parameters have a monotonic increase or decrease in their values and approach a limiting value as the grid number increases. An examination of the parameters indicates that the grid size of Nx  Ny ¼ 55  400, is sufficiently fine to obtain results that are essentially grid independent. The local temperature, the axial velocity, the average Nusselt number and the average friction factor have a relative deviation of 1.11%, 0.39%, 0.005% and 0.11%, respectively, compared to the finest grid size. Therefore, taking into account the computational time and the accuracy of the results, the grid size of 55  400 is retained for l* ¼ 0.5. Similarly, the grid independence studies for the cases with other wall thicknesses were carried out, and the grid sizes were adopted as Nx ¼ (l*  60 þ 25), Ny ¼ 400 in the computation. 4. Results and discussion 4.1. Effect of solid-fluid thermal conductivity ratio As is well known, for mixed convection flows in a vertical channel, the fluids near the walls are accelerated by the assisted buoyancy force, while those near the centerline of the channel are decelerated due to the conservation of the mass flow rate. If the buoyancy effect is large enough, the centerline axial velocity even becomes negative, which is regarded as “flow reversal”. Such phenomenon is also exemplified in Fig. 3(a)e(d), which presents the contours of axial velocity (Vy ¼ 0, 0.5, 1, 1.5 and 2) in the fluid zones and isotherms (q ¼ 0.1, 0.3, 0.5, 0.7 and 0.9) in the channel

Fig. 4. The critical buoyancy parameter for the occurrence of the flow reversal for (a) various k* at l* ¼ 1 and (b) various 1/k* at 0.5  l*  2, Re ¼ 300 and Pr ¼ 0.7.

250

G. Yang, J.Y. Wu / International Journal of Thermal Sciences 98 (2015) 245e254

monotonously from 3540 to 625, which indicates an easier occurrence of the flow reversal. It should be mentioned that when the k* increases to infinitely large, which corresponds to the case where there is no wall resistance, the critical Gr/Re is found to be 596, which is very close to the 600 reported by Desrayaud and Lauriat [11]. It can also be inferred that if k* decreases to 0, the (Gr/Re)cr will becomes infinite large, due to the infinite wall thermal resistance. As such variation trend is similar to inverse proportional function, the variation of (Gr/Re)cr with thermal conductivity ratios is plotted with 1/k* in Fig. 4(b), together with those for other wall thickness ratios of l* ¼ 0.5, 1.5 and 2. Fig. 4(b) clearly shows that, for a constant wall thickness ratio, the (Gr/Re)cr increases linearly with 1/k*. Such linearly variation of (Gr/Re)cr with 1/k* has never been reported before. 4.2. Effect of wall thickness Fig. 5(a)e(d) presents the contour lines of axial velocity (Vy ¼ 0, 0.5, 1, 1.5 and 2) in the fluid zones and isotherms (q ¼ 0.1, 0.3, 0.5, 0.7 and 0.9) in the channel with thick walls for various length ratios of l* ¼ 0.5, 1, 1.5 and 2 at fixed values of k* ¼ 50, Re ¼ 300, Gr ¼ 1.98  105 (Gr/Re ¼ 660 or Ri ¼ 2.2) and Pr ¼ 0.7. It is seen from the distribution of axial velocity that, for each case, the minimum axial velocity is located at the centerline of the channel and the maximum axial velocity is located near the channel walls, due to the effect of assisted buoyancy force. The increase of l* is found to restrain the flow reversal phenomenon, as the flow reversal occurs at l*  1, and disappears at l*  1.5 for k* ¼ 50, Re ¼ 300 and Gr ¼ 1.98  105. The temperature field in Fig. 5 shows that though the thermal resistance increases with the increasing wall thickness, the shape of isotherms in the solid zone is almost the same for a certain thermal conductivity ratio. The variation of (Gr/Re)cr with l* for constant thermal conductivity ratios of k* ¼ 1, 2, 5, 10, 20 and 50, at Re ¼ 300, Gr ¼ 1.98  105

and Pr ¼ 0.7 is presented in Fig. 6. It's seen that, for a constant thermal conductivity ratio, the (Gr/Re)cr increases almost linearly with l* from (Gr/Re)cr ¼ 596 at l* ¼ 0. Based on the results in Figs. 4 and 6, it can be inferred that the (Gr/Re)cr in conjugate mixed convection can be obtained by linearly fitting the data of l* and 1/k*. As the (Gr/Re)cr is constant at the minimum value for the case of either 1/k* ¼ 0 or l* ¼ 0, the fitted correlation for (Gr/Re)cr with 1/k* and l* is presented as the following form:

 ðGr=ReÞcr ¼ 2934:40l* k* þ 595:13

(17)

The Eq. (17) is fitted by 49 sets of data from 0  l*  5 and 0 < 1/ k*  1 at Re ¼ 300 and Pr ¼ 0.7, and the R-squared is 0.999. As the l*/ k* is the lumped wall thermal resistance in radial direction (X), see Eq. (13), Eq. (17) indicates that the wall conduction in the axial direction (Y) has a negligible effect on the occurrence of the flow reversal in the fluid zone. In engineering design or calculation procedure, ignoring the wall axial conduction (i.e. consider one dimensional wall conduction) can greatly reduce the computational cost, as the heat conduction in the solid wall does not need to be solved and just with Eq. (13) imposed as the boundary condition [35]. In order to further investigate the relative effect of the wall axial conduction, Fig. 7 presents the distribution of temperature at the fluidesolid interface along the channel height considering 1-D wall conduction as well as the results for a complete 2-D simulation for l*/k* ¼ 1/40, l* ¼ 1, 2 and 5 and l*/k* ¼ 1/5, l* ¼ 1 and 2 at Gr ¼ 1.5  105, Re ¼ 300 and Pr ¼ 0.7. The temperature at the fluidesolid interface is selected as the parameter for comparison, as it is directly affected by the wall thermal conduction, which thereafter influence the flow and heat transfer in the fluid zone. It can be seen from the curves under the same l*/k*, that the effect of the wall axial conduction on the overall temperature distribution at the fluidesolid clet number interface is not obvious, which is due to the large Pe (Pe ¼ Re  Pr ¼ 210) selected in this figure [19]. However, its effect

Fig. 5. Contour lines of axial velocity (left) in the fluid zone and isotherms (right) in the channel with thick walls for (a) l* ¼ 0.5, (b) l* ¼ 1, (c) l* ¼ 1.5 and (d) l* ¼ 2 at k* ¼ 50, Re ¼ 300, Gr ¼ 1.98  105 (Gr/Re ¼ 660 or Ri ¼ 2.2) and Pr ¼ 0.7.

G. Yang, J.Y. Wu / International Journal of Thermal Sciences 98 (2015) 245e254

Fig. 6. The critical buoyancy parameter for the occurrence of the flow reversal for various l* at 1  k*  50, Re ¼ 300 and Pr ¼ 0.7.

251

the assisted buoyancy (Fig. 4), which increases the heat transfer efficiency. Furthermore, at a sufficiently large buoyancy parameter, in which case the reversed flow occurs (Gr/Re  1000 in Fig. 8(a)), the Nub shows a local maximum at the location corresponds to the reversed flow regions as the flow is mostly intensified near the wall. For conjugate mixed convection flow, Fig. 8(b) presents the axial development of bulk Nusselt number at the fluidesolid interface for various thermal conductivity ratios of k*  2, Gr/Re ¼ 660, Re ¼ 300 and Pr ¼ 0.7. It is seen that, increasing k*/l* can always increase the Nub. For the case of l* ¼ 0 (k*/l* ¼ ∞), the Nub decreases to the asymptote of Nub ¼ 3.77 along the channel length, which is the same as that by Desrayaud and Lauriat [11]. Similar to the results in Fig. 8(a), for a thermal conductivity ratio of k*  50, a local maximum of Nub occurs at the flow reversal zone. Fig. 9 presents the axial development of ambient Nusselt number, based on the wall-to-ambient temperature difference, for various thermal conductivity ratios of k*/l*  2 at Gr/Re ¼ 660 and 0, Re ¼ 300 and Pr ¼ 0.7. For a constant k*/l*, Fig. 9 shows the Nu∞ decreases monotonously along channel height, and the variation becomes gentler as k*/l* decreases. So, the Nu∞ is higher for a larger k*/l* near the inlet, and becomes larger for a smaller k*/l* from some

on the temperature distribution near the channel inlet is prominent. For example, the relative deviation between the 1-D and 2-D simulation results for l* ¼ 5 at Y ¼ 0.1 for k* ¼ 40l*, Gr ¼ 1.5  105, Re ¼ 300 and Pr ¼ 0.7 is about 32%. Such temperature deviation decreases sharply towards downstream, which becomes not obvious when Y  4. 4.3. Heat transfer Fig. 8(a) presents the axial development of local bulk Nusselt number at the fluidesolid interface of the channel with l* ¼ 1 and k* ¼ 10 for various buoyancy parameters of 0  Gr/Re  2000 at Re ¼ 300 and Pr ¼ 0.7. For Gr/Re ¼ 0 (forced convection), the Nub decreases monotonously along the channel height until it reaches to the minimum value of Nub ¼ 2.75, when the flow becomes thermal fully developed. For mixed convection, the bulk Nusselt number is found to be increased by the assisted buoyancy force. The reason is that the fluid flow near the channel wall is accelerated by

Fig. 7. Streamwise variation of the temperature at the solid-fluid interface for k*/ l* ¼ 40, l* ¼ 1, 2 and 5 and k*/l* ¼ 5, l* ¼ 1 and 2 at Gr ¼ 1.5  105, Re ¼ 300 and Pr ¼ 0.7.

Fig. 8. Streamwise variation of the local bulk Nusselt number at the fluidesolid interface for (a) various Gr/Re at l* ¼ 1, k* ¼ 10 and (b) various k*/l* at Gr/Re ¼ 660, Re ¼ 300 and Pr ¼ 0.7.

252

G. Yang, J.Y. Wu / International Journal of Thermal Sciences 98 (2015) 245e254

complete 2-D simulation in Fig. 9, where a slight deviation can be found in the entrance region, which is due to the effect of wall conduction in the axial direction. 4.4. Effect of Reynolds number

Fig. 9. Streamwise variation of the local ambient Nusselt number at the fluidesolid interface for various k*/l* at Gr/Re ¼ 0 and 660, Re ¼ 300 and Pr ¼ 0.7.

downstream regions. Similarly, compared to that in forced convection case (curve 40 ), when the buoyancy assisted the flow (curve 4), higher Nub is attained near the inlet and lower Nub is attained at the downstream regions. Such phenomenon is caused by the increasing convective heat transfer rates that results from the accelerating flow near the walls at the inlet. This, in turn, diminishes the magnitude of the wall-to-fluid temperature difference and causes the lowering of local mixed convection Nusselt numbers in the downstream regions. The values of Nu∞ considering 1-D wall conduction at k*/l* ¼ 10 is also compared with the results of

Fig. 10 presents the effect of Reynolds number (Re ¼ 200, 300, 500 and 700) on the flow and temperature patterns in the channel for l* ¼ 1, k* ¼ 20, Gr ¼ 2.25  105 and Pr ¼ 0.7. It clearly shows that the increased Reynolds number increases the relative effect of inertial force to buoyancy force, which makes the axial velocity in the fluid zone more uniform. Therefore, the size of the flow reversal cell diminishes with the increasing Reynolds number, until it disappears at Re  700. The temperature in the solid zone also shows to decrease with the increasing Reynolds number, due to the increase of heat transfer between the fluids and the solid wall. Then, the critical Gr/Re for the occurrence of the flow reversal for 200  Re  1000, in the channel of l* ¼ 1 and k* ¼ 20 at Pr ¼ 0.7 is presented in Fig. 11, in order to obtain the effect of Reynolds number on the flow reversal in the vertical channel with solid walls. The solid dots are for flows with flow reversal while the hollow ones are for the flow without flow reversal. The dashed line in Fig. 11, indicating the (Gr/Re)cr, is independent of the Reynolds number. Such effect of Re has also been checked for 0  l*  5 and 0 < 1/k*  1, and the same result is found. So the Eq. (17), which is obtained from Re ¼ 300, can be further extended to 200  Re  1000. Fig. 12 presents the streamwise development of ambient Nusselt number at the fluidesolid interface for Re ¼ 200, 300, 500 and 700, Gr ¼ 2.25  105 and Pr ¼ 0.7, in the channel of l* ¼ 1 and k* ¼ 20. Increasing Reynolds is found to increase the heat transfer at the channel walls, due to the increase of fluid velocity. The enhancement of heat transfer caused by Reynolds number is even obvious in the region of Y  5, where the effect of the assisted buoyancy on

Fig. 10. Contour lines of axial velocity (left) in the fluid zone and isotherms (right) in the channel with thick walls for (a) Re ¼ 200, (b) Re ¼ 300, (c) Re ¼ 500 and (d) Re ¼ 700 at l* ¼ 1, k* ¼ 20, Gr ¼ 2.25  105 and Pr ¼ 0.7.

G. Yang, J.Y. Wu / International Journal of Thermal Sciences 98 (2015) 245e254

253

walls, the assisted buoyancy induced flow reversal tends to occur with the increasing k* or decreasing l*, and the critical buoyancy parameter for the occurrence of the flow reversal increases linearly with the increasing l*/k*, and it can be obtained by (Gr/ Re)cr ¼ 2934.40l*/k* þ 595.13. The effect of wall axial conduction on the temperature distribution near the channel inlet is notable, but its effects on the occurrence of the flow reversal and overall heat transfer are negligible for the ranges of parameters concerned. The bulk Nusselt number at the solid-fluid interface increases with the increasing assisted buoyancy force and the increasing k*/ l*, especially at the locations where the flow reversal occurs. The ambient Nusselt number decreases monotonously along the channel height, and the variation becomes gentler as the k*/l* decreases. The occurrence of the assisted buoyancy induced flow reversal in the vertical channel with thick walls is independent of the Reynolds number for 200  Re  1000. Increasing Reynolds number increases the heat transfer between the fluids and the wall, which is even more prominent in the regions of Y  5. Fig. 11. Regime of reversed flow occurrence for 200  Re  1000, in the channel of l* ¼ 1 and k* ¼ 20 at Pr ¼ 0.7.

the heat transfer is less strong. It's also seen from Fig. 12 that the simplify of 1-D wall conduction overestimates the local Nusselt number near the inlet region of Y  4, while the average Nusselt number along the channel has a relative deviation of 0.83% and 0.72% for Re ¼ 200 and 700, respectively. 5. Conclusions A detailed numerical analysis of buoyancy assisted mixed convection in the entrance region of a symmetrically heated vertical channel with thermal conducting walls is performed in the present study. Results are obtained for the solid-to-fluid thermal conductivity ratios of 1  k* < ∞, the wall thickness-to-channel length ratios of 0  l*  5, the Reynolds numbers of 200  Re  1000, the Prandtl number of 0.7 and for various Grashof numbers. This work confirms the effect of assisted buoyancy force on the flow reversal phenomenon, which occurs at the centerline of the vertical channel. For mixed convection flow in a channel with thick

Fig. 12. Streamwise development of ambient Nusselt number at the fluidesolid interface for various Reynolds numbers at l* ¼ 1, k* ¼ 20, Gr ¼ 2.25  105 and Pr ¼ 0.7.

Acknowledgments This work was supported by the National Natural Science Foundation of China (NSFC) under the contract No. 51476096.

References [1] B. Metais, E.R.G. Eckert, Forced, mixed, and free convection regimes, J. Heat. Transf. 86 (1964) 295e296, http://dx.doi.org/10.1115/1.3687128. [2] J.D. Jackson, M.A. Cotton, B.P. Axcell, Studies of mixed convection in vertical tubes, Int. J. Heat. Fluid Flow. 10 (1989) 2e15, http://dx.doi.org/10.1016/0142727X(89)90049-0. [3] B. Zeldin, F.W. Schmidt, Developing flow with combined forced-free convection in an isothermal vertical tube, J. Heat. Transf. 94 (1972) 211e223. [4] W. Aung, G. Worku, Theory of fully developed, combined convection including flow reversal, J. Heat. Transf. 108 (1986) 485e488. http://www.scopus.com/ inward/record.url?eid¼2-s2.00022711595&partnerID¼40&md5¼c0254f5adbdf37de63c4a822858cafe8. [5] D. Ingham, D. Keen, P. Heggs, Two dimensional combined convection in vertical parallel plate ducts, including situations of flow reversal, Int. J. Numer. Methods Fluids 26 (1988) 1645e1664. http://onlinelibrary.wiley.com/doi/10. 1002/nme.1620260713/abstract (accessed 16.05.13). [6] M.A. Bernier, B.R. Baliga, Visualization of upward mixed-convection flows in vertical pipes using a thin semitransparent gold-film heater and dye injection, Int. J. Heat. Fluid Flow. 13 (1992) 241e249, http://dx.doi.org/10.1016/0142727X(92)90037-A. [7] A. Moutsoglou, Y.D. Kwon, Laminar mixed convection flow in a vertical tube, J. Thermophys. Heat. Transf. 7 (1993) 361e368. http://arc.aiaa.org/doi/pdf/10. 2514/3.428. [8] M. Wang, T. Tsuji, Y. Nagano, Mixed convection with flow reversal in the thermal entrance region of horizontal and vertical pipes, Int. J. Heat. Mass Transf. 37 (1994) 2305e2319, http://dx.doi.org/10.1016/0017-9310(94) 90372-7. [9] C.H. Cheng, C.J. Weng, W. Aung, Buoyancy-assisted flow reversal and convective heat transfer in entrance region of a vertical rectangular duct, Int. J. Heat. Fluid Flow. 21 (2000) 403e411. http://www.sciencedirect.com/science/ article/pii/S0142727X00000059 (accessed 15.04.13). [10] M. Zghal, N. Galanis, C. Tam, Developing mixed convection with aiding buoyancy in vertical tubes: a numerical investigation of different flow regimes, Int. J. Therm. Sci. 40 (2001) 816e824. [11] G. Desrayaud, G. Lauriat, Flow reversal of laminar mixed convection in the entry region of symmetrically heated, vertical plate channels, Int. J. Therm. Sci. 48 (2009) 2036e2045, http://dx.doi.org/10.1016/j.ijthermalsci.2009.03.002. [12] G. Yang, J.Y. Wu, L. Yan, Flow reversal and entropy generation due to buoyancy assisted mixed convection in the entrance region of a three dimensional vertical rectangular duct, Int. J. Heat. Mass Transf. 67 (2013) 741e751. [13] G. Yang, J.Y. Wu, Effect of aspect ratio and assisted buoyancy on flow reversal for mixed convection with imposed flow rate in a vertical three dimensional rectangular duct, Int. J. Heat. Mass Transf. 77 (2014) 335e343, http:// dx.doi.org/10.1016/j.ijheatmasstransfer.2014.05.027. [14] R.K. Shah, A.L. London, Laminar Flow Forced Convection in Ducts, Academic Press, New York, 1978. [15] E.J. Davis, W.N. Gill, The effects of axial conduction in the wall on heat transfer with laminar flow, Int. J. Heat. Mass Transf. 13 (1970) 459e470, http:// dx.doi.org/10.1016/0017-9310(70)90143-2.

254

G. Yang, J.Y. Wu / International Journal of Thermal Sciences 98 (2015) 245e254

[16] S. Mori, M. Sakakibara, A. Tanimoto, Steady heat transfer to laminar flow in a circular tube with conduction in tube wall, Heat. Transf. Jpn. Res. 3 (1974) 37e46. [17] S. Mori, T. Shinke, M. Sakakibara, A. Tanimoto, Steady heat transfer to laminar flow between parallel plates with conduction in wall, Heat. Transf. Jpn. Res. 5 (1976) 17e25. [18] M. Faghri, E.M. Sparrow, Simultaneous Wall and fluid axial conduction in laminar pipe-flow heat transfer, ASME Trans. J. Heat. Transf. 102 (1980) 58e63, http://dx.doi.org/10.1115/1.3244249. [19] S. Bilir, Laminar flow heat transfer in pipes including two-dimensional wall and fluid axial conduction, Int. J. Heat. Mass Transf. 38 (1995) 1619e1625. http://www.sciencedirect.com/science/article/pii/0017931094002692 (accessed 05.03.14). [20] A.O. Adelaja, J. Dirker, J.P. Meyer, Effects of the thick walled pipes with convective boundaries on laminar flow heat transfer, Appl. Energy (2014) 1e8, http://dx.doi.org/10.1016/j.apenergy.2014.01.072. [21] F.C. Chou, W.Y. Lien, Effect of wall heat conduction on laminar mixed convection in the thermal entrance region of horizontal rectangular channels, €rme Und Stoffübertragung 26 (1991) 121e127, http://dx.doi.org/10.1007/ Wa BF01590110. [22] M.A. Bernier, B.R. Baliga, Conjugate conduction and laminar mixed convection in vertical pipes for upward flow and uniform wall heat flux, Numer. Heat. Transf. Part A Appl. 21 (1992) 313e332, http://dx.doi.org/10.1080/ 10407789208944879. [23] G. Laplante, M. Bernier, Opposing mixed convection in vertical pipes with Wall heat-conduction, Int. J. Heat. Mass Transf. 40 (1997) 3527e3536. http:// www.cheric.org/research/tech/periodicals/view.php?seq¼96174. [24] A. Omara, S. Abboudi, Numerical analysis of transient conjugated downward laminar mixed convection in a vertical pipe, Numer. Heat. Transf. Part A Appl. (2007) 225e247, http://dx.doi.org/10.1080/10407780600762681. [25] A. Omara, S. Abboudi, Transient conjugated upward and downward mixed convection in a partially heated thick pipe, Comput. Therm. Sci. 3 (2011)

[26] [27]

[28] [29]

[30]

[31]

[32] [33] [34]

[35]

1e14. http://www.dl.begellhouse.com/journals/ 648192910890cd0e,1e01f3fc0d41604f,48c726e809a445b4.html (accessed 16.04.14). J.P. Van Doormaal, D.G. Raithby, Enhancements of the SIMPLE method for predicting incompressible fluid flows, Numer. Heat. Transf. 7 (1984) 147e163. B.P. Leonard, A stable and accurate convective modelling procedure based on quadratic upstream interpolation, Comput. Methods Appl. Mech. Eng. 19 (1979) 59e98, http://dx.doi.org/10.1016/0045-7825(79)90034-3. A. Toselli, O. Widlund, Domain Decomposition Method e Algorithms and Theory, Springer-Verlag, Berlin, 2005. K.D. Carlson, W. Lin, C.-J. Chen, Pressure boundary conditions of incompressible flows with conjugate heat transfer on nonstaggered grids. Part II: applications, Numer. Heat. Transf. Part A Appl. 32 (1997) 481e501, http:// dx.doi.org/10.1080/10407789708913902. C.-H. Zhang, W. Zhang, G. Xi, A pseudospectral multidomain method for conjugate conduction-convection in enclosures, Numer. Heat. Transf. Part B Fundam. 57 (2010) 260e282, http://dx.doi.org/10.1080/10407790.2010. 489880. W.H. Press, S.A. Teukolsky, W.T. Vetterling, B.P. Flannery, Numerical Recipes: the Art of Scientific Computing, third ed., Cambridge University Press, New York, 2007. W.-Q. Tao, Numerical Heat Transfer, second ed., Xi’an Jiaotong University Press, Xi'an, China, 2001. J.F. Thompton, S. Warsi, C.W. Mastin, Numerical Grid Generation-Foundations and Applications, 1986. M.K. Das, K.S.K. Reddy, Conjugate natural convection heat transfer in an inclined square cavity containing a conducting block, Int. J. Heat. Mass Transf. 49 (2006) 4987e5000, http://dx.doi.org/10.1016/j.ijheatmasstransfer.2006. 05.041. D.A. Kaminski, C. Prakash, Conjugate natural convection in a square enclosure: effect of conduction in one of the vertical walls, Int. J. Heat. Mass Transf. 29 (1986) 1979e1988, http://dx.doi.org/10.1016/0017-9310(86)90017-7.