- Email: [email protected]
Lattice-Boltzmann computation of natural convection in a partitioned enclosure with inclined partitions attached to its hot wall
ARTICLE IN PRESS
Physica A 368 (2006) 481–494 www.elsevier.com/locate/physa
Lattice-Boltzmann computation of natural convection in a partitioned enclosure with inclined partitions attached to its hot wall Mohammed Jamia, Ahmed Mezrhaba,, M’hamed Bouzidib, Pierre Lallemandc a
De´partement de Physique, Faculte´ des Sciences, Universite´ Mohamed I, Oujda, Morocco b Universite´ Clermont-2, I.U.T. Ave A. Briand, 03100 Montluc- on, France c LIMSI, Baˆtiment 508, Universite´ Paris-Sud, 91405 Orsay, France Received 7 October 2005 Available online 13 January 2006
Abstract This study presents a numerical investigation of laminar natural convection heat transfer in an inclined enclosure, differentially heated, with inclined partitions attached to its hot wall. The top and bottom walls of the cavity are insulated. Numerical solutions are obtained by using the hybrid lattice-Boltzmann method finite-difference. Effects of partition length, partition inclination angle, partitions number and aspect ratio cavity are studied for Prandtl number 0.71 and cavity inclination angles of 45 and 90 . Rayleigh numbers considered are 5 105 for a square cavity and 104 for a cavity with an aspect ratio equal to 6. Results are reported in terms of isotherms, streamlines, local and average Nusselt numbers. r 2006 Elsevier B.V. All rights reserved.
1. Introduction Natural convection in a partitioned cavity has been extensively studied using numerical simulations and experiments because of its interest and importance in industrial applications. Some applications are solar collectors, fire research, electronic cooling, aeronautics, chemical apparatus and building engineering. Most of the papers in this field are substantially oriented toward the study of natural convection in enclosed square or rectangular cavities [1]. Among these works, we mention the one of Mezrhab and Bchir [2] who used the finitevolume method to study the effect of adding a thick partition located vertically close to the hot wall of a differentially heated cavity, forming a narrow vertical channel in which the flow is controlled by vents at the bottom and top of the partition. They have shown that radiation has a significant influence on the flow and heat transfer in the channel. Frederick [3] has studied numerically natural convection in an air-filled, differentially heated, inclined square cavity, with a single partition attached to its cold wall, at Rayleigh numbers of 103 –105 . He showed that the partition leads to the suppression of convection, and reduces the heat transfer by up to 47% in comparison to the empty cavity at the same Rayleigh number. He used Corresponding author. Tel.: +212 5650 0601; fax: +212 5650 0603.
E-mail addresses: [email protected] (M. Jami), [email protected], [email protected] (A. Mezrhab), [email protected] (M. Bouzidi), [email protected] (P. Lallemand). 0378-4371/$ - see front matter r 2006 Elsevier B.V. All rights reserved. doi:10.1016/j.physa.2005.12.029
ARTICLE IN PRESS M. Jami et al. / Physica A 368 (2006) 481–494
482
Nomenclature A g H k L Lp Lp Np Num Nu Pr Ra Tc Th T0 u v x; y X; Y
cavity aspect ratio, A ¼ H=L acceleration of gravity enclosure height thermal conductivity of the fluid enclosure width partition length dimensionless partition length, Lp ¼ Lp =L number of partitions RH average Nusselt number ð1=AÞ 0 qy=qx wall dy local Nusselt number, L qy=qx wall Prandtl number, Pr ¼ n=a Rayleigh number, gb DTL3 =na temperature of the cold wall temperature of the hot wall average fluid temperature ðT h þ T c Þ=2 x velocity component y velocity component Cartesian coordinates dimensionless Cartesian coordinates, X ¼ x=L, Y ¼ y=L
Greek symbols a thermal diffusivity of the fluid b volumetric expansion coefficient d partition thickness d dimensionless partition thickness, d ¼ d=L DT maximal difference temperature, T h T c DX ; DY dimensionless grid spacings in X and Y directions g partition inclination n kinematic viscosity of the fluid f inclination angle of the enclosure with respect to the horizontal c stream function C dimensionless stream function, c=a y dimensionless temperature ðT T 0 Þ=DT Subscripts max pþ p
maximal value right surface of the partition left surface of the partition
a finite-difference over-relaxation procedure for the solution of the mass, momentum and energy transfer governing equations. Let us mention that studies on natural convection in rectangular enclosures heated from below and cooled along a single side or both sides have been carried out respectively by Anderson and Lauriat [4] and by Ganzarolli and Milanez [5]. More recently, the case of heating from one side and cooling from the top has been analyzed by Aydin et al. [6] who investigated the effect of Prandtl number upon heat and momentum transfer inside square cavities. They also studied the influence of the aspect ratio for air-filled, rectangular enclosures [7]. Acharya and Jetli [8] have studied numerically the heat transfer and the fluid flow
ARTICLE IN PRESS M. Jami et al. / Physica A 368 (2006) 481–494
483
within a square cavity divided by a single obstacle, attached to the ceiling or the floor. They considered different positions and heights of the partition. They showed that the heat transfer is strongly influenced by the height of the partition; nevertheless, its position has a rather weak effect on the total heat transfer. Finally, we mention that in a similar study of natural convection in an inclined square enclosure with partitions attached to its cold wall, Mezrhab et al. [9] concluded that at lower Rayleigh number Ra (Rap105 ), the average Nusselt number Num is higher in inclined cavities than in vertical ones; while at larger Ra (Ra ¼ 106 ), the opposite occurs. In vertical cavities, Num decreases with increasing the dimensionless length of the partition Lp . However in inclined ones, this trend is found only at lower Ra, whereas at higher Ra, Num shows a maximum at about Lp ¼ 0:5. Num also decreases when the number of partitions ðN p Þ attached to the cold wall of the enclosure increases; however, beyond a certain value of ðN p Þ, which depends of Ra, Num becomes almost constant. The objective of this paper is to study numerically the heat transfer in an inclined square enclosure with inclined partitions attached to its hot wall. The numerical method proposed is based on the coupling between the lattice-Boltzmann method and finite-difference. 2. Numerical technique and mathematical formulation 2.1. Physical system Fig. 1 presents the geometry and the boundary conditions of a two-dimensional square enclosure placed in a gravity field. The laminar air flow can be described by assuming that the fluid is Newtonian and satisfies the Boussinesq approximation, it is initially stationary at a uniform temperature T 0 . We consider a cavity with two walls maintained respectively at T c and T h , making an angle f with the horizontal. The other two walls are adiabatic, midway between them, an inclined partition of length Lp and thickness d is fitted to the hot wall. Concerning applications of this study to systems such as solar collectors, the partition should be very thin and it is taken as a single line of mesh points. Thus, the dimensionless partition thickness d was fixed to one grid spacing. 2.2. Lattice-Boltzmann equation The simulations have been performed using a combination of lattice-Boltzmann equation (LBE) and finitedifference scheme for the temperature with suitable coupling between them to describe thermo-convective effects. LBE has been developed over the last 20 years as a simple means to simulate fluid flows [10]. It is
Fig. 1. Square cavity simulated.
ARTICLE IN PRESS M. Jami et al. / Physica A 368 (2006) 481–494
484
inspired by the kinetic theory of gases which considers a function f ðr; v; tÞ that describes the population of particles of the gas at location r, with velocity v at time t. In LBE space and time are discretized on simple uniform grids, and the possible particle velocities are limited to a small number chosen such that particles move exactly from one node of the grid to a neighbouring one during one time step. In the usual LBE models, the space grid is taken square or cubic, respectively for 2-d and 3-d problems. The number of velocities is chosen to include no move, move to the nearest neighbour and move to the next nearest neighbour. As we describe here results for 2-d situations, we take the simple D2Q9 model with 9 velocities (example: node located at the center of the Fig. 2). The variables are thus f i ðrj ; tl Þ;
i 2 f0; . . . ; 8g,
(1)
which are symbolically represented by the vector F, with rj a node of a square lattice of unit size DX , and tl ¼ l dt. All quantities are normalized in terms of the proper combination of DX and dt. The equation of motion for the quantities f i is given by f i ðrj þ ci ; t þ 1Þ f i ðrj ; tÞ ¼ OF ,
(2)
where the left-hand side represents the motion of particles (treated exactly with the particular choice of elementary velocities ci ) and the right-hand side describes particle redistribution at each lattice node due to collisions. The operator O is highly simplified compared to what would be required to take into account real collisions but it can be chosen in such a way that the large-scale behaviour of the model is similar to that of a real viscous fluid following the Navier–Stokes equations. In this work we use the technique proposed by d’Humie`res [11] (through the computation of moments) nowadays known as the MRT–LBE technique (multirelaxation time) as it is numerically more stable than the simple LBE–BGK technique and can be tuned to give more satisfactory behaviour at boundaries. Compared to a standard finite-difference code one may say that LBE uses more degrees of freedom (here 9 per node). The operator O is such that some combinations of f i are conserved (mass and linear momentum) whereas others follow some relaxation equation. It is the relaxation of these combinations that is responsible for the actual values of the speed of sound and of the transport coefficients. Apart from its simplicity, LBE allows to implement boundaries without much difficulty.
Fig. 2. Grid distributions of computational domain.
ARTICLE IN PRESS M. Jami et al. / Physica A 368 (2006) 481–494
485
The simplest boundary conditions so-called ‘‘bounce-back’’ allow to implement Dirichlet-type boundaries located to a good approximation half-way between the last ‘‘fluid node’’ and the first ‘‘solid node’’ [12,13]. To represent boundaries that cross a link between fluid and solid that is located at a relative distance q (see Fig. 3), one can use a variety of approximations. Here, we consider a simple combination of ‘‘bounce-back’’ together with quadratic interpolation for each particle link that crosses the interface [14,15]. The LBE part of our code yields values of the fluid velocity uðtÞ, vðtÞ at the nodes of the lattice. It is straightforward to add the effect of a buoyancy force through modification of the collision process [16]: instead of conserving linear momentum, we apply a change of the components of the linear momentum proportional to the buoyancy gbðT T 0 Þ. The dynamics of the temperature is computed using the advection–diffusion equation with finite-difference stencils to compute the space derivatives. Time evolution is done with a simple second-order predictor–corrector scheme. As we use the same thermal diffusivity for the fluid and the material used to make the partitions, there is no special treatment of the partitions (LBE automatically gives 0 for the fluid velocity inside the partition). In practice, we take the external walls of the cavity parallel to the axis of the LBE-grid and apply gravity at an angle chosen to satisfy the geometry of Fig. 1. This allows us to use simple ‘‘bounce-back’’ on the external boundary of the fluid, therefore, the fluid velocity is 0 at mid-distance between the last fluid node and the first solid node, that is along the axis Ox at x ¼ 12 and x ¼ N x þ 12 and at y ¼ 12 and y ¼ N y þ 12 along the axis Oy as shown in Fig. 2. On the interface between the fluid and the partition we have used a combination of the bounce-back scheme and interpolations because there are some particles that are going to be in interstitial positions after being rebounded. This technique is graphically depicted in Fig. 3. We have used the following interpolation formulas [14]: f ¯i ðxf ; tÞ ¼ qð1 þ 2qÞf i ðxf þ ci ; tÞ þ ð1 4q2 Þf i ðxf ; tÞ qð1 2qÞf i ðxf ci ; tÞ;
f ¯i ðxf ; tÞ ¼
qo1=2
1 ð2q 1Þ ð2q 1Þ f i ðxf þ ci ; tÞ þ f ¯i ðxf ci ; tÞ f ¯ ðxf 2ci ; tÞ; qð1 þ 2qÞ q ð2q þ 1Þ i
qX1=2
(3)
(4)
where f ¯i is the distribution function of the velocity c¯i ci . The fluid velocity on the enclosure walls and the partition is taken equal to 0.
/ / (a)
/ (b)
/ (c)
Fig. 3. Illustration of the boundary conditions for a rigid wall located between two grid sites in the direction of the motion of the particle that is rebounded.
ARTICLE IN PRESS M. Jami et al. / Physica A 368 (2006) 481–494
486
2.3. Equation of temperature The energy equation is solved by the finite-difference, in order to provide us with the field of temperature inside the cavity. The energy equation for the present system can be expressed by the following 2-d equation: 2 qy qy qy q y q2 y þu þv ¼a þ . (5) qt qx qy qx2 qy2 The coupling between the flow and the temperature field is done as follows: dimensionless temperature y provides a buoyancy force in the x and y directions because the square cavity is inclined by an angle f with respect to the horizontal. Therefore, f x ¼ gb DTyðr; tÞ sin f and f y ¼ gb DTyðr; tÞ cos f where yðr; tÞ is the local dimensionless temperature at the position of the air node r at time t. Let us note that the terms introduced to couple flow and temperature neglect compressibility. This procedure gives very good results in the simulation of the convective flows as that was shown by Mezrhab et al. [17]. As a matter of fact, many attempts have been made to simulate the thermal effects with LBM techniques [18]. However, this path presented difficulties as far as isotropy and stability [19,20] are concerned. Thus, it is more convenient to treat separately density and momentum (athermal fluid) and temperature. As the temperature grid nodes are the same as the fluid nodes, this leads to consider ‘‘ghost’’ nodes outside the fluid at x ¼ 0 and N x þ 1 and at y ¼ 0 and N y þ 1 with a temperature computed by extrapolating the information on the solid boundary (Dirichlet or Neuman) and the first and second point inside the fluid. The thermal boundary conditions at the walls are as follows: yð1=2; yÞ ¼ yc ¼ 0:5 yðN x þ 1=2; yÞ ¼ yh ¼ 0:5 qy ðx; 1=2Þ ¼ 0 qy qy ðx; N y þ 1=2Þ ¼ 0 qy
1 1 pypN y þ ; 2 2 1 1 for pypN y þ ; 2 2 1 1 for pxpN x þ ; 2 2 1 1 for pxpN x þ : 2 2 for
At the partition, the conditions of continuity of temperature and temperature gradient are imposed. These conditions are written as ypþ ¼ yp ;
qy qy ¼ , qzp qzpþ
(6)
where z is the normal direction to the partition. 3. Results and discussions The parameters governing the problem are A; Lp ; N p ; Pr; Ra; f and g. Except for the Prandtl number that was fixed at Pr ¼ 0:71, the other parameters were varied with the aim of studying their effect on the thermal transfer and the flow in the cavity. A preliminary study was carried out to determine the optimum uniform grid (i.e. the best compromise between accuracy and computational costs). In the case of a square cavity with a single partition attached to its hot wall, Table 1 presents the average Nusselt number
Table 1 Variation of Num and Cmax according to the mesh size N x N y Nx Ny
41 41
61 61
81 81
101 101
Num Cmax
5.759 24.90
5.898 25.18
5.917 25.24
5.917 25.25
ARTICLE IN PRESS M. Jami et al. / Physica A 368 (2006) 481–494
487
pffiffiffi and maximal dimensionless stream function as a function of the mesh size N x N y for A ¼ 1, Lp ¼ 2=4, Ra ¼ 5 105 , f ¼ 45 and g ¼ 45 . As can be seen, the solutions become grid independent at Nx Ny ¼ 81 81. Consequently, for the computations reported in the case of a square cavity, a grid with 81 81 points was chosen to optimize the relation between accuracy and computing time. The code was extensively exercised on benchmark problems to check its validity. In the cases of an empty cavity, a square cavity with a hot circle inside, and concentric cylinders, the results that we obtain with our simulation technique were checked for accuracy against the earlier published numerical and experimental results reported by different authors, and the agreement between the present and previous results was very good in Ref. [17]. The code was also validated for cavities divided by a thin partition as shown in Ref. [9]. It was concluded that the largest discrepancies between our and published results can be estimated to be less than 1%. For this reason, and for the sake of brevity it is not repeated here. Based on the above studies, it was concluded that the code could be reliably applied to the problem considered here. 3.1. Square cavity with a single partition 3.1.1. Isotherms and streamlines pffiffiffi Figs. 4–6 show the streamlines and isotherms for f ¼ 90 , Lp ¼ 2=4, g ¼ 90 , 45 and 135 , respectively. Let us note that near the partition and the cold wall, the flow is intense, the fluid rotates counter-clockwise owing to the positions of the hot and cold walls. Whatever the partition orientation, the isotherms and streamlines are considerably altered, in comparison to those corresponding to the empty cavity at the same Rayleigh number. In fact, when the partition is inclined by 45 , the major part of the fluid is heated by the lower part of the hot wall. Then it is reflected toward a region close to the cold wall, the rest of the fluid recirculates in the upper right corner of the enclosure, following the same direction of the main circulation. In this case the circulation is stronger than for the two others cases of g. For the partition inclination (g ¼ 135 ), the major part of the fluid circulates near the walls. The circulation and the heat transfer decrease and the point of maximum stream function is displaced towards the left part of the cavity. Independently of the value of g, from the temperature contours, it is apparent that thermal boundary layers form along the hot and cold walls. Far from the partition, the isotherms are nearly horizontal and resemble to those obtained in the case of an empty cavity, so the temperature stratification persists; however, around the partition the isotherms become more affected. We note that the streamlines and isotherms corresponding to g ¼ 45 and 135 are very different from those corresponding to g ¼ 90 . pffiffiffi The flow and temperature profiles, for f ¼ 45 , Lp ¼ 2=4, g ¼ 90 , 45 and 135 , respectively, are presented in Figs. 7–9. Concerning the streamlines, the same remarks for f ¼ 90 remain valid in this orientation (f ¼ 45 ). In addition, we note a strong circulating flow in the cavity. The point of maximum stream function is located at the center of the enclosure for g ¼ 45 and 135 , whereas it is displaced towards
0.46
0.5 0.3
0.2 7.17
0.4
8.51 0.1 0.0 -0.1
13.88
-0.2 -0.3 -0.4 -0.5
(a)
(b) p ffiffi ffi Fig. 4. Isograms at Ra ¼ 5 105 , Lp ¼ 2=4, f ¼ 90 , g ¼ 90 : (a) streamlines; (b) isotherms.
ARTICLE IN PRESS M. Jami et al. / Physica A 368 (2006) 481–494
488
0.5
0.2 3.73
0.4
0.1
0.3
14.92 0.0 -0.1 -0.2 7.93
-0.3
0.93
-0.5
-0.4
(a)
(b) p ffiffi ffi Fig. 5. Isograms at Ra ¼ 5 105 , Lp ¼ 2=4, f ¼ 90 , g ¼ 45 : (a) streamlines; (b) isotherms. 0.70
0.4
0.3
6.09
0.5
0.2 0.1 0.0 -0.1
11.48 -0.2 -0.3 5.01 -0.5
(a)
-0.4
(b) pffiffiffi 2=4, f ¼ 90 , g ¼ 135 : (a) streamlines; (b) isotherms.
-0
.2
10
.4
9
0.
1
-2
.1
-0
0
.1
0.
2
23
.0
8
0.
0
Fig. 6. Isograms at Ra ¼ 5 105 , Lp ¼
0.
(a)
42
45o
45 o
(b) p ffiffi ffi Fig. 7. Isograms at Ra ¼ 5 105 , Lp ¼ 2=4, f ¼ 45 , g ¼ 90 : (a) streamlines; (b) isotherms.
the cold wall for g ¼ 90 . In this last case, there is formation of a co-rotating vortex near the partition. The isotherms plots, near the hot and cold walls, are very similar for those of g ¼ 90 , however, they are very different in the center of the cavity.
ARTICLE IN PRESS M. Jami et al. / Physica A 368 (2006) 481–494
.5
5
2
0.
0.
3 0.
4
-0
0.
.1
25
-0
.2
4
4.
27
0.
0
0.
1
489
.75
1.
45 o
45 o
(b) p ffiffi ffi Fig. 8. Isograms at Ra ¼ 5 105 , Lp ¼ 2=4, f ¼ 45 , g ¼ 45 : (a) streamlines; (b) isotherms.
0.
5
-0
11 .0
9
23
.1
.5
1
0.
0
0.
1
(a)
-0
65
.2
14
1.
47
8.
68
-0
.2
.1 -0
(a)
45 o
45 o
(b) p ffiffi ffi Fig. 9. Isograms at Ra ¼ 5 105 , Lp ¼ 2=4, f ¼ 45 , g ¼ 135 : (a) streamlines; (b) isotherms.
3.1.2. Average Nusselt number The distributions of the average Nusselt number Num pffiffiffi vs Lp for g ¼ 90 , 45 , 135 , f ¼ 90 and Ra ¼ 5 5 p10 ffiffiffi are shown in Fig. 10(a). Lp is varied from 0 to 2=2 (maximal dimensionless partition length) by step of 2=8 for ga90 , and it is varied from 0 to 1 (maximal dimensionless partition length) for g ¼ 90p . The ffiffiffi average Nusselt number decreasespwith increasing L and it becomes almost constant when L exceeds 2=8 p p ffiffi ffi for g ¼ 90 and when it exceeds 3 2=8 for g ¼ 45 or g ¼ 135 , one notices also that the decrease of Num is more pronounced for g ¼ 135 . The profiles of the average Nusselt number Num with respect to Lp for f ¼ 45 are plotted in Fig. 10(b). Num has the same behaviour when the partition is inclined as for f ¼ 90 . When g ¼ 90 , Num presents pffiffiffi pffiffiffi a low value at Lp ¼ 2=8, then it increases up to its maximum value that corresponds to Lp ¼ 3 2=8. From its maximum pffiffiffivalue, Num decreases considerably. In addition, for this last case, the local Nusselt number for Lp ¼ 2=8, as shown in Fig. pffiffiffi 11, is the smallest for small values of Y and along the top part of the hot wall ðY X0:5Þ. For Lp ¼ 3 2=8, the opposite occurs. For each value pffiffiffi of g the local Nusselt number has two local maxima. An example can be seen in Fig 7(a), for Lp ¼ 2=4 where a second vortex appears in the upper half of the cavity and a transition to bicellular flow occurs. In fact, the part of the fluid heated by the hot wall is diverted towards the cold one, producing a local increase in Nu from Y ¼ 0:5 downwards.
ARTICLE IN PRESS 490
M. Jami et al. / Physica A 368 (2006) 481–494
Fig. 10. Variation of the average Nusselt number with Lp for Ra ¼ 5 105 : (a) f ¼ 90 ; (b) f ¼ 45 .
Fig. 11. Local Nusselt number on hot wall at Ra ¼ 5 105 , g ¼ 90 and f ¼ 45 .
3.2. Effect of number of partitions and aspect ratio The results presented here are limited to the inclination angle equal to 45 and cavity aspect ratio equal to 6. After having studied the grid sensitivity test for A ¼ 6 as seen in Table 2, we have estimated that the solution becomes grid independent at Nx Ny ¼ 31 191. Ra is chosen equal to 104 because the solution becomes time dependent for Ra4104 . As a matter of fact, transition to time dependence of the solution depends upon the angle g, the number of partitions N p and the dimensionless partition length Lp . For example, we found 4 that Ra pffiffiffi¼ 10 is the limit for which the calculation reaches a steady state when g ¼ 45 , N p ¼ 5 and Lp ¼ 2=4. 3.2.1. Isotherms and streamlines Figs. 12–15 show the streamlines and isotherms for g ¼ 45 and 135 in the case of N p ¼ 1 and 5, respectively. For the two cases of g, the increase in the number of cells occurring as the number of partitions increases may be explained through the progressive breakdown of the density stratification in the fluid layers adjacent to the active walls that brings to the formation of alternating ‘‘lines’’ of hot and cold fluid circulating from right to left and from left to right across the cavity, with direct effects upon the temperature distribution. One observes a spreading and a compression of the isotherms, of an alternate manner, close to each active wall. It was also noted that around each partition, a significant temperature gradient is established. The isotherms and streamlines are affected in the same manner for g ¼ 45 that for g ¼ 135 . In the first case when
ARTICLE IN PRESS M. Jami et al. / Physica A 368 (2006) 481–494
491
Table 2 pffiffiffi Variation of Num and Cmax according to the mesh size N x N y for Ra ¼ 104 , g ¼ 45 , N p ¼ 5 and Lp ¼ 2=8. 21 131
31 191
41 251
Num Cmax
1.669 9.85
1.6 10.37
1.605 10.4
-0
.5
-0
.1
13
0
.0
6
0.
1
-0
.1
6
13
.0
6
0.
5
Nx Ny
45o
45o
(b)
(a)
pffiffiffi 2=4, f ¼ 45 , g ¼ 45 , N p ¼ 1: (a) streamlines; (b) isotherms.
-0
10
.5
.3
8
-0
.1
0
9.
40
0.
1
0.
5
Fig. 12. Isograms at Ra ¼ 104 , Lp ¼
45o
45o
(a)
(b)
pffiffiffi Fig. 13. Isograms at Ra ¼ 104 , Lp ¼ 2=4, f ¼ 45 , g ¼ 45 , N p ¼ 5: (a) streamlines; (b) isotherms.
N p increases, the circulation becomes intense in the lower part of the cavity, however, in the second case it is intense in the upper part of the cavity. For only one partition a small negative recirculation cell appears towards the center of the cavity. It intensifies when g ¼ 135 , but it remains always secondary with respect to the positive cell. This is due to the secondary buoyancy forces generated by the partition. 3.2.2. Average Nusselt number 4 The results obtained pffiffiffi in this section at f ¼ 45 and Ra ¼ 10 are shown in Fig. 16. For only one partition of small length (Lp p 2=2), the average Nusselt number increases with increasing Lp . This is related with the secondary buoyancy forces generated by the partition. When the partition is hot, these forces cause a local increase of circulation above pffiffiitffi (Fig. 14(a)), with a corresponding increase in the heat transfer. The opposite behaviour occurs for Lp 4 2=2, because the partition comes closer to the cold part of the cavity, so the secondary buoyancy force becomes weak. This behaviour becomes complex in the case of N p ¼ 5. For short
ARTICLE IN PRESS M. Jami et al. / Physica A 368 (2006) 481–494
-0
.5
12
.6
0
-0
.9
8
-0
.1
12
0
.6
0
0.
1
0.
5
492
45o
45o
(b)
(a)
1 0. -0
.5
6. 73
-0
.1
0
9. 03
0.
5
pffiffiffi 2=4, f ¼ 45 , g ¼ 135 , N p ¼ 1: (a) streamlines; (b) isotherms.
8. 26
Fig. 14. Isograms at Ra ¼ 104 , Lp ¼
45o
45o
(a) 4
Fig. 15. Isograms at Ra ¼ 10 ,
Lp
(b) pffiffiffi ¼ 2=4, f ¼ 45 , g ¼ 135 , N p ¼ 5: (a) streamlines; (b) isotherms.
Fig. 16. Variation of the average Nusselt number with Lp at Ra ¼ 104 , f ¼ 45 .
ARTICLE IN PRESS M. Jami et al. / Physica A 368 (2006) 481–494
493
pffiffiffi pffiffiffi (Lp p 2=2) and long (Lp X 2=4) partition, Num decreases with increasing Lp , but it increases between these two values of Lp . 4. Conclusion We have performed 2-d computations of natural convection in a differentially heated inclined enclosure which is fitted with inclined partitions attached to the hot wall. Effects of partition length, aspect ratio of the enclosure, number of partitions, and inclination angles of the partitions and the enclosure have been studied. The main results obtained in the present work may be summarized as: In the case of A ¼ 1: 1. For Ra ¼ 5 105 and for any value of g, the heat transfer is higher in vertical cavities than in inclined ones. 2. The heat transfer reduction increases with increasing Lp when the partition is inclined and whatever the inclination angle pffiffiffi of the cavity. 3. From Lp ¼ 2=8, f ¼ 45 and g ¼ 90 , the average Nusselt number increases up to a maximum value, then it decreases owing to the effect of the secondary buoyancy force generated by the partition. 4. Independently of the value pffiffiffiof f, when the partition is inclined, the average Nusselt number remains almost constant beyond Lp ¼ 3 2=8. In the case of A ¼ 6: 5. At f ¼ 45 , the average Nusselt number decreases with increasing the number of partitions. 6. The decrease of the average Nusselt number Num with increasing Lp is more pronounced for g ¼ 135 than for g ¼ 45 . Finally, we may indicate that extending this work to 3-d situations mainly requires access to more computing power as the present simulation technique has already been tested [17]. In addition consideration of systems with complex geometry can be done with satisfactory accuracy using the boundary conditions recalled in Eqs. (3) and (4).
References [1] J.A. Khan, G.F. Yaho, Comparison of natural convection of water and air in a partitioned rectangular enclosure, Int. J. Heat Mass Transfer 36 (12) (1993) 3107–3117. [2] A. Mezrhab, L. Bchir, Radiation–natural convection interactions in partitioned cavities, Int. J. Numer. Methods Heat Fluid Flow 9 (2) (1999) 186–203. [3] R.L. Frederick, Natural convection in an inclined square enclosure with a partition attached to its cold wall, Int J. Heat Mass Transfer 32 (1989) 87–94. [4] R. Anderson, G. Lauriat, The horizontal natural convection boundary layer regime in a closed cavity, in: Proceedings of the Eighth International Heat Transfer Conference, vol. 4, San Francisco, CA, 1986, pp. 1453–1458. [5] M.M. Ganzarolli, L.F. Milanez, Natural convection in rectangular enclosures heated from below and symmetrically cooled from the sides, J. Heat Mass Transfer 38 (1995) 1063–1073. [6] O. Aydin, A. Unal, T. Ayhan, Numerical solutions for buoyancy-driven flow in a 2-D square enclosure heated from one side and cooled from above, in: Proceedings of the Advances in Computational Heat Transfer Symposium, Begell House, New York, 1997, pp. 387–394. [7] O. Aydin, A. Unal, T. Ayhan, Natural convection in rectangular enclosures heated from one side and cooled from the ceiling, Int. J. Heat Mass Transfer 42 (1999) 2345–2355. [8] S. Acharya, R. Jetli, Heat transfer due to buoyancy in a partially divided square box, Int. J. Heat Mass Transfer 33 (5) (1990) 931–942. [9] A. Mezrhab, M. Jami, C. Abid, M. Bouzidi, P. Lallemand, Lattice-Boltzmann modelling of natural convection in an inclined square enclosure with partitions attached to its cold wall, Int. J. Heat Fluid Flow, 2005, in press. [10] S. Chen, G.D. Doolen, Lattice Boltzmann method for fluid flows, Ann. Rev. Fluid Mech. 30 (1998) 329–364. [11] D. d’Humie`res, Generalized lattice-Boltzmann equations, in: B.D. Shizgal, D.P. Weaver (Eds.), Rarefied Gas Dynamics: Theory and Simulations, Progress in Astronautics and Aeronautics, vol. 159, AIAA, Washington, DC, 1992, pp. 450–458. [12] L. Zhenyu, L. Ying, Q. Dongying, J.B. McLaughlin, J.J. Derksen, K. Kontomaris, Large eddy simulations of a stirred tank using the lattice Boltzmann method on a nonuniform grid, J. Comput. Phys. 181 (2002) 675–704.
ARTICLE IN PRESS 494
M. Jami et al. / Physica A 368 (2006) 481–494
[13] I. Ginzbourg, P.M. Adler, Boundary flow condition analysis for the three-dimensional lattice Boltzmann model, J. Phys. II France 4 (1994) 191–214. [14] P. Lallemand, L.S. Luo, Lattice Boltzmann method for moving boundaries, J. Comput. Phys. 184 (2003) 406–421. [15] M. Bouzidi, M. Firdaouss, P. Lallemand, Momentum transfer of a Boltzmann lattice fluid with boundaries, Phys. Fluids 13 (2001) 3452–3459. [16] D.R. Noble, S. Chen, J.G. Georgiadis, R.O. Buckius, A consistent hydrodynamic boundary condition for the lattice Boltzmann method, Phys. Fluids 7 (1995) 203–209. [17] A. Mezrhab, M. Bouzidi, P. Lallemand, Hybrid lattice Boltzmann finite-difference simulation of convective flows, Comput. Fluids 33 (2004) 623–641. [18] F.J. Alexander, S. Chen, J.D. Sterling, Lattice Boltzmann thermodynamics, Phys. Rev. (E47) (1993) R2249–R2252. [19] P. Pavlo, G. Vahala, L. Vahala, M. Soe, Linear stability analysis of thermo-lattice Boltzmann models, J. Comput. Phys. 139 (1998) 79–91. [20] P. Lallemand, L.S. Luo, Theory of the lattice Boltzmann method: acoustic and thermal properties in two and three dimensions, Phys. Rev. (E68) (2003) 036706-1–0367025-1.
Physica A 368 (2006) 481–494 www.elsevier.com/locate/physa
Lattice-Boltzmann computation of natural convection in a partitioned enclosure with inclined partitions attached to its hot wall Mohammed Jamia, Ahmed Mezrhaba,, M’hamed Bouzidib, Pierre Lallemandc a
De´partement de Physique, Faculte´ des Sciences, Universite´ Mohamed I, Oujda, Morocco b Universite´ Clermont-2, I.U.T. Ave A. Briand, 03100 Montluc- on, France c LIMSI, Baˆtiment 508, Universite´ Paris-Sud, 91405 Orsay, France Received 7 October 2005 Available online 13 January 2006
Abstract This study presents a numerical investigation of laminar natural convection heat transfer in an inclined enclosure, differentially heated, with inclined partitions attached to its hot wall. The top and bottom walls of the cavity are insulated. Numerical solutions are obtained by using the hybrid lattice-Boltzmann method finite-difference. Effects of partition length, partition inclination angle, partitions number and aspect ratio cavity are studied for Prandtl number 0.71 and cavity inclination angles of 45 and 90 . Rayleigh numbers considered are 5 105 for a square cavity and 104 for a cavity with an aspect ratio equal to 6. Results are reported in terms of isotherms, streamlines, local and average Nusselt numbers. r 2006 Elsevier B.V. All rights reserved.
1. Introduction Natural convection in a partitioned cavity has been extensively studied using numerical simulations and experiments because of its interest and importance in industrial applications. Some applications are solar collectors, fire research, electronic cooling, aeronautics, chemical apparatus and building engineering. Most of the papers in this field are substantially oriented toward the study of natural convection in enclosed square or rectangular cavities [1]. Among these works, we mention the one of Mezrhab and Bchir [2] who used the finitevolume method to study the effect of adding a thick partition located vertically close to the hot wall of a differentially heated cavity, forming a narrow vertical channel in which the flow is controlled by vents at the bottom and top of the partition. They have shown that radiation has a significant influence on the flow and heat transfer in the channel. Frederick [3] has studied numerically natural convection in an air-filled, differentially heated, inclined square cavity, with a single partition attached to its cold wall, at Rayleigh numbers of 103 –105 . He showed that the partition leads to the suppression of convection, and reduces the heat transfer by up to 47% in comparison to the empty cavity at the same Rayleigh number. He used Corresponding author. Tel.: +212 5650 0601; fax: +212 5650 0603.
E-mail addresses: [email protected] (M. Jami), [email protected], [email protected] (A. Mezrhab), [email protected] (M. Bouzidi), [email protected] (P. Lallemand). 0378-4371/$ - see front matter r 2006 Elsevier B.V. All rights reserved. doi:10.1016/j.physa.2005.12.029
ARTICLE IN PRESS M. Jami et al. / Physica A 368 (2006) 481–494
482
Nomenclature A g H k L Lp Lp Np Num Nu Pr Ra Tc Th T0 u v x; y X; Y
cavity aspect ratio, A ¼ H=L acceleration of gravity enclosure height thermal conductivity of the fluid enclosure width partition length dimensionless partition length, Lp ¼ Lp =L number of partitions RH average Nusselt number ð1=AÞ 0 qy=qx wall dy local Nusselt number, L qy=qx wall Prandtl number, Pr ¼ n=a Rayleigh number, gb DTL3 =na temperature of the cold wall temperature of the hot wall average fluid temperature ðT h þ T c Þ=2 x velocity component y velocity component Cartesian coordinates dimensionless Cartesian coordinates, X ¼ x=L, Y ¼ y=L
Greek symbols a thermal diffusivity of the fluid b volumetric expansion coefficient d partition thickness d dimensionless partition thickness, d ¼ d=L DT maximal difference temperature, T h T c DX ; DY dimensionless grid spacings in X and Y directions g partition inclination n kinematic viscosity of the fluid f inclination angle of the enclosure with respect to the horizontal c stream function C dimensionless stream function, c=a y dimensionless temperature ðT T 0 Þ=DT Subscripts max pþ p
maximal value right surface of the partition left surface of the partition
a finite-difference over-relaxation procedure for the solution of the mass, momentum and energy transfer governing equations. Let us mention that studies on natural convection in rectangular enclosures heated from below and cooled along a single side or both sides have been carried out respectively by Anderson and Lauriat [4] and by Ganzarolli and Milanez [5]. More recently, the case of heating from one side and cooling from the top has been analyzed by Aydin et al. [6] who investigated the effect of Prandtl number upon heat and momentum transfer inside square cavities. They also studied the influence of the aspect ratio for air-filled, rectangular enclosures [7]. Acharya and Jetli [8] have studied numerically the heat transfer and the fluid flow
ARTICLE IN PRESS M. Jami et al. / Physica A 368 (2006) 481–494
483
within a square cavity divided by a single obstacle, attached to the ceiling or the floor. They considered different positions and heights of the partition. They showed that the heat transfer is strongly influenced by the height of the partition; nevertheless, its position has a rather weak effect on the total heat transfer. Finally, we mention that in a similar study of natural convection in an inclined square enclosure with partitions attached to its cold wall, Mezrhab et al. [9] concluded that at lower Rayleigh number Ra (Rap105 ), the average Nusselt number Num is higher in inclined cavities than in vertical ones; while at larger Ra (Ra ¼ 106 ), the opposite occurs. In vertical cavities, Num decreases with increasing the dimensionless length of the partition Lp . However in inclined ones, this trend is found only at lower Ra, whereas at higher Ra, Num shows a maximum at about Lp ¼ 0:5. Num also decreases when the number of partitions ðN p Þ attached to the cold wall of the enclosure increases; however, beyond a certain value of ðN p Þ, which depends of Ra, Num becomes almost constant. The objective of this paper is to study numerically the heat transfer in an inclined square enclosure with inclined partitions attached to its hot wall. The numerical method proposed is based on the coupling between the lattice-Boltzmann method and finite-difference. 2. Numerical technique and mathematical formulation 2.1. Physical system Fig. 1 presents the geometry and the boundary conditions of a two-dimensional square enclosure placed in a gravity field. The laminar air flow can be described by assuming that the fluid is Newtonian and satisfies the Boussinesq approximation, it is initially stationary at a uniform temperature T 0 . We consider a cavity with two walls maintained respectively at T c and T h , making an angle f with the horizontal. The other two walls are adiabatic, midway between them, an inclined partition of length Lp and thickness d is fitted to the hot wall. Concerning applications of this study to systems such as solar collectors, the partition should be very thin and it is taken as a single line of mesh points. Thus, the dimensionless partition thickness d was fixed to one grid spacing. 2.2. Lattice-Boltzmann equation The simulations have been performed using a combination of lattice-Boltzmann equation (LBE) and finitedifference scheme for the temperature with suitable coupling between them to describe thermo-convective effects. LBE has been developed over the last 20 years as a simple means to simulate fluid flows [10]. It is
Fig. 1. Square cavity simulated.
ARTICLE IN PRESS M. Jami et al. / Physica A 368 (2006) 481–494
484
inspired by the kinetic theory of gases which considers a function f ðr; v; tÞ that describes the population of particles of the gas at location r, with velocity v at time t. In LBE space and time are discretized on simple uniform grids, and the possible particle velocities are limited to a small number chosen such that particles move exactly from one node of the grid to a neighbouring one during one time step. In the usual LBE models, the space grid is taken square or cubic, respectively for 2-d and 3-d problems. The number of velocities is chosen to include no move, move to the nearest neighbour and move to the next nearest neighbour. As we describe here results for 2-d situations, we take the simple D2Q9 model with 9 velocities (example: node located at the center of the Fig. 2). The variables are thus f i ðrj ; tl Þ;
i 2 f0; . . . ; 8g,
(1)
which are symbolically represented by the vector F, with rj a node of a square lattice of unit size DX , and tl ¼ l dt. All quantities are normalized in terms of the proper combination of DX and dt. The equation of motion for the quantities f i is given by f i ðrj þ ci ; t þ 1Þ f i ðrj ; tÞ ¼ OF ,
(2)
where the left-hand side represents the motion of particles (treated exactly with the particular choice of elementary velocities ci ) and the right-hand side describes particle redistribution at each lattice node due to collisions. The operator O is highly simplified compared to what would be required to take into account real collisions but it can be chosen in such a way that the large-scale behaviour of the model is similar to that of a real viscous fluid following the Navier–Stokes equations. In this work we use the technique proposed by d’Humie`res [11] (through the computation of moments) nowadays known as the MRT–LBE technique (multirelaxation time) as it is numerically more stable than the simple LBE–BGK technique and can be tuned to give more satisfactory behaviour at boundaries. Compared to a standard finite-difference code one may say that LBE uses more degrees of freedom (here 9 per node). The operator O is such that some combinations of f i are conserved (mass and linear momentum) whereas others follow some relaxation equation. It is the relaxation of these combinations that is responsible for the actual values of the speed of sound and of the transport coefficients. Apart from its simplicity, LBE allows to implement boundaries without much difficulty.
Fig. 2. Grid distributions of computational domain.
ARTICLE IN PRESS M. Jami et al. / Physica A 368 (2006) 481–494
485
The simplest boundary conditions so-called ‘‘bounce-back’’ allow to implement Dirichlet-type boundaries located to a good approximation half-way between the last ‘‘fluid node’’ and the first ‘‘solid node’’ [12,13]. To represent boundaries that cross a link between fluid and solid that is located at a relative distance q (see Fig. 3), one can use a variety of approximations. Here, we consider a simple combination of ‘‘bounce-back’’ together with quadratic interpolation for each particle link that crosses the interface [14,15]. The LBE part of our code yields values of the fluid velocity uðtÞ, vðtÞ at the nodes of the lattice. It is straightforward to add the effect of a buoyancy force through modification of the collision process [16]: instead of conserving linear momentum, we apply a change of the components of the linear momentum proportional to the buoyancy gbðT T 0 Þ. The dynamics of the temperature is computed using the advection–diffusion equation with finite-difference stencils to compute the space derivatives. Time evolution is done with a simple second-order predictor–corrector scheme. As we use the same thermal diffusivity for the fluid and the material used to make the partitions, there is no special treatment of the partitions (LBE automatically gives 0 for the fluid velocity inside the partition). In practice, we take the external walls of the cavity parallel to the axis of the LBE-grid and apply gravity at an angle chosen to satisfy the geometry of Fig. 1. This allows us to use simple ‘‘bounce-back’’ on the external boundary of the fluid, therefore, the fluid velocity is 0 at mid-distance between the last fluid node and the first solid node, that is along the axis Ox at x ¼ 12 and x ¼ N x þ 12 and at y ¼ 12 and y ¼ N y þ 12 along the axis Oy as shown in Fig. 2. On the interface between the fluid and the partition we have used a combination of the bounce-back scheme and interpolations because there are some particles that are going to be in interstitial positions after being rebounded. This technique is graphically depicted in Fig. 3. We have used the following interpolation formulas [14]: f ¯i ðxf ; tÞ ¼ qð1 þ 2qÞf i ðxf þ ci ; tÞ þ ð1 4q2 Þf i ðxf ; tÞ qð1 2qÞf i ðxf ci ; tÞ;
f ¯i ðxf ; tÞ ¼
qo1=2
1 ð2q 1Þ ð2q 1Þ f i ðxf þ ci ; tÞ þ f ¯i ðxf ci ; tÞ f ¯ ðxf 2ci ; tÞ; qð1 þ 2qÞ q ð2q þ 1Þ i
qX1=2
(3)
(4)
where f ¯i is the distribution function of the velocity c¯i ci . The fluid velocity on the enclosure walls and the partition is taken equal to 0.
/ / (a)
/ (b)
/ (c)
Fig. 3. Illustration of the boundary conditions for a rigid wall located between two grid sites in the direction of the motion of the particle that is rebounded.
ARTICLE IN PRESS M. Jami et al. / Physica A 368 (2006) 481–494
486
2.3. Equation of temperature The energy equation is solved by the finite-difference, in order to provide us with the field of temperature inside the cavity. The energy equation for the present system can be expressed by the following 2-d equation: 2 qy qy qy q y q2 y þu þv ¼a þ . (5) qt qx qy qx2 qy2 The coupling between the flow and the temperature field is done as follows: dimensionless temperature y provides a buoyancy force in the x and y directions because the square cavity is inclined by an angle f with respect to the horizontal. Therefore, f x ¼ gb DTyðr; tÞ sin f and f y ¼ gb DTyðr; tÞ cos f where yðr; tÞ is the local dimensionless temperature at the position of the air node r at time t. Let us note that the terms introduced to couple flow and temperature neglect compressibility. This procedure gives very good results in the simulation of the convective flows as that was shown by Mezrhab et al. [17]. As a matter of fact, many attempts have been made to simulate the thermal effects with LBM techniques [18]. However, this path presented difficulties as far as isotropy and stability [19,20] are concerned. Thus, it is more convenient to treat separately density and momentum (athermal fluid) and temperature. As the temperature grid nodes are the same as the fluid nodes, this leads to consider ‘‘ghost’’ nodes outside the fluid at x ¼ 0 and N x þ 1 and at y ¼ 0 and N y þ 1 with a temperature computed by extrapolating the information on the solid boundary (Dirichlet or Neuman) and the first and second point inside the fluid. The thermal boundary conditions at the walls are as follows: yð1=2; yÞ ¼ yc ¼ 0:5 yðN x þ 1=2; yÞ ¼ yh ¼ 0:5 qy ðx; 1=2Þ ¼ 0 qy qy ðx; N y þ 1=2Þ ¼ 0 qy
1 1 pypN y þ ; 2 2 1 1 for pypN y þ ; 2 2 1 1 for pxpN x þ ; 2 2 1 1 for pxpN x þ : 2 2 for
At the partition, the conditions of continuity of temperature and temperature gradient are imposed. These conditions are written as ypþ ¼ yp ;
qy qy ¼ , qzp qzpþ
(6)
where z is the normal direction to the partition. 3. Results and discussions The parameters governing the problem are A; Lp ; N p ; Pr; Ra; f and g. Except for the Prandtl number that was fixed at Pr ¼ 0:71, the other parameters were varied with the aim of studying their effect on the thermal transfer and the flow in the cavity. A preliminary study was carried out to determine the optimum uniform grid (i.e. the best compromise between accuracy and computational costs). In the case of a square cavity with a single partition attached to its hot wall, Table 1 presents the average Nusselt number
Table 1 Variation of Num and Cmax according to the mesh size N x N y Nx Ny
41 41
61 61
81 81
101 101
Num Cmax
5.759 24.90
5.898 25.18
5.917 25.24
5.917 25.25
ARTICLE IN PRESS M. Jami et al. / Physica A 368 (2006) 481–494
487
pffiffiffi and maximal dimensionless stream function as a function of the mesh size N x N y for A ¼ 1, Lp ¼ 2=4, Ra ¼ 5 105 , f ¼ 45 and g ¼ 45 . As can be seen, the solutions become grid independent at Nx Ny ¼ 81 81. Consequently, for the computations reported in the case of a square cavity, a grid with 81 81 points was chosen to optimize the relation between accuracy and computing time. The code was extensively exercised on benchmark problems to check its validity. In the cases of an empty cavity, a square cavity with a hot circle inside, and concentric cylinders, the results that we obtain with our simulation technique were checked for accuracy against the earlier published numerical and experimental results reported by different authors, and the agreement between the present and previous results was very good in Ref. [17]. The code was also validated for cavities divided by a thin partition as shown in Ref. [9]. It was concluded that the largest discrepancies between our and published results can be estimated to be less than 1%. For this reason, and for the sake of brevity it is not repeated here. Based on the above studies, it was concluded that the code could be reliably applied to the problem considered here. 3.1. Square cavity with a single partition 3.1.1. Isotherms and streamlines pffiffiffi Figs. 4–6 show the streamlines and isotherms for f ¼ 90 , Lp ¼ 2=4, g ¼ 90 , 45 and 135 , respectively. Let us note that near the partition and the cold wall, the flow is intense, the fluid rotates counter-clockwise owing to the positions of the hot and cold walls. Whatever the partition orientation, the isotherms and streamlines are considerably altered, in comparison to those corresponding to the empty cavity at the same Rayleigh number. In fact, when the partition is inclined by 45 , the major part of the fluid is heated by the lower part of the hot wall. Then it is reflected toward a region close to the cold wall, the rest of the fluid recirculates in the upper right corner of the enclosure, following the same direction of the main circulation. In this case the circulation is stronger than for the two others cases of g. For the partition inclination (g ¼ 135 ), the major part of the fluid circulates near the walls. The circulation and the heat transfer decrease and the point of maximum stream function is displaced towards the left part of the cavity. Independently of the value of g, from the temperature contours, it is apparent that thermal boundary layers form along the hot and cold walls. Far from the partition, the isotherms are nearly horizontal and resemble to those obtained in the case of an empty cavity, so the temperature stratification persists; however, around the partition the isotherms become more affected. We note that the streamlines and isotherms corresponding to g ¼ 45 and 135 are very different from those corresponding to g ¼ 90 . pffiffiffi The flow and temperature profiles, for f ¼ 45 , Lp ¼ 2=4, g ¼ 90 , 45 and 135 , respectively, are presented in Figs. 7–9. Concerning the streamlines, the same remarks for f ¼ 90 remain valid in this orientation (f ¼ 45 ). In addition, we note a strong circulating flow in the cavity. The point of maximum stream function is located at the center of the enclosure for g ¼ 45 and 135 , whereas it is displaced towards
0.46
0.5 0.3
0.2 7.17
0.4
8.51 0.1 0.0 -0.1
13.88
-0.2 -0.3 -0.4 -0.5
(a)
(b) p ffiffi ffi Fig. 4. Isograms at Ra ¼ 5 105 , Lp ¼ 2=4, f ¼ 90 , g ¼ 90 : (a) streamlines; (b) isotherms.
ARTICLE IN PRESS M. Jami et al. / Physica A 368 (2006) 481–494
488
0.5
0.2 3.73
0.4
0.1
0.3
14.92 0.0 -0.1 -0.2 7.93
-0.3
0.93
-0.5
-0.4
(a)
(b) p ffiffi ffi Fig. 5. Isograms at Ra ¼ 5 105 , Lp ¼ 2=4, f ¼ 90 , g ¼ 45 : (a) streamlines; (b) isotherms. 0.70
0.4
0.3
6.09
0.5
0.2 0.1 0.0 -0.1
11.48 -0.2 -0.3 5.01 -0.5
(a)
-0.4
(b) pffiffiffi 2=4, f ¼ 90 , g ¼ 135 : (a) streamlines; (b) isotherms.
-0
.2
10
.4
9
0.
1
-2
.1
-0
0
.1
0.
2
23
.0
8
0.
0
Fig. 6. Isograms at Ra ¼ 5 105 , Lp ¼
0.
(a)
42
45o
45 o
(b) p ffiffi ffi Fig. 7. Isograms at Ra ¼ 5 105 , Lp ¼ 2=4, f ¼ 45 , g ¼ 90 : (a) streamlines; (b) isotherms.
the cold wall for g ¼ 90 . In this last case, there is formation of a co-rotating vortex near the partition. The isotherms plots, near the hot and cold walls, are very similar for those of g ¼ 90 , however, they are very different in the center of the cavity.
ARTICLE IN PRESS M. Jami et al. / Physica A 368 (2006) 481–494
.5
5
2
0.
0.
3 0.
4
-0
0.
.1
25
-0
.2
4
4.
27
0.
0
0.
1
489
.75
1.
45 o
45 o
(b) p ffiffi ffi Fig. 8. Isograms at Ra ¼ 5 105 , Lp ¼ 2=4, f ¼ 45 , g ¼ 45 : (a) streamlines; (b) isotherms.
0.
5
-0
11 .0
9
23
.1
.5
1
0.
0
0.
1
(a)
-0
65
.2
14
1.
47
8.
68
-0
.2
.1 -0
(a)
45 o
45 o
(b) p ffiffi ffi Fig. 9. Isograms at Ra ¼ 5 105 , Lp ¼ 2=4, f ¼ 45 , g ¼ 135 : (a) streamlines; (b) isotherms.
3.1.2. Average Nusselt number The distributions of the average Nusselt number Num pffiffiffi vs Lp for g ¼ 90 , 45 , 135 , f ¼ 90 and Ra ¼ 5 5 p10 ffiffiffi are shown in Fig. 10(a). Lp is varied from 0 to 2=2 (maximal dimensionless partition length) by step of 2=8 for ga90 , and it is varied from 0 to 1 (maximal dimensionless partition length) for g ¼ 90p . The ffiffiffi average Nusselt number decreasespwith increasing L and it becomes almost constant when L exceeds 2=8 p p ffiffi ffi for g ¼ 90 and when it exceeds 3 2=8 for g ¼ 45 or g ¼ 135 , one notices also that the decrease of Num is more pronounced for g ¼ 135 . The profiles of the average Nusselt number Num with respect to Lp for f ¼ 45 are plotted in Fig. 10(b). Num has the same behaviour when the partition is inclined as for f ¼ 90 . When g ¼ 90 , Num presents pffiffiffi pffiffiffi a low value at Lp ¼ 2=8, then it increases up to its maximum value that corresponds to Lp ¼ 3 2=8. From its maximum pffiffiffivalue, Num decreases considerably. In addition, for this last case, the local Nusselt number for Lp ¼ 2=8, as shown in Fig. pffiffiffi 11, is the smallest for small values of Y and along the top part of the hot wall ðY X0:5Þ. For Lp ¼ 3 2=8, the opposite occurs. For each value pffiffiffi of g the local Nusselt number has two local maxima. An example can be seen in Fig 7(a), for Lp ¼ 2=4 where a second vortex appears in the upper half of the cavity and a transition to bicellular flow occurs. In fact, the part of the fluid heated by the hot wall is diverted towards the cold one, producing a local increase in Nu from Y ¼ 0:5 downwards.
ARTICLE IN PRESS 490
M. Jami et al. / Physica A 368 (2006) 481–494
Fig. 10. Variation of the average Nusselt number with Lp for Ra ¼ 5 105 : (a) f ¼ 90 ; (b) f ¼ 45 .
Fig. 11. Local Nusselt number on hot wall at Ra ¼ 5 105 , g ¼ 90 and f ¼ 45 .
3.2. Effect of number of partitions and aspect ratio The results presented here are limited to the inclination angle equal to 45 and cavity aspect ratio equal to 6. After having studied the grid sensitivity test for A ¼ 6 as seen in Table 2, we have estimated that the solution becomes grid independent at Nx Ny ¼ 31 191. Ra is chosen equal to 104 because the solution becomes time dependent for Ra4104 . As a matter of fact, transition to time dependence of the solution depends upon the angle g, the number of partitions N p and the dimensionless partition length Lp . For example, we found 4 that Ra pffiffiffi¼ 10 is the limit for which the calculation reaches a steady state when g ¼ 45 , N p ¼ 5 and Lp ¼ 2=4. 3.2.1. Isotherms and streamlines Figs. 12–15 show the streamlines and isotherms for g ¼ 45 and 135 in the case of N p ¼ 1 and 5, respectively. For the two cases of g, the increase in the number of cells occurring as the number of partitions increases may be explained through the progressive breakdown of the density stratification in the fluid layers adjacent to the active walls that brings to the formation of alternating ‘‘lines’’ of hot and cold fluid circulating from right to left and from left to right across the cavity, with direct effects upon the temperature distribution. One observes a spreading and a compression of the isotherms, of an alternate manner, close to each active wall. It was also noted that around each partition, a significant temperature gradient is established. The isotherms and streamlines are affected in the same manner for g ¼ 45 that for g ¼ 135 . In the first case when
ARTICLE IN PRESS M. Jami et al. / Physica A 368 (2006) 481–494
491
Table 2 pffiffiffi Variation of Num and Cmax according to the mesh size N x N y for Ra ¼ 104 , g ¼ 45 , N p ¼ 5 and Lp ¼ 2=8. 21 131
31 191
41 251
Num Cmax
1.669 9.85
1.6 10.37
1.605 10.4
-0
.5
-0
.1
13
0
.0
6
0.
1
-0
.1
6
13
.0
6
0.
5
Nx Ny
45o
45o
(b)
(a)
pffiffiffi 2=4, f ¼ 45 , g ¼ 45 , N p ¼ 1: (a) streamlines; (b) isotherms.
-0
10
.5
.3
8
-0
.1
0
9.
40
0.
1
0.
5
Fig. 12. Isograms at Ra ¼ 104 , Lp ¼
45o
45o
(a)
(b)
pffiffiffi Fig. 13. Isograms at Ra ¼ 104 , Lp ¼ 2=4, f ¼ 45 , g ¼ 45 , N p ¼ 5: (a) streamlines; (b) isotherms.
N p increases, the circulation becomes intense in the lower part of the cavity, however, in the second case it is intense in the upper part of the cavity. For only one partition a small negative recirculation cell appears towards the center of the cavity. It intensifies when g ¼ 135 , but it remains always secondary with respect to the positive cell. This is due to the secondary buoyancy forces generated by the partition. 3.2.2. Average Nusselt number 4 The results obtained pffiffiffi in this section at f ¼ 45 and Ra ¼ 10 are shown in Fig. 16. For only one partition of small length (Lp p 2=2), the average Nusselt number increases with increasing Lp . This is related with the secondary buoyancy forces generated by the partition. When the partition is hot, these forces cause a local increase of circulation above pffiffiitffi (Fig. 14(a)), with a corresponding increase in the heat transfer. The opposite behaviour occurs for Lp 4 2=2, because the partition comes closer to the cold part of the cavity, so the secondary buoyancy force becomes weak. This behaviour becomes complex in the case of N p ¼ 5. For short
ARTICLE IN PRESS M. Jami et al. / Physica A 368 (2006) 481–494
-0
.5
12
.6
0
-0
.9
8
-0
.1
12
0
.6
0
0.
1
0.
5
492
45o
45o
(b)
(a)
1 0. -0
.5
6. 73
-0
.1
0
9. 03
0.
5
pffiffiffi 2=4, f ¼ 45 , g ¼ 135 , N p ¼ 1: (a) streamlines; (b) isotherms.
8. 26
Fig. 14. Isograms at Ra ¼ 104 , Lp ¼
45o
45o
(a) 4
Fig. 15. Isograms at Ra ¼ 10 ,
Lp
(b) pffiffiffi ¼ 2=4, f ¼ 45 , g ¼ 135 , N p ¼ 5: (a) streamlines; (b) isotherms.
Fig. 16. Variation of the average Nusselt number with Lp at Ra ¼ 104 , f ¼ 45 .
ARTICLE IN PRESS M. Jami et al. / Physica A 368 (2006) 481–494
493
pffiffiffi pffiffiffi (Lp p 2=2) and long (Lp X 2=4) partition, Num decreases with increasing Lp , but it increases between these two values of Lp . 4. Conclusion We have performed 2-d computations of natural convection in a differentially heated inclined enclosure which is fitted with inclined partitions attached to the hot wall. Effects of partition length, aspect ratio of the enclosure, number of partitions, and inclination angles of the partitions and the enclosure have been studied. The main results obtained in the present work may be summarized as: In the case of A ¼ 1: 1. For Ra ¼ 5 105 and for any value of g, the heat transfer is higher in vertical cavities than in inclined ones. 2. The heat transfer reduction increases with increasing Lp when the partition is inclined and whatever the inclination angle pffiffiffi of the cavity. 3. From Lp ¼ 2=8, f ¼ 45 and g ¼ 90 , the average Nusselt number increases up to a maximum value, then it decreases owing to the effect of the secondary buoyancy force generated by the partition. 4. Independently of the value pffiffiffiof f, when the partition is inclined, the average Nusselt number remains almost constant beyond Lp ¼ 3 2=8. In the case of A ¼ 6: 5. At f ¼ 45 , the average Nusselt number decreases with increasing the number of partitions. 6. The decrease of the average Nusselt number Num with increasing Lp is more pronounced for g ¼ 135 than for g ¼ 45 . Finally, we may indicate that extending this work to 3-d situations mainly requires access to more computing power as the present simulation technique has already been tested [17]. In addition consideration of systems with complex geometry can be done with satisfactory accuracy using the boundary conditions recalled in Eqs. (3) and (4).
References [1] J.A. Khan, G.F. Yaho, Comparison of natural convection of water and air in a partitioned rectangular enclosure, Int. J. Heat Mass Transfer 36 (12) (1993) 3107–3117. [2] A. Mezrhab, L. Bchir, Radiation–natural convection interactions in partitioned cavities, Int. J. Numer. Methods Heat Fluid Flow 9 (2) (1999) 186–203. [3] R.L. Frederick, Natural convection in an inclined square enclosure with a partition attached to its cold wall, Int J. Heat Mass Transfer 32 (1989) 87–94. [4] R. Anderson, G. Lauriat, The horizontal natural convection boundary layer regime in a closed cavity, in: Proceedings of the Eighth International Heat Transfer Conference, vol. 4, San Francisco, CA, 1986, pp. 1453–1458. [5] M.M. Ganzarolli, L.F. Milanez, Natural convection in rectangular enclosures heated from below and symmetrically cooled from the sides, J. Heat Mass Transfer 38 (1995) 1063–1073. [6] O. Aydin, A. Unal, T. Ayhan, Numerical solutions for buoyancy-driven flow in a 2-D square enclosure heated from one side and cooled from above, in: Proceedings of the Advances in Computational Heat Transfer Symposium, Begell House, New York, 1997, pp. 387–394. [7] O. Aydin, A. Unal, T. Ayhan, Natural convection in rectangular enclosures heated from one side and cooled from the ceiling, Int. J. Heat Mass Transfer 42 (1999) 2345–2355. [8] S. Acharya, R. Jetli, Heat transfer due to buoyancy in a partially divided square box, Int. J. Heat Mass Transfer 33 (5) (1990) 931–942. [9] A. Mezrhab, M. Jami, C. Abid, M. Bouzidi, P. Lallemand, Lattice-Boltzmann modelling of natural convection in an inclined square enclosure with partitions attached to its cold wall, Int. J. Heat Fluid Flow, 2005, in press. [10] S. Chen, G.D. Doolen, Lattice Boltzmann method for fluid flows, Ann. Rev. Fluid Mech. 30 (1998) 329–364. [11] D. d’Humie`res, Generalized lattice-Boltzmann equations, in: B.D. Shizgal, D.P. Weaver (Eds.), Rarefied Gas Dynamics: Theory and Simulations, Progress in Astronautics and Aeronautics, vol. 159, AIAA, Washington, DC, 1992, pp. 450–458. [12] L. Zhenyu, L. Ying, Q. Dongying, J.B. McLaughlin, J.J. Derksen, K. Kontomaris, Large eddy simulations of a stirred tank using the lattice Boltzmann method on a nonuniform grid, J. Comput. Phys. 181 (2002) 675–704.
ARTICLE IN PRESS 494
M. Jami et al. / Physica A 368 (2006) 481–494
[13] I. Ginzbourg, P.M. Adler, Boundary flow condition analysis for the three-dimensional lattice Boltzmann model, J. Phys. II France 4 (1994) 191–214. [14] P. Lallemand, L.S. Luo, Lattice Boltzmann method for moving boundaries, J. Comput. Phys. 184 (2003) 406–421. [15] M. Bouzidi, M. Firdaouss, P. Lallemand, Momentum transfer of a Boltzmann lattice fluid with boundaries, Phys. Fluids 13 (2001) 3452–3459. [16] D.R. Noble, S. Chen, J.G. Georgiadis, R.O. Buckius, A consistent hydrodynamic boundary condition for the lattice Boltzmann method, Phys. Fluids 7 (1995) 203–209. [17] A. Mezrhab, M. Bouzidi, P. Lallemand, Hybrid lattice Boltzmann finite-difference simulation of convective flows, Comput. Fluids 33 (2004) 623–641. [18] F.J. Alexander, S. Chen, J.D. Sterling, Lattice Boltzmann thermodynamics, Phys. Rev. (E47) (1993) R2249–R2252. [19] P. Pavlo, G. Vahala, L. Vahala, M. Soe, Linear stability analysis of thermo-lattice Boltzmann models, J. Comput. Phys. 139 (1998) 79–91. [20] P. Lallemand, L.S. Luo, Theory of the lattice Boltzmann method: acoustic and thermal properties in two and three dimensions, Phys. Rev. (E68) (2003) 036706-1–0367025-1.