Abstract
Natural convection of gas flow (air) confined within an enclosed square-section cavity is investigated numerically using the lattice Boltzmann method (LBM). The right (left) side of the enclosure is partially heated (cooled) by a hot (cold) chip, while the left (right) one is completely kept at cold (hot) temperature. However, the horizontal walls and vertical parts near the chip are kept adiabatic. The buoyancy effect induced by the gravity acceleration, related to the convection force, is evaluated through the Rayleigh number in the range of (laminar regime). The wall heating-ratio effect on the flow properties such as temperature and velocity profiles was examined. The heat transfer is analyzed through the Nusselt number for different chip lengths. Results show that the wall heat ratio has an interesting effect on the flow behavior. Results show good agreement with those of full natural convection in the literature, experimental, and simulation data.
1. Introduction
Modern framework-based systems including gas flows at various sizes of uses are in growth development. In the recent few decades, the lattice Boltzmann method was successfully used as a mesoscopic approach based on both continuum and kinetic description of fluid flows. Hence, regarding its extended applicability, LBM has become a recent powerful computational fluid dynamics (CFD) tool. Its popularity is gained in light of the fact that the technique can be stretched out to different applications, unlike solvers based on the Navier-Stokes equations which are limited to the pure macroscopic description of flows. Using the LBM models, the simulation of incompressible flows and heat transfer problems can be typically carried out. Both governing equations of density and temperature distribution functions are solved according to a specified arrangement in lattices. Depending on the algorithm and/or boundary conditions used, different varieties of LBM can be found in the literature [1–4]. In the hybrid schemes, besides the classical LBM used to describe the hydrodynamic flow properties, the finite difference method (FDM) can be used to calculate the flow temperature distribution [5]. In addition, recent works are focused on the finite element method (FEM) based solvers to simulate transient heat transfer in composite materials where the thermal conductivity is also dependent on the temperature [6], diffusion [7], and time-dependent value problems [8].
For macrosize systems in which the characteristic length is quite larger than the mean free path of gas, the flow is governed by the classical equations of Navier-Stokes and Fourier. However, at the micron level, a kinetic-based approach is needed to capture nonequilibrium effects described by Boltzmann equation [9]. In this context, Tang et al. [10] proposed a thermal boundary condition for TLBM. The results obtained are in good agreement with finite volume method (FVM) and kinetic solutions for Poiseuille and Couette flows. For convection problems, Zhang and Sun [11] have used the Chapman-Enskog expansion in the LBM method, to prove its consistency with macroscopic conservation laws. They simulated Rayleigh-Benard convection as a benchmark problem. Their results are interesting, reasonable, and meet well experimental data [12]. The stability of this scheme is also proved through different cases in a large range of Rayleigh numbers, until . Unlike the classical CFD approaches, regarding its kinetics nature, LBM describes the flow according to the Boltzmann equation. Therefore, the method ensures a link between the macro- and microscopic views of flow. Thus, LBM becomes a promising and powerful CFD tool as well as for microflow applications [13–15] usually reserved for pure kinetic approaches.
In the present work, the thermal and hydrodynamic behaviors of flow are modeled by means of thermal Lattice Boltzmann (TLBM) [2]. In this model, the temperature and density field are modeled independently using two distribution functions which allow us to simulate easily the heat transfer convection processes [3]. The LBM approach is used to study the classical benchmark problem of natural or free convection within a two-dimensional cavity, usually encountered as a basic geometry. In the most previous studies, of this problem, the wall (right or left) is fully heated. However, few studies treat a system with a partially heated (or cooled) wall which is the goal of this paper. Three heating ratios, length of a chip, are evaluated in this study for different Rayleigh numbers in the laminar regime range of flow. Note that the criterion for switching between the laminar () and turbulent regimes () are given by [16].
2. Problem Statement
A two-dimensional cavity enclosing a gas (air) is studied numerically using the lattice Boltzmann method. The enclosure considered has a long dimension in the -direction. Thus, a two-dimensional () and square, , cross-section geometry is studied. Natural convection flow is assured by the buoyancy effect, gravity, and thermal gradient imposed between the right wall kept at cold temperature , located at , and heated chip kept at in the left side located at (left case) or inversely (right case) (Figure 1), while the horizontal walls (located at and ) and apart from the chip are insulated (adiabatic).
3. Lattice Boltzmann Method
At the kinetic level, the gas flow behavior is governed by the Boltzmann equation given, in presence of external force per unit of mass , by where , and are, respectively, the time, particle velocity, particle position, velocity distribution function, and the collision operator. In the lattice Boltzmann method, the above equation is discretized and assumed it is valid along specific directions, linkages. Hence, the discrete Boltzmann equation can be written along a specified direction as,
Note that the collision operator represents the number of molecules rate of change in a controlled volume due to the particle collision.
In the Bhatnagar-Gross-Krook (BGK) approximation, it is related to the relaxation of distribution function towards equilibrium state by
In which represents the relaxation time of , to the equilibrium distribution function , related to the kinematic viscosity and the collision frequency by and , where is the velocity of sound related, in lattice, to the molecular velocity by
Using the finite difference scheme, the left side of the equation (2) can be written as, where . Using , the above equation becomes, where represents the buoyancy-driven force given in the Boussinesq approximation by [17]. where , , and are, respectively, the temperature at the reference state, the thermal expansion coefficient, and the weight factors given below.
The Maxwell distribution function can be written in the Taylor expansion as a function of density , macroscopic , and lattice velocities as
Note that for the discrete scheme D2Q9 (see Figure 2) used in this study, the lattices velocities and weights are given, respectively, by and .
Similarly, the following equation governs the internal energy distribution function variation: where the thermal relaxation time is related to the thermal diffusion coefficient by the relation . The first order approximation of is:
Using the following relations, the main macroscopic properties of flow such as density , velocity u, and internal energy by the unit of mass are calculated as where (two-dimensional flow) in which R is the universal gas constant and T is the gas temperature. In the following, we use the normalized variables:
Therefore, the hot temperature becomes , and the cold one will be . At the walls, no-slip velocity boundary conditions are used.
4. Results and Discussion
4.1. Validation Test
In this work, air natural convection flow within a partially heated cavity is investigated numerically using the LBM approach. A FORTRAN code is adopted to simulate the gas flow inside the enclosure. The classical problem of, fully heated wall, natural convection is used as a validation test case for the present code with the literature results. Figure 3 shows the temperature and the -velocity component profiles along the horizontal heated chip (left case) middle line. Good agreement is observed between the present simulation results and those obtained by Lai and Yang [12], on the one hand, and with the experimental data obtained by Krane and Jessee [18], on the other hand.
4.2. Mesh Size Independence
To test the grid size effect on the results, the average Nusselt number variation is calculated for different mesh sizes. Table 1 lists the average Nusselt number values for in the case of . It is shown that the grid size of is enough to give a good agreement with the experimental data [12] with an error less than 1%. To save the computation cost, time of calculation, and required memory, this mesh size was chosen for all the following simulations.
4.3. Flow Property Discussion
In this study, we evaluate the effect of Rayleigh number on the temperature and velocity fields. Figures 4 and 5 show streamlines overlaid on the temperature contours as a function of for the heated chip (left case) and cooled chip (right case) for , respectively. It is shown that the convection force , related to the Rayleigh number through the term (Eq. (16)), has a significant effect on both heat and viscous flow behaviors. At the low values of ( and ), the gas undergoes a circular motion of convection which leads to one primary vortex centered near the cavity center ().
(a)
(b)
(c)
(d)
(a)
(b)
(c)
(d)
According to the temperature of the walls, full or partially heated (or cooled), the temperature isotherm shape changes from the vertical walls forward the cavity bulk. This indicates that most of the heat transfer near the vertical walls is mainly assured by the conduction process. The single vortex observed is induced by the horizontal temperature gradient, , across the section. This gradient is always negative; therefore, clockwise motion is observed. Comparing both cases of hot and cold chips, it is shown that circulation of streamlines is not affected for the lower values of . However, by increasing value, the effect becomes more pronounced on both fields of velocity and temperature. The stronger convection () given by the warm-to-cold motion of flow gives birth to double stagnation points, secondary elliptic vortexes, in both cases. Their symmetry and shape break as the convection force increases (). The following table, Table 2, shows the locations of the vortex center for and . It is shown that the vortex center is affected by the Rayleigh number value. For the hot chip case, the vortex moves toward the left-top corner as increases. However, the stagnation point motion is toward the right-bottom corner in the cold chip case. This is due to the temperature distribution which is pushed under the convection force effect, the warm region at the top while the cold one is located at the bottom. Note that, for , the turbulence effect becomes more pronounced, and therefore, the flow behavior is more affected. This leads to a significant change in the vortex shape. Such behavior is not observed in the fully heated-wall case, [12]. Therefore, the chip length affects the flow characteristics. This effect is also illustrated in the heat transfer rate inside the cavity. Table 3 shows the average Nusselt number variation as a function of for different chip lengths.
The results show that the increases with and with the length of the chip. This agrees well with the definition which corresponds to the ratio of convection and conduction heat transfers (Eq. (16)). To illustrate more the flow circulation induced by the velocity field, the stream function () is evaluated for , , and , for both cases (see Figure 6). The plots show the sensitivity of stream function to the wall heated/cooled portion. According to the flow velocity, changes from a stronger, absolute, value near the vortex center () to lower values () at the walls. To show the characteristic variation across the cavity the figures (see Figures 7 and 8) show the normalized - and -components of velocity profiles along the vertical and horizontal center lines, respectively, as a function of for a hot and cold chip length of . It is shown that the velocity profile is sensitive to the Rayleigh number value. As shown in the 2D plots, the absolute value of velocity magnitude increases with . However, at the bulk of flow, the velocity is almost equal to zero, stagnation points.
(a)
(b)
(c)
(d)
(a)
(b)
(a)
(b)
By comparing profiles for both cases, it is shown, on the one hand, that as decreases the velocity magnitude undergoes a significant diminution from for to for near the left and right walls. On the second hand, in the cold chip case, this magnitude decreases with respect to the hot case, for . This effect of fast and slow vertical circulation in the hot and cold chip cases, respectively, induced by the temperature gradient vanishes far from the walls .
To more illustrate the stream function variation, the plots in Figure 8 show the vertical profiles of the velocity -component. In the hot chip case, a stronger horizontal motion () near the top (bottom) wall for hot (cold) chip case and the magnitude decreases with the decreasing. Comparing with the classical natural convection of fully heated walls, it is shown that the chip length has a significant effect on the viscous and thermal flow characteristics. Figures 9 and 10 plots show the velocity of -component and temperature profiles for different chip lengths. It is shown that the magnitude decreases with the chip length decreasing. Unlike the component, profiles have a similar shape in both cases (hot and cold chips). However, in the cold chip situation, the velocity has a lower magnitude compared with those of hot case. This means that the flow motion is also affected in the presence of partially cooled walls not only partially heated, usually investigated in the literature.
(a)
(b)
(a)
(b)
Besides the viscous flow behavior, the rate of heat transfer across the cavity is one of the most important characteristics targeted by engineers. To further clarify this effect, Figure 10 shows the temperature profiles along the horizontal center line for and different hot and cold chip lengths. It is shown that the temperature profile is also sensitive to the chip lengths. Comparing the profiles in both cases, the isotherm lines are similar and symmetric with respect to those of the fully heated wall situation for which the bulk temperature is almost equal to the average of the vertical wall temperatures. The temperature increases with the hot chip length and inversely with those of the cold one. The rate of heat convection with respect to the conduction process is estimated by means of local Nusselt number along with the left wall (); its profile is given by the following equation: where is the temperature at the bulk of flow which is calculated as the average temperature in the enclosure.
The average Nusselt number along the heated wall is given using the normalized parameters as where is the number of nodes along the vertical direction. Figure 11 shows the variation of as a function of for three chip lengths ( and ). The figure shows the relationship between heat transfer () and the Rayleigh number (). With the increasing of , the convection force increases according to Eq. (7) in which the thermal expansion coefficient is related to the Rayleigh number by the relation: where and are, respectively, the gravity acceleration and the temperature gap between the vertical walls. Note that the coefficients of viscosity and thermal diffusion are related to the Prandtl number by (air).
(a)
(b)
Along the chip, decreases until a constant value. At the hot chip extremities, discontinuous points, stronger values of are observed elsewhere corresponding to adiabatic wall portions. In the cold chip case, the value decreases significantly, and as in the full convection, undergoes two excitations, respectively, near the bottom () and top () walls. In addition, a supplementary discontinuous point is observed in the bulk. Its position changes according to the Rayleigh number and moves toward the bottom wall as increases, in agreement with the right secondary vortex observed in Figure 4.
5. Conclusion
In this work, the natural convection problem of gas flow within a square cavity partially heated or cooled cavity is studied numerically using the thermal lattice Boltzmann method (TLBM). According to the Rayleigh number, different hot chip lengths are compared with those of cold chip based on the profiles of temperature, velocity components, and Nusselt number. The results obtained in the full convection problem are promising and in good agreement with the experimental data of the literature. It is shown that the chip length has a significant effect on the thermal and hydrodynamic behaviors of flow. For higher value of , the flow stagnation point is pushed in the cold chip case toward the lower corner. However, in the hot chip one, the vortex is observed at the top left corner. In addition, the flow characteristics are strongly affected by the convection strength which is related to the Rayleigh number. It is shown that the rate of heated or cooled wall affects the heat transfer process and the viscous flow behavior with respect to the usually treated case of fully heated natural convection. To sum up, the flow behavior is sensitive to the heated or cooled rate of the wall.
Nomenclature and Symbols
| : | Aspect ratio, |
| : | Particle velocity |
| : | Sound velocity |
| : | Internal energy density |
| : | Velocity distribution function |
| : | Velocity and internal energy distribution functions in direction |
| : | External force, |
| : | Gravity acceleration vector, |
| : | Cavity height, |
| : | Chip length, |
| : | Heat transfer convection coefficient () |
| : | Thermal conductivity, () |
| : | Cavity width, |
| : | Local Nusselt number, |
| : | Average Nusselt number |
| : | Prandtl number, () |
| : | Rayleigh number, |
| : | Universal gases constant |
| : | Time, () |
| : | Temperature at the bulk of flow |
| : | Cold temperature, () |
| : | Hot temperature, () |
| : | Temperature, () |
| : | Macroscopic velocity vector |
| : | Macroscopic velocity components, () |
| : | Normalized velocity components |
| (x, y): | Spatial coordinates |
| : | Weight factors |
| (): | Normalized coordinates |
| : | Thermal diffusivity coefficient, () |
| : | Thermal expansion coefficient, () |
| : | Dynamic viscosity, () |
| : | Density, () |
| : | Kinematic viscosity, () |
| : | Gap of temperature, |
| : | Normalized temperature |
| : | Stream function |
| : | Particle collision operator |
| : | Relaxation times of distribution functions and to their corresponding ones and at the equilibrium state |
Data Availability
The data used to support the findings of this study are available from the corresponding author upon request.
Disclosure
This article is a revised and an expanded version of a paper presented at the International Conference on Intelligent Systems and Advanced Computing Sciences (ISACS’19), 26-27 December 2019, Taza-Morocco [19].
Conflicts of Interest
The authors declare that they have no conflicts of interest.