A Study on Natural Convection of Air in a Square Cavity with Partially Thermally Active Side Walls

Mame Khady Kane, Cheikh Mbow, Mamadou Lamine Sow, Joseph Sarr

Department of Physics, Cheikh Anta Diop University of Dakar, Dakar, Senegal.

DOI: 10.4236/ojfd.2017.74041                

Abstract

In this present work, we study heat transfer in a confined environment. We have to determine the thermal and dynamics fields of the cavity while observing the effect of the Rayleigh number which depends on the characteristics of the fluid and the temperatures imposed. The behavior of boundary layers in natural convection is analyzed along this square cavity. The central halves of its vertical walls are heated at different temperatures. The left active part is at a higher temperature than the one on the right wall. The remaining inactive parts and the horizontal walls (upper and lower) are adiabatic. The thermal and dynamic modeling of two-dimensional problem was done using a calculation code Fortran 90 and a visualization software ParaView based on the finite volume method. The equations governing this phenomenon of unsteady flow have thus been solved. This allows the modeling of both air flow and heat transfer with a numerical stabilization of the solution. So, we have presented our results of numerical simulations using a visualization tool. The results show the different velocity and temperature curves, velocity vectors and isotherms in laminar flow regime.

Keywords

Heat Transfer, Natural Convection, Square Cavity, Laminar Flow Regime

1. Introduction

Convection is the transfer of thermal energy due to the motion of fluids. Natural convection is defined as the fluid flow and heat transport arising due to density and temperature differences in the fluid. This is different from forced convection where the transport is because of an external energy source like a fan or a pump. Natural convection is an important topic of scientific research with applications in nature and engineering. For example in nature, a better understanding of natural convection will help us to predict the flow of magma in the earth’s core and mantle [1] . In engineering, natural convection is having important applications in cooling of electronics [2] , the design of efficient mixers and buildings [3] , nano-fluid devices [4] , etc.

There are mainly two configurations, the first being that of an enclosure containing a fluid and subjected to a vertical temperature gradient (for example Rayleigh-Benard convection), the second being that of a cavity with a horizontal temperature gradient. In this study, we focus on the case of a cavity with a horizontal temperature gradient. The case of vertical walls remains highly coveted in the industrial field because of its application to double skin facades, maintenance of equipment, electronic compounds, etc. Various geometric configurations have been studied theoretically, numerically or in the experimental design. The cavity studied has its thermally active vertical walls [5] at their centers. As a result, the temperature gradient is horizontal. It was interesting to conduct research in this field of study. Based on the approximations of Boussinesq [6] , we have developed a mathematical model which describes the problem and we were able to study numerically the heat exchange in a confined environment.

2. Materials and Methods

2.1. Basic Equations

The chosen geometry is a square cavity differentially heated with a side of 0.05 m. It is formed by two vertical walls, each of which is heated at its center. The left active part has a higher temperature $\left(Tr+\Delta T/2\right)$ than that of the right active part $\left(Tr-\Delta T/2\right)$ as shown in Figure 1. The boundary conditions are of Dirichlet type (imposed temperatures). The direction x is normal to the vertical walls, the gravity is carried by the axis y.

Figure 1. Description of the geometry of the model studied.

In Figure 1: s = 0.025 m, represents the height of each heated zone; A/4 = l = 0.0125 m―the height of the adiabatic zone; A = 0.05 m―the side of the square cavity; Tg and Td―the temperatures imposed respectively on the central halves of the left and right vertical walls. Convective heat transfer is a very complex phenomenon, and it depends on the boundary conditions (dynamic and thermal) of the system under consideration.

The dynamic boundary conditions are:

- to the walls, u = 0 m/s and v = 0 m/s.

The thermal boundary conditions give:

- for heated areas, we have:

for Ra = 1000 (Tg = 295.027 K and Td = 294.973 K),

for Ra = 10,000 (Tg = 295.27 K and Td = 294.73 K),

for Ra = 100,000 (Tg = 297.7 K and Td = 292.3 K),

for Ra = 1,000,000 (Tg = 322 K and Td = 268 K).

- for the adiabatic zone, φ = 0 W/m².

The flow is governed by the continuity equation, the equations of motion and the energy equation under the Boussinesq hypothesis, with: $U_r=\alpha /A$ for velocity, A for length, $A/U_r=A^2/\alpha$ for time and $\theta =\left(T-Tr\right)/\Delta T$ corresponds to the dimensionless temperature.

$\frac{\partial U}{\partial X}+\frac{\partial V}{\partial Y}=0$ (1)

$\frac{\partial U}{\partial \tau }+U\frac{\partial U}{\partial X}+V\frac{\partial U}{\partial Y}=-\frac{\partial P}{\partial X}+Pr\left(\frac{{\partial }^{2}U}{\partial X^2}+\frac{{\partial }^{2}U}{\partial Y^2}\right)$ (2)

$\frac{\partial V}{\partial \tau }+U\frac{\partial V}{\partial X}+V\frac{\partial V}{\partial Y}=-\frac{\partial P}{\partial Y}+Pr\left(\frac{{\partial }^{2}V}{\partial X^2}+\frac{{\partial }^{2}V}{\partial Y^2}\right)+Ra\cdot Pr\cdot \theta$ (3)

$\frac{\partial \theta }{\partial \tau }+U\frac{\partial \theta }{\partial X}+V\frac{\partial \theta }{\partial Y}=\frac{{\partial }^2\theta }{\partial X^2}+\frac{{\partial }^2\theta }{\partial Y^2}$ (4)

We have the following boundary conditions (hydrodynamic and thermal) at the walls:

$X=0$ , (Heated part of the left wall), $U=V=0$ , $\theta =1.$ (5)

$X=1$ , (Heated part of the right wall), $U=V=0$ , $\theta =0.$ (6)

$X=0$ and $X=1$ , (Adiabatic parts of left and right walls), $U=V=0.$ (7)

$\frac{\partial \theta }{\partial X}=0.$ (8)

$Y=0$ and $Y=1$ , (Lower and upper walls), $U=V=0.$ (9)

$\frac{\partial \theta }{\partial Y}=0.$ (10)

The dimensionless numbers used are:

The Prandtl number is physically defined as ratio of momentum diffusivity to thermal diffusivity. It provides a measure of the efficiency of diffusion transport through the velocity boundary layer and the thermal boundary layer.

$Pr=\frac{\nu }{\alpha }$ (11)

The Rayleigh number is physically the ratio of buoyancy and viscosity forces multiplied by the ratio of momentum and thermal diffusivities. Below the critical value of the Rayleigh number, heat transfer is primarily due to conduction. Above this critical value, heat transfer is due to convection.

$Ra=\frac{g\beta \Delta T{A}^3}{\nu \alpha }$ (12)

2.2. Validation of the Code

The results produced by our calculation code were compared with some results available in the literature. Let’s consider the same model as the benchmark [7] . Table 1 summarizes the maximum horizontal velocity on the vertical mid-plane of the cavity (together with its location) and the maximum vertical velocity on the horizontal mid-plane of the cavity (together with its location). We have also the maximum value of the local Nusselt number on the boundary at X = 0 (together with its location) and the minimum value of the local Nusselt number on the boundary at X = 0 (together with its location). They are compared with those obtained by different authors G. De Vahl Davis [7] (1983) and P. M. Gresho (1980) [8] for a Rayleigh number equal to 100,000. For this Rayleigh number, we used a mesh 81 × 81.

For the thermal and dynamic modeling of our problem, we first discretized the domain of study by a staggered cartesian mesh; the control volume of the component u is staggered in the direction x relative to the main control volume, that of the component v is staggered in the direction of y. Moreover, the staggered mesh technique consists in taking the vector and scalar variables in different nodes. This mesh is adapted to our choice of the distribution of the calculation points. Those of the temperature and the pressure are at the center of the meshes while those of the components u and v of the velocity are located on the faces of the control volume. We use the finite volume method which allows to

Table 1. Maximum velocities dimensionless and local Nusselt number.

$U_{\max }$ : the maximum horizontal velocity, $V_{\max }$ : the maximum vertical velocity, $N{u}_{\max }$ : the maximum value of the local Nusselt number and $N{u}_{\min }$ : the minimum value of the local Nusselt number.

model both the flow of a fluid and the transfer of heat. This method relies on a discretization technique, which converts the differential equations to the partial derivatives into nonlinear algebraic equations, which can then be solved numerically. It was described in 1971 by Patankar and Spalding and published in 1980 by Patankar [9] . Close to the active walls, the gradient of temperature and velocity is strong. For a better representation to the required accuracy of the reality of thermal transfers especially near the walls, the dimensions of the mesh are variable as shown in Figure 2 (41 × 41 for Ra = 1000, 10,000 to 81 × 81 (mesh finer)

Figure 2. Horizontal velocity profile for different meshes along the vertical centerline: (a) Ra = 1000, (b) Ra = 10000, (c) Ra = 100,000, (d) Ra = 1,000,000.

for Ra = 100,000, 1,000,000).

There is a change in the horizontal velocity curves from 41 × 41, 61 × 61 to 81 × 81. Hence considering the accuracy of the results required and computational time involved, a 41 × 41 grid size is chosen for computations with the two lower Rayleigh numbers. An 81 × 81 grid size is chosen also for the two upper Rayleigh numbers. The theoretical analysis of Rayleigh [10] was done at the beginning of the 20th century. Natural convective heat transfer is became an interesting field of study phenomenon [11] .

In the literature, some papers describe natural convection in rectangular enclosures [12] [13] [14] [15] [16] . However, in other studies, authors investigate the natural convection heat transfer in square cavities [17] . The study of heat transfer by natural convection in enclosures is an active field of research [18] and the fluid studied in the cavity may be air, water [19] [20] , etc.

2.3. Numerical Procedure

The differential Equations (1)-(4) governing the physical situation are translated into algebraic equations using the finite volume scheme. The system of dimensionless algebraic equations with boundary conditions associated (5)-(10) is iteratively solved. For accurate numerical simulation, a mesh size of 41 × 41 was selected for computations with Rayleigh numbers Ra = (1000, 10,000) and a mesh size of 81 × 81 for Ra = (100,000, 1,000,000). We have developed a numerical code with FORTRAN 90. The governing equations are discretized with the Central Difference Scheme (CDS). This scheme (order 2) generates a linear interpolation to the problem boundary. An iterative process is employed to find the velocity and temperature fields. The process is repeated until the following convergence criterion for velocity and temperature is met. The calculation stopped when the follows inequalities were satisfied:

$\text{Max}\left[\frac{{\displaystyle {\sum }_{i,j}\left|{\varphi }_{i,j}^{n+1}-{\varphi }_{i,j}^n\right|}}{{\displaystyle {\sum }_{i,j}\left|{\varphi }_{i,j}^{n+1}\right|}}\right]\le {10}^{-5}.$ (13)

The same convergence criterion is imposed in terms of relative error for velocity and temperature. In the above expression n is any time level and ϕ = U, V, θ. In this study, we used the finite volume method with quadrilateral control volumes and a staggered mesh. The latter is the subdivision of the field of study into longitudinal and transverse grids whose intersection represents a node, where the variables P and θ are located while the components U and V of the velocity vector are in the middle of the segments connecting two nodes adjacent. After discretization of the differential transport equations we obtain a system of nonlinear algebraic equations, these equations describe the discrete properties of the fluid at the nodes in the solution domain. For temporal discretization, we have used a numerical method of solving numerical differential equations, based on a multi-step method. The temporal integration was carried out on a mesh staggered according to an Adams-Bashforth scheme (order 2).

$\frac{\text{d}Y}{\text{d}\tau }=f\left(\tau ,Y\right)$ , with $\tau ={\tau }_0$ , $Y\left({\tau }_0\right)=Y_0.$ (14)

$\{\begin{array}{l}{Y}_0\left(\text{given}\right),\\ Y_1\text{approached}\left(one\text{-}stepmethod\right),\\ Y_{n+1}=Y_n+\frac{\Delta \tau }{2}\left(3f\left({\tau }_n,Y_n\right)-f\left({\tau }_{n-1},Y_{n-1}\right)\right).\end{array}$ (15)

3. Results and Discussion

We performed our simulations by varying the Rayleigh number. We have analyzed below the velocity and temperature curves in the cavity for four different cases. Natural convection heat transfer is studied in a square enclosure for different thermally active locations with a Prandtl number Pr = 0.72 (air).

3.1. For Ra = 1000

In the first case, we have the following solutions for different thermally active locations. Figure 3 shows the profile for the different locations of the thermally active walls. It is observed that the fluid particle moves with greater velocity for the middle active locations and the velocity is less for the top/bottom active locations. Figure 3 shows that the variation of the vertical velocity reaches a maximum in the vicinity of the left heated wall and decreases as one moves away in order to be canceled towards x = 0.025 m. It descends by forming a kind of central symmetry and reaches a minimum in the vicinity of the right heated wall. The evolution of the vertical velocity in the cavity along the ordinate shows that the amplitude of this velocity is greater when the height is between the active

Figure 3. Vertical velocity profile along the ordinate for Ra = 1000.

parts than on those which are inactive (adiabatic). In Figure 4, we have represented the temperature profiles along the ordinate. Temperature decreases almost linearly with abscissa. Temperatures increase in the vicinity of the heated left part and decrease near the active right part.

3.2. For Ra = 10,000

In a second case, we have the following results for different thermally active locations. Figure 5 shows the profile for the different locations of the thermally

Figure 4. Temperature profile along the ordinate for Ra = 1000.

Figure 5. Vertical velocity profile along the ordinate for Ra = 10,000.

active walls. It is observed that the fluid particle moves with greater velocity for the middle active locations and the velocity is less for the top/bottom active locations. With a larger temperature variation and therefore a higher Rayleigh number, we observe in Figure 5 that the amplitudes of velocities increased in the vicinity of the heated left wall by being canceled towards x = 0.025 m and decreasing near the active right wall. The evolution of the vertical velocity in the cavity along the ordinate has become greater than in the preceding case, with the velocity amplitude being greater when the height approaches the center of the active part (y = 0.025 m). In Figure 6, temperature decreases with abscissa and the convection is accentuated because we have increased the variation of temperatures between active parts. Temperatures increase again in the vicinity of the heated left part and decrease near the active right part.

3.3. For Ra = 100,000

In the third case, we have the following findings for different thermally active locations. Figure 7 shows the profile for the different locations of the thermally active walls. It is observed that the fluid particle moves with greater velocity for the middle active locations and the velocity is less for the top/bottom active locations. For an even larger Rayleigh number, the flow regime remains laminar with greater velocity amplitudes in Figure 7. Temperatures have decreased again with abscissa in Figure 8.

3.4. For Ra = 1,000,000

In the fourth case, we have the following results for different thermally active locations. Figure 9 shows the profile for the different locations of the thermally

Figure 6. Temperature profile along the ordinate for Ra = 10,000.

Figure 7. Vertical velocity profile along the ordinate for Ra = 100,000.

Figure 8. Temperature profile along the ordinate for Ra = 100,000.

active walls. It is observed that the fluid particle moves with greater velocity for the middle active locations and the velocity is less for the top/bottom active locations. The flow regime always remains laminar with higher velocity and temperature variations (Figure 9 and Figure 10).

Figure 9. Vertical velocity profile along the ordinate for Ra = 1,000,000.

Figure 10. Temperature profile along the ordinate for Ra = 1,000,000.

3.5. Effects of Rayleigh Number

Figures 11-14 show the impact of variation in the Rayleigh number. With an increase in the Rayleigh number from 1000 to 1,000,000, we have some noticeable changes. There are fluctuations in the velocity profiles and in particular for

Figure 11. Variation of horizontal velocity at mid-cavity (vertical line) for different Rayleigh numbers.

Figure 12. Variation of vertical velocity at mid-cavity (horizontal line) for different Rayleigh numbers.

the maximum values u_max and v_max. The velocity components (Figure 11 and Figure 12) increase progressively in amplitude and the maximum of v_max moves to the vertical walls with the increase of the Rayleigh number. The centro-symmetry’s property is retained throughout the range of the Rayleigh number studied. We observe a symmetry of the curves of velocity with respect to the

Figure 13. Variation of temperature at mid-cavity (vertical line) for different Rayleigh numbers.

Figure 14. Variation of temperature at mid-cavity (horizontal line) for different Rayleigh numbers.

center of the cavity. For our model, the variation of the Rayleigh number represents a variation proportional to the temperature difference between the right and left active parts. In the zone 0.0125 m ≤ y ≤ 0.0375 m, the temperature increases progressively (Figure 13). This is where the heated parts are located. Fluid particles move with greater velocity due to an increase in buoyancy force for large values of the Rayleigh number. On the other hand, in the zone 0.0125 m ≤ x ≤ 0.375 m, the temperature remains almost constant (Figure 14). For a small Rayleigh number, of the order of 1000, we have noticed the dominance of the heat transfer mode by conduction. Beyond this value, convection dominates and appears more clearly for Ra = 1,000,000. Temperature and velocity are influenced by the increasing the Rayleigh number, which increases them at the same time. The local temperature profile certifies a centro-symmetry property of the local temperature field. The determination of these parameters allowed us to analyze the behavior of the fluid inside the cavity.

3.6. The Velocity Vector Maps and the Isotherms

The results are displayed graphically in terms of velocity vectors and isotherms. In Figure 15, we show the velocity vectors corresponding to the end of the simulation (final time). The circulation of the fluid inside the cavity is shown by the velocity vectors. In Figure 15(a), we can observe the existence of a stagnant zone in the center of the cavity, which means that the heat exchange takes place

(a) (b) (c) (d)

Figure 15. Velocity vector maps: (a) [Ra = 1000; contours at (0 - 0.001166) (m/s)]; (b) [Ra = 10,000; contours at (0 - 0.006481) (m/s)]; (c) [Ra = 100,000; contours at (0 - 0.026042) (m/s)]; (d) [Ra = 1,000,000; contours at (0 - 0.087263) (m/s)].

in an intense way at the corners of the cavity. In Figure 15(b), the stagnant zone at the center is always present inside the cavity. In the case of a low Rayleigh number, the fluid rotates around the cavity midpoint. In Figure 15(c) & Figure 15(d), we noticed the existence of two cells of recirculation inside the cavity, one to the upper left side and the other to the lower right side. The two remaining corners are less active. The flow leaving an active wall region is free to move either upwards or downwards. The formation of the two recirculation cells shows that the solution has a symmetrical flow. Fluid circulation therefore fills the entire cavity even at a very high Rayleigh number.

In Figure 16(a), the isotherms are almost parallel to the active parts. It is indicating that only the conduction mode of heat transfer prevails in this stage. For a low value of the Rayleigh number equal to 1000, the dominance of the conduction heat transfer mode was noted. In Figure 16(a) & Figure 16(b), isotherms show a monocellular flow represented by a sort of vortex rotating in the sense of a clock. The flow of the fluid takes place in the form of a large recirculation loop.

(a) (b) (c) (d)

Figure 16. Contour maps of temperature T: (a) [Ra = 1000; contours at (294.9738 - 295.0262) (K)]; (b) [Ra = 10,000; contours at (294.7406 - 295.2587) (K)]; (c) [Ra = 100,000; contours at (292.3952 - 297.5999) (K)]; (d) [Ra = 1,000,000; contours at (269.6427 - 320.2993) (K)].

As it approaches the hot wall, the particle gains heat (its kinetic energy), heats up (becomes lighter) and rises along the hot wall with a trajectory that will be further modified by the presence of the upper horizontal wall. Along the latter, the fluid particle cools as it approaches the cold wall: we have the formation of a circulation cell. The warm air rises and the cold air descends. In Figure 16(c) & Figure 16(d), we notice that the isotherms are close to one another and condense near the active situations, where the temperature gradients are high, whereas they are negligible for the rest of the walls of the cavity. This reflects the existence of the thermally strong boundary layers close to these zones. We show that there is the phenomenon of convection in the vicinity of the active places, in particular for a high Rayleigh number, of the order of 100,000. Isotherms are crowded around the active locations on the left and right sides of the enclosure. The heat transfer is high because this situation gives the minimum distance between points on the hot and cold surfaces.

4. Conclusions

We have studied the behavior of boundary layers in natural convection along a cavity. The central parts of the vertical walls are at imposed temperatures (horizontal gradient). We have seen that from Ra = 1000 to Ra = 1,000,000, the flow regime is laminar and evolves remarkably. Analysis of the results obtained with different ordinates showed that for different heights the values of velocity and temperature fluctuated due to the convection which creates upward and downward movements of the air. Moreover, we have shown that there are recirculations and a property of centro-symmetry in the flow. The vertical velocity profiles show a correspondence with the results obtained in the form of velocity vectors.

We know that in physics, particularly in solar energy, for collectors due to shading it is only the unshaded part of the wall that is thermally active. Buoyancy-driven flow in a square cavity with vertical sides which are differentially and partially heated is a very important way to understand the thermal and dynamic effects in many practical problems.

Acknowledgements

We thank all the contributing members of SOLMATS group for their support, their assistance and their encouragement.

Competing Interests

We author(s) declare that we have no competing interests.

Nomenclature

A: Side of the cavity m

g: Gravity intensity m/s²

l: Length of adiabatic parts m

P: Dimensionless pressure

Pr: Prandtl Number

Ra: Rayleigh Number

s: Length of active parts m

T: Dimensional temperature K

Td: Dimensional temperature of right active part K

Tg: Dimensional temperature of left active part K

Tr: Reference temperature K

u: Dimensional horizontal velocity m/s

U: Dimensionless horizontal velocity

U_r: Reference velocity m/s

v: Dimensional vertical velocity m/s

V: Dimensionless vertical velocity

x: Dimensional abscissa m

X: Dimensionless abscissa

y: Dimensional ordinate m

Y: Dimensionless ordinate

Greek Letters

$\alpha$ : Thermal diffusivity m²/s

$\beta$ : Coefficient of thermal expansion 1/K

$\Delta$ : Variation

$\theta$ : Dimensionless temperature

$\nu$ : Kinematic viscosity m²/s

$\tau$ : Dimensionless time

$\varphi$ : Dimensionless variable

$\phi$ : Heat flux density W/m²

Subscript

Max or max Maximum value

min Minimum value

Conflicts of Interest

The authors declare no conflicts of interest.

References

[1]	Mark Jellinek, A., Kerr, R.C. and Griffiths, R.W. (1999) Mixing and Compositional Stratification Produced by Natural Convection 1. Experiments and their Application to Earth’s Core and Mantle. Journal of Geophysical Research, 104, 7183-7201. https://doi.org/10.1029/1998JB900117 https://doi.org/10.1029/1998JB900116
[2]	Florio, L.A. and Harnoy, A. (2007) Combination Technique for Improving Natural Convection Cooling in Electronics. International Journal of Thermal Sciences, 46, 76-92. https://doi.org/10.1016/j.ijthermalsci.2006.03.007
[3]	Mishra, A.A., Hasan, N., Sanghi, S. and Kumar, R. (2008) Two-Dimensional Buoyancy Driven Thermal Mixing in a Horizontally Partitioned Adiabatic Enclosure. Physics of Fluids, 20, 063601. https://doi.org/10.1063/1.2931567
[4]	Putra, N., Roetzel, W. and Das, S.K. (2003) Natural Convection of Nano-Fluids. Heat and Mass Transfer, 39, 775-784. https://doi.org/10.1007/s00231-002-0382-z
[5]	Valencia, A. and Frederick, R.L. (1989) Heat Transfer in Square Cavities with Partially Active Vertical Walls. International Journal of Heat and Mass Transfer, 32, 1567-1574. https://doi.org/10.1016/0017-9310(89)90078-1
[6]	Bae, J.H. and Hyun, J.M. (2004) Time-Dependent Buoyant Convection in an En-closure with Discrete Heat Sources. International Journal of Thermal Sciences, 43, 3-11. https://doi.org/10.1016/S1290-0729(03)00102-9
[7]	De Vahl Davis, G. (1983) Natural Convection of Air in a Square Cavity: a Bench Mark Numerical Solution. International Journal for Numerical Methods in Fluids, 3, 249-264. https://doi.org/10.1002/fld.1650030305 https://doi.org/10.1002/fld.1650030304
[8]	Gresho, P.M., Upson, C.D. and Lee, R.L. (1980) Finite Element Simulations of Thermally Induced Convection in an Enclosed Cavity. Lawrence Livermore Laboratory Report UCID-18602.
[9]	Patankar, S.V. (1980) Numerical Heat Transfer and Fluid Flow. Hemisphere Publishing Corporation, Taylor and Francis Group, New York.
[10]	Rayleigh, L. (1916) On Convection Currents in a Horizontal Layer of Fluid When the Higher Temperature is on the Under Side. Philosophical Magazine, 32, 529-546. https://doi.org/10.1080/14786441608635602
[11]	Calcagni, B., Marsili, F. and Paroncini, M. (2005) Natural Convective Heat Transfer in Square Enclosures Heated From Below. Applied Thermal Engineering, 25, 2522-2531. https://doi.org/10.1016/j.applthermaleng.2004.11.032
[12]	November, M. and Nansteel, M.W. (1987) Natural Convection in Rectangular En-closures Heated from Below and Cooled along One Side. International Journal of Heat and Mass Transfer, 30, 2433-2440. https://doi.org/10.1016/0017-9310(87)90233-X
[13]	Ganzarolli, M.M. and Milanez, L.F. (1995) Natural Convection in Rectangular Enclosures Heated from Below and Symmetrically Cooled from the Sides. International Journal of Heat and Mass Transfer, 38, 1063-1073. https://doi.org/10.1016/0017-9310(94)00217-J
[14]	Hasnaoui, M., Bilgen, E. and Vasseur, P. (1992) Natural Convection Heat Transfer in Rectangular Cavities Partially Heated from Below. Journal of Thermophysics and Heat Transfer, 6, 255-264. https://doi.org/10.2514/3.353
[15]	Ho, C.J. and Chang, J.Y. (1994) A Study of Natural Convection Heat Transfer in a Vertical Rectangular Enclosure with Two-Dimensional Discrete Heating: Effect of Aspect Ratio. International Journal of Heat and Mass Transfer, 37, 917-925. https://doi.org/10.1016/0017-9310(94)90217-8
[16]	Deng, Q.-H., Tang, G.-F. and Li, Y. (2002) A Combined Temperature Scale for Analyzing Natural Convection in Rectangular Enclosures with Discrete Wall Heat Sources. International Journal of Heat and Mass Transfer, 45, 3437-3446. https://doi.org/10.1016/S0017-9310(02)00060-1
[17]	Nithyadevi, N., Kandaswamy, P. and Sivasankaran, S. (2006) Natural Convection in a Square Cavity with Partially Active Vertical Walls: Time-Periodic Boundary Condition. Mathematical Problems in Engineering, 1-16. https://doi.org/10.1155/MPE/2006/23425
[18]	Jaluria, Y. (1980) Natural Convection Heat and Mass Transfer. Pergamon Press, Oxford.
[19]	Tong, W. (1999) Aspect Ratio Effect on Natural Convection in Water near Its Density Maximum Temperature. International Journal of Heat and Fluid Flow, 20, 624-633. https://doi.org/10.1016/S0142-727X(99)00027-2
[20]	Kandaswamy, P., Sivasankaran, S. and Nithyadevi, N. (2007) Buoyancy-Driven Convection of Water near Its Density Maximum with Partially Active Vertical Walls. International Journal of Heat and Mass Transfer, 50, 942-948. https://doi.org/10.1016/j.ijheatmasstransfer.2006.08.013