Simulation of Thermal Radiation and Turbulent Free Convection in an Enclosure with a Glass Wall and a Local Heater

Abstract

In this study, a numerical modelling of thermal radiation and turbulent thermogravitational convection in a large-scale chamber containing a thermally-generating element is conducted. The lower border of the cabinet is maintained under adiabatic conditions, while on the other walls the convective boundary conditions (Robin boundary condition) are used. The managing equations with corresponding restrictions are transformed using the stream function–vorticity formulation and then solved by employing a finite difference method. The influence of both the height and wall emissivity of the heated source on fluid motion and the heat transmission in a large-scale chamber is investigated. Our results of the calculations on non-uniform grids with algebraic transformation are in excellent agreement with other available experimental and numerical outcomes for turbulent thermal convection in enclosures. The computations indicate that the average total Nusselt number is enhanced up to 2 times with an increase in the heater height. The results show that the surface emissivity of the heat source has a great influence on the total thermal transference coefficient. Furthermore, a growth of the heater surface emissivity has no significant effect on the flow structure.

1. Introduction

The study of convective–radiative heat transfer is extremely important for various fields of modern science and technology, such as mechanical engineering, microelectronics, building design, and nuclear technologies. The significant interest in this type of research is due, first of all, to its primary importance for understanding the physical processes taking place in modern energy and in technological structures of various kinds. In this regard, the combination of thermogravitational convection and surface radiation in chambers is the subject of many theoretical and experimental works, the number of which is rapidly increasing due to the development of methods of mathematical modeling and improvement of computing technology, along with improved visualization and methods for measuring experimental data.

In many works devoted to convective heat transfer in closed cavities, the influence of the radiative mechanism of energy transfer has not been taken into account due to the difficulties associated with modeling radiation or due to its insignificant influence on the circulation structure and energy transfer in some engineering applications [1,2,3,4,5]. For many gases (at low and moderate temperatures), such as nitrogen, hydrogen, ozone, oxygen, and helium, their own radiation is quite small, and the gases are essentially transparent to radiation. For example, for turbulent free convection in a cubic chamber, it has been shown earlier [6] that gas radiation has an extremely insignificant influence on the mass and energy transport inside the considered volume. In particular, the average total Nusselt number changes by only 1.5%. However, for the considered class of problems, it was found that one should not neglect the radiation between the bounding surfaces. If the temperature of such a surface is high enough or if the fluid flow occurs at a relatively low velocity (for example, thermogravitational convection), then the contribution of the radiation component can be quite significant.

In the modern construction industry, the problem of heat and energy conservation is very relevant. Minimizing energy losses allows us not only to make housing more comfortable and of high quality for living, but it also leads to enormous savings. Mathematical modeling can help to solve the problem of improving energy efficiency, as well as calculate the operational characteristics of buildings depending on the ambient temperature and the required indoor climate. For a correct description of heat transfer processes in building structures, it is necessary to take into account their unsteady and turbulent nature, as well as the influence of external media and internal heat sources. Wang et al. [7] studied the impact of convection with radiation interaction in an industrial building with a high-temperature heater, where the high-temperature heater was located in the central part of the floor. The obtained results were reported for various values of Grashof number (Gr) (from 10⁷ to 10¹¹) and surface emissivity (0.2–0.8). They found that the ratio between convection and radiation is a raising function of Gr and a diminution function of the value of wall emissivity. Rahimi and Sabernaeemi [8] experimentally investigated the radiant ceiling heating system in a large-scale room. The main goal of their work was to study the proportion of heat transmitted by radiation. Their results clearly demonstrated that only 10% of the heat was transferred by the convection (or conduction). Ibrahim et al. [9] reported detailed studies of the impact of radiation on thermogravitational convection in enclosures. The obtained results were reported for a fixed value of Rayleigh number, Ra = 1.5 × 10⁹. The radiation model took into account the radiation of the gas and walls. Four cases were considered, namely, combined gas and wall radiation, neglecting radiation, surface radiation only, and gas radiation only. It was shown that gas radiation has a minor influence on the circulation nature.

A detailed experimental analysis of the turbulent modes of free convection in a closed parallelepiped having differentially warmed side borders was carried out by Tian and Karayiannis [10,11], who described the same experiment but were devoted to different problems. In the first study [10], the authors obtained the distributions of the average temperature and average velocity, and they also determined the circulation nature and the thickness of the temperature and velocity boundary layers. In the second part [11], the main attention was focused on the study of turbulent quantities (pulsation values of temperature and velocity, as well as Reynolds stresses). The evolution of a thermal plume in a closed cavity near an isothermal wall by the digital tracer visualization method was carried out by Caudwell et al. [12]. A year earlier [13], a group of Chinese colleagues from Tsinghua University conducted a complicated study of turbulent free convection in a rectangular region with dimensions of 1960 × 980 × 800 mm³ with a local source of energy release with the characteristic dimensions of 200 × 300 mm². Air at normal pressure was used as a working medium, and the Rayleigh number was equal to 10⁹. It was shown that the results of numerical simulation using the standard k–ε (SKE) model and the renormalization k–ε (RNGKE) turbulence model are in good agreement with the experimental outcomes. The obtained results can help us to better understand the flow behavior.

A numerical estimate of the efficiency of six semiempirical turbulence models as applied to problems of free convection in closed 2D and 3D domains was carried out by Altac and Ugurlubilek [14]. The authors considered the following RANS models: a standard k-ω, a Reynolds stress model, a standard k-ε, shear stress transport k-ω, realizable k-ε, and renormalization group k-ε. The Rayleigh number varied from 10⁸ to 10¹³. The obtained simulation results were compared in detail with more than fifteen experimental and numerical works of other authors. As a result of the studies, the possibility of using the data of the two-dimensional approximation to describe the processes of heat transfer in spatial objects at Rayleigh numbers up to 10¹⁰ was established.

Despite the constantly increasing computing power, the numerical modeling of turbulent thermogravitational convection combined with radiation is still a difficult task from the point of view of numerical implementation. To perform mathematical simulations within a reasonable amount of time, a number of assumptions are often used. The most common assumptions are the following: Radiation does not depend on the wavelength and is diffuse, that is, the radiation intensity at any point on the surface is evenly distributed in all directions. However, regardless of the assumptions used, the radiation intensity is still dependent on the geometrical parameters and the thermal attributes of the interacting walls. A detailed review of different analytical and numerical techniques with their strengths and weaknesses was prepared by Mishra and Prasad [15].

Shati et al. [16] reported a computational study of the turbulent modes of free convection and surface radiation in a closed chamber that has differentially warmed vertical borders. In the framework of computational experiments, the aspect ratio varied between 0.0625 and 16; for the square chamber, the length of the sides varied from 40 to 240 cm, and the temperature ratio of the cold T_c and hot T_h walls varied from 1.02 to 2.61. The authors obtained correlations for calculating the average Nu as a function of Gr and Pr. It was also found that the ratio between the strength of energy transference owing to convection and the strength of energy transference owing to radiation increases with an increase in $T_h/T_c$ from 1.02 to 1.2, and decreases at $T_h/T_c>1.2$. Later, the authors experimentally confirmed the obtained numerical results by performing a series of experiments [17] under conditions similar to [16]. Convective–conductive–radiative energy transference in a closed cube with a heater along the lower wall was investigated by Yang and Wu [18]. Correlations were obtained for the average Nusselt numbers on the surface of the energy source, and a multivariable analysis of the effect of emissivity and thermal conductivity coefficient on solid walls for heat transfer and flow structure was carried out. Moreover, it was found that by increasing the aspect ratio of the cubic cavity, one can significantly increase the intensity of convective energy transference. An extensive bibliography on turbulent free convection under a possible influence of radiation in enclosures can be found in the review prepared by Miroshnichenko and Sheremet [19]. An analytical solution of two-dimensional Navier–Stokes equation with a non-Newtonian type of viscosity was carried out by Barna et al. [20]. They concluded that the velocity field of the fluid—in contrast to a Newtonian type of viscosity—has a compact support. Bognar and Hriczo [21] investigated thermal Marangoni convection; in this work the influence of the power law exponent and Prandtl number on the velocity and temperature profiles were studied.

Previous research demonstrated that a large-scale enclosure heated by a local heater with a fixed temperature has been widely investigated. The results presented in this work correspond to the conditions of constant heat generation inside the power source. These conditions are of maximum interest from the point of view of practical applications. This is due to the fact that in practice, publicly available information is precisely the power of the device, not the various kinds of temperature dependence on the walls of the warmed elements. The main goal of this research is to investigate the influence of both height and wall emissivity of the heat source on the fluid circulation and energy transference in a large-scale chamber.

2. Governing Equations and Numerical Method

Figure 1 shows the system configuration, including the coordinate system and computational domain. The system consisted of a square enclosure with massive walls, one of which is glass, and a local heat-generating source. The air is inside the enclosure (with the Prandtl number being 0.71). To some extent, the system configuration can be viewed as a small room with a panoramic window. It was assumed that the outer edge of the bottom wall was ideally insulated. Heat transfer with the environment was considered on other outside surfaces of solid borders of finite thickness. The local heater of height l was situated near the glass wall of the large-scale enclosure and had an inner volumetric heat generation. The inner surfaces of the enclosure were considered to be both gray and opaque emitters. The flow inside the enclosure was considered turbulent. To describe the change in density with temperature, the Boussinesq approximation was used. In the present work, a two-dimensional formulation was considered (the middle cross section of the enclosure along the axis Z).

The set of unsteady Reynolds-averaged Navier–Stokes equations (URANS) was applied for a mathematical statement of the problem of convective–radiative heat transfer inside the enclosure [22]. The use of URANS model is due to the rational requirements of this model for computing resources, as well as its ability to provide quite acceptable accuracy in calculating various types of turbulent flows.

Taking into account the above assumptions, for the mathematical description of the considered problem in the analyzed domain (Figure 1), the governing equations for the conservation of mass, momentum, and energy in primitive variables were defined in the following form:

The equation for the air cavity:

$$\frac{\partial u}{\partial x}+\frac{\partial v}{\partial y}=0$$

(1)

$$\frac{\partial u}{\partial t}+u\frac{\partial u}{\partial x}+v\frac{\partial u}{\partial y}=-\frac{1}{\mathsf{\rho }}\frac{\partial p}{\partial x}+2\frac{\partial }{\partial x}\left[\left(\mathsf{\nu }+{\mathsf{\nu }}_t\right)\frac{\partial u}{\partial x}\right]+\frac{\partial }{\partial y}\left[\left(\mathsf{\nu }+{\mathsf{\nu }}_t\right)\left(\frac{\partial u}{\partial y}+\frac{\partial v}{\partial x}\right)\right]$$

(2)

$$\frac{\partial v}{\partial t}+u\frac{\partial v}{\partial x}+v\frac{\partial v}{\partial y}=-\frac{1}{\mathsf{\rho }}\frac{\partial p}{\partial y}+\frac{\partial }{\partial x}\left[\left(\mathsf{\nu }+{\mathsf{\nu }}_t\right)\left(\frac{\partial u}{\partial y}+\frac{\partial v}{\partial x}\right)\right]+2\frac{\partial }{\partial y}\left[\left(\mathsf{\nu }+{\mathsf{\nu }}_t\right)\frac{\partial v}{\partial y}\right]+g\mathsf{\beta }\left(T-T_0\right)$$

(3)

$$\frac{\partial T}{\partial t}+u\frac{\partial T}{\partial x}+v\frac{\partial T}{\partial y}=\frac{\partial }{\partial x}\left[\left({\mathsf{\alpha }}_{air}+\frac{{\mathsf{\nu }}_t}{{\mathrm{Pr}}_t}\right)\frac{\partial T}{\partial x}\right]+\frac{\partial }{\partial y}\left[\left({\mathsf{\alpha }}_{air}+\frac{{\mathsf{\nu }}_t}{{\mathrm{Pr}}_t}\right)\frac{\partial T}{\partial y}\right]$$

(4)

$$\frac{\partial k}{\partial t}+u\frac{\partial k}{\partial x}+v\frac{\partial k}{\partial y}=\frac{\partial }{\partial x}\left[\left(\mathsf{\nu }+\frac{{\mathsf{\nu }}_t}{{\sigma }_k}\right)\frac{\partial k}{\partial x}\right]+\frac{\partial }{\partial y}\left[\left(\mathsf{\nu }+\frac{{\mathsf{\nu }}_t}{{\sigma }_k}\right)\frac{\partial k}{\partial y}\right]+P_k+G_k-\mathsf{\epsilon }$$

(5)

$$\frac{\partial \mathsf{\epsilon }}{\partial t}+u\frac{\partial \mathsf{\epsilon }}{\partial x}+v\frac{\partial \mathsf{\epsilon }}{\partial y}=\frac{\partial }{\partial x}\left[\left(\mathsf{\nu }+\frac{{\mathsf{\nu }}_t}{{\sigma }_{\mathsf{\epsilon }}}\right)\frac{\partial \mathsf{\epsilon }}{\partial x}\right]+\frac{\partial }{\partial y}\left[\left(\mathsf{\nu }+\frac{{\mathsf{\nu }}_t}{{\sigma }_{\mathsf{\epsilon }}}\right)\frac{\partial \mathsf{\epsilon }}{\partial y}\right]+\left(c_{1\mathsf{\epsilon }}\left(P_k+c_{3\mathsf{\epsilon }}{G}_k\right)-c_{2\mathsf{\epsilon }}\mathsf{\epsilon }\right)\frac{\mathsf{\epsilon }}{k}$$

(6)

The equation for the solid walls:

$$\frac{\partial T}{\partial t}={\mathsf{\alpha }}_w\left(\frac{{\partial }^{2}T}{\partial x^2}+\frac{{\partial }^{2}T}{\partial y^2}\right)$$

(7)

The equation for the glass door:

$$\frac{\partial T}{\partial t}={\mathsf{\alpha }}_g\left(\frac{{\partial }^{2}T}{\partial x^2}+\frac{{\partial }^{2}T}{\partial y^2}\right)$$

(8)

The equation for the heat-generating source:

$$\frac{\partial T}{\partial t}={\mathsf{\alpha }}_{hs}\left(\frac{{\partial }^{2}T}{\partial x^2}+\frac{{\partial }^{2}T}{\partial y^2}\right)+\frac{q_v}{\mathsf{\rho }{c}_p}.$$

(9)

The air flow inside the enclosure was examined in terms of streamlines. In this connection, the vorticity $\left(\omega =\frac{\partial v}{\partial x}-\frac{\partial u}{\partial y}\right)$ and stream function $\left(u=\frac{\partial \mathsf{\psi }}{\partial y}, v=-\frac{\partial \mathsf{\psi }}{\partial x}\right)$ variables were used.

The coordinate transformation is used to thicken the computational mesh near heater and solid walls (Figure 2). Transformation of coordinate is considered in the following form [23]:

$$\begin{array}{l}\xi =a+\frac{b-a}{2}\left\{1+\mathrm{tg}\left[\frac{\pi \kappa }{b-a}\left(x-\frac{a+b}{2}\right)\right]/\mathrm{tg}\left[\frac{\pi }{2}\kappa \right]\right\}, \\ \mathsf{\eta }=a+\frac{b-a}{2}\left\{1+\mathrm{tg}\left[\frac{\pi \kappa }{b-a}\left(y-\frac{a+b}{2}\right)\right]/\mathrm{tg}\left[\frac{\pi }{2}\kappa \right]\right\}.\end{array}$$

where κ is a compaction parameter, and a, b are the geometrical characteristics.

Taking into account the algebraic coordinate transformation noted above, the transport equations in the non-dimensional non-primitive variables of stream function–vorticity can be shown in the following form:

Inside the enclosure:

$$\frac{d^2\xi }{d{X}^2}\frac{\partial \Psi }{\partial \xi }+{\left(\frac{d\xi }{dX}\right)}^2\frac{{\partial }^2\Psi }{\partial {\xi }^2}+\frac{d^2\mathsf{\eta }}{d{Y}^2}\frac{\partial \Psi }{\partial \mathsf{\eta }}+{\left(\frac{d\mathsf{\eta }}{dY}\right)}^2\frac{{\partial }^2\Psi }{\partial {\mathsf{\eta }}^2}=-\Omega ,$$

(10)

$$\begin{array}{l}\frac{\partial \Omega }{\partial \tau }+\left(U-\frac{d\xi }{dX}\frac{\partial {\mathsf{\nu }}_t}{\partial \xi }\right)\frac{d\xi }{dX}\frac{\partial \Omega }{\partial \xi }+\left(V-\frac{d\mathsf{\eta }}{dY}\frac{\partial {\mathsf{\nu }}_t}{\partial \mathsf{\eta }}\right)\frac{d\mathsf{\eta }}{dY}\frac{\partial \Omega }{\partial \mathsf{\eta }}=\frac{d\xi }{dX}\frac{\partial }{\partial \xi }\left[\left(\sqrt{\frac{\mathrm{Pr}}{\mathrm{Ra}}}+{\mathsf{\nu }}_t\right)\frac{d\xi }{dX}\frac{\partial \Omega }{\partial \xi }\right]+\\ +\frac{d\mathsf{\eta }}{dY}\frac{\partial }{\partial \mathsf{\eta }}\left[\left(\sqrt{\frac{\mathrm{Pr}}{\mathrm{Ra}}}+{\mathsf{\nu }}_t\right)\frac{d\mathsf{\eta }}{dY}\frac{\partial \Omega }{\partial \mathsf{\eta }}\right]+\left(\frac{d^2\xi }{d{X}^2}\frac{\partial {\mathsf{\nu }}_t}{\partial \xi }+{\left(\frac{d\xi }{dX}\right)}^2\frac{{\partial }^2{\mathsf{\nu }}_t}{\partial {\xi }^2}-\frac{d^2\mathsf{\eta }}{d{Y}^2}\frac{\partial {\mathsf{\nu }}_t}{\partial \mathsf{\eta }}-{\left(\frac{d\mathsf{\eta }}{dY}\right)}^2\frac{{\partial }^2{\mathsf{\nu }}_t}{\partial {\mathsf{\eta }}^2}\right)\times \\ \times \left(\Omega +2\frac{d\mathsf{\eta }}{dY}\frac{\partial U}{\partial \mathsf{\eta }}\right)+4\frac{d\xi }{dX}{\left(\frac{d\mathsf{\eta }}{dY}\right)}^2\frac{{\partial }^2{\mathsf{\nu }}_t}{\partial \xi \partial \mathsf{\eta }}\frac{\partial V}{\partial \mathsf{\eta }}+\frac{d\xi }{dX}\frac{\partial \Theta }{\partial \xi },\end{array}$$

(11)

$$\frac{\partial \Theta }{\partial \tau }+U\frac{d\xi }{dX}\frac{\partial \Theta }{\partial \xi }+V\frac{d\mathsf{\eta }}{dY}\frac{\partial \Theta }{\partial \mathsf{\eta }}=\frac{d\xi }{dX}\frac{\partial }{\partial \xi }\left[\left(\frac{1}{\sqrt{\mathrm{Ra}\cdot \mathrm{Pr}}}+\frac{{\mathsf{\nu }}_t}{{\mathrm{Pr}}_t}\right)\frac{d\xi }{dX}\frac{\partial \Theta }{\partial \xi }\right]+\frac{d\mathsf{\eta }}{dY}\frac{\partial }{\partial \mathsf{\eta }}\left[\left(\frac{1}{\sqrt{\mathrm{Ra}\cdot \mathrm{Pr}}}+\frac{{\mathsf{\nu }}_t}{{\mathrm{Pr}}_t}\right)\frac{d\mathsf{\eta }}{dY}\frac{\partial \Theta }{\partial \mathsf{\eta }}\right],$$

(12)

$$\begin{array}{l}\frac{\partial K}{\partial \tau }+U\frac{d\xi }{dX}\frac{\partial K}{\partial \xi }+V\frac{d\mathsf{\eta }}{dY}\frac{\partial K}{\partial \mathsf{\eta }}=\frac{d\xi }{dX}\frac{\partial }{\partial \xi }\left[\left(\sqrt{\frac{\mathrm{Pr}}{\mathrm{Ra}}}+\frac{{\mathsf{\nu }}_t}{{\sigma }_k}\right)\frac{d\xi }{dX}\frac{\partial K}{\partial \xi }\right]+\frac{d\mathsf{\eta }}{dY}\frac{\partial }{\partial \mathsf{\eta }}\left[\left(\sqrt{\frac{\mathrm{Pr}}{\mathrm{Ra}}}+\frac{{\mathsf{\nu }}_t}{{\sigma }_k}\right)\frac{d\mathsf{\eta }}{dY}\frac{\partial K}{\partial \mathsf{\eta }}\right]+\\ +P_k+G_k-E,\end{array}$$

(13)

$$\begin{array}{l}\frac{\partial E}{\partial \tau }+U\frac{d\xi }{dX}\frac{\partial E}{\partial \xi }+V\frac{d\mathsf{\eta }}{dY}\frac{\partial E}{\partial \mathsf{\eta }}=\frac{d\xi }{dX}\frac{\partial }{\partial \xi }\left[\left(\sqrt{\frac{\mathrm{Pr}}{\mathrm{Ra}}}+\frac{{\mathsf{\nu }}_t}{{\sigma }_{\mathsf{\epsilon }}}\right)\frac{d\xi }{dX}\frac{\partial E}{\partial \xi }\right]+\frac{d\mathsf{\eta }}{dY}\frac{\partial }{\partial \mathsf{\eta }}\left[\left(\sqrt{\frac{\mathrm{Pr}}{\mathrm{Ra}}}+\frac{{\mathsf{\nu }}_t}{{\sigma }_{\mathsf{\epsilon }}}\right)\frac{d\mathsf{\eta }}{dY}\frac{\partial E}{\partial \mathsf{\eta }}\right]+\\ +c_{1\mathsf{\epsilon }}\left(P_k+c_{3\mathsf{\epsilon }}{G}_k\right)\frac{E}{K}-c_{2\mathsf{\epsilon }}\frac{E^2}{K},\end{array}$$

(14)

Inside the solid walls:

$$\frac{\partial \Theta }{\partial \tau }=\frac{{\mathsf{\alpha }}_{w,air}}{\sqrt{Ra\cdot \mathrm{Pr}}}\left(\frac{d^2\xi }{d{X}^2}\frac{\partial \Theta }{\partial \xi }+{\left(\frac{d\xi }{dX}\right)}^2\frac{{\partial }^2\Theta }{\partial {\xi }^2}+\frac{d^2\mathsf{\eta }}{d{Y}^2}\frac{\partial \Theta }{\partial \mathsf{\eta }}+{\left(\frac{d\mathsf{\eta }}{dY}\right)}^2\frac{{\partial }^2\Theta }{\partial {\mathsf{\eta }}^2}\right),$$

(15)

Inside the glass wall:

$$\frac{\partial \Theta }{\partial \tau }=\frac{{\mathsf{\alpha }}_{g,air}}{\sqrt{Ra\cdot \mathrm{Pr}}}\left(\frac{d^2\xi }{d{X}^2}\frac{\partial \Theta }{\partial \xi }+{\left(\frac{d\xi }{dX}\right)}^2\frac{{\partial }^2\Theta }{\partial {\xi }^2}+\frac{d^2\mathsf{\eta }}{d{Y}^2}\frac{\partial \Theta }{\partial \mathsf{\eta }}+{\left(\frac{d\mathsf{\eta }}{dY}\right)}^2\frac{{\partial }^2\Theta }{\partial {\mathsf{\eta }}^2}\right),$$

(16)

Inside the heat-generating source:

$$\frac{\partial \Theta }{\partial \tau }=\frac{{\mathsf{\alpha }}_{hs,air}}{\sqrt{Ra\cdot \mathrm{Pr}}}\left(\frac{d^2\xi }{d{X}^2}\frac{\partial \Theta }{\partial \xi }+{\left(\frac{d\xi }{dX}\right)}^2\frac{{\partial }^2\Theta }{\partial {\xi }^2}+\frac{d^2\mathsf{\eta }}{d{Y}^2}\frac{\partial \Theta }{\partial \mathsf{\eta }}+{\left(\frac{d\mathsf{\eta }}{dY}\right)}^2\frac{{\partial }^2\Theta }{\partial {\mathsf{\eta }}^2}+P{o}_v\right),$$

(17)

where $G_k=-\frac{{\mathsf{\nu }}_t}{{\mathrm{Pr}}_t}\frac{d\mathsf{\eta }}{dY}\frac{\partial \Theta }{\partial \mathsf{\eta }}$, $P_k={\mathsf{\nu }}_t\left[2{\left(\frac{d\xi }{dX}\frac{\partial U}{\partial \xi }\right)}^2+2{\left(\frac{d\mathsf{\eta }}{dY}\frac{\partial V}{\partial \mathsf{\eta }}\right)}^2+{\left(\frac{d\mathsf{\eta }}{dY}\frac{\partial U}{\partial \mathsf{\eta }}+\frac{d\xi }{dX}\frac{\partial V}{\partial \xi }\right)}^2\right]$, and ${\mathsf{\nu }}_t=c_{\mu }\frac{K^2}{E}$.

To determine the turbulent viscosity ${\mathsf{\nu }}_t$, the Kolmogorov–Prandtl formula ${\mathsf{\nu }}_t=c_{\mu }{k}^2/\mathsf{\epsilon }$ is utilized.

The values of the constants for the k-ε turbulence model are presented in

Table 1. Constants of the turbulence model.

Parameters	cμ	c1ε	c2ε	c3ε	σk	σε	Prt
Values	0.09	1.44	1.92	0.8	1.0	1.3	1.0

. The hydrodynamic and heat transfer are defined by the following dimensionless numbers: the Rayleigh number (Ra) and the Pomerantsev heat number (Po_v). These similarity numbers are defined as follows:

$$P{o}_v=\frac{q_v{L}^2}{{\lambda }_{hs}\Delta T}, Ra=\frac{g\mathsf{\beta }\Delta T{L}^3}{\mathsf{\nu }{\mathsf{\alpha }}_{air}}.$$

When calculating the heat transfer by radiation between the surfaces for each of the surfaces, it is necessary to fully take into account all radiation incidents on it, coming from all possible directions in space. All surfaces are assumed to be diffuse emitters. This means that the intensity of radiation emitted by a given isothermal surface does not depend on the direction. The geometry of the system is taken into account by introducing the view factors $F_{k-i}$. Dimensionless radiation flux $Q_{rad,k}$ is determined by solving the following system:

$$\begin{gathered} Q_{rad,k}=R_k-{\displaystyle \sum _{i=1}^N{F}_{k-i}{R}_i}, \\ {R}_k=\left(1-{\tilde{\mathsf{\epsilon }}}_k\right){\displaystyle \sum _{i=1}^N{F}_{k-i}{R}_i}+{\tilde{\mathsf{\epsilon }}}_k{\left(1-\zeta \right)}^4{\left({\Theta }_k+0.5\frac{1+\zeta }{1-\zeta }\right)}^4. \end{gathered}$$

The dimensionless boundary and initial conditions, with a special algebraic transformation of coordinates being taken into account, are determined in the following form:

at τ = 0

$$\Psi \left(\xi ,\mathsf{\eta },0\right)=\Omega \left(\xi ,\mathsf{\eta },0\right)=K\left(\xi ,\mathsf{\eta },0\right)=E\left(\xi ,\mathsf{\eta },0\right)=0, \Theta \left(\xi ,\mathsf{\eta },0\right)=0.5,{\Theta }_h\left(\xi ,\mathsf{\eta },0\right)=1,$$

at τ > 0

• at the boundary η = 0: $\frac{\partial \Theta }{\partial \mathsf{\eta }}=0$;
• at the boundaries ξ = 0 and $\xi =1+2h/L$: $\frac{\partial \xi }{\partial X}\frac{\partial \Theta }{\partial \xi }=Bi\cdot \Theta$;
• at the boundary $\mathsf{\eta }=1+2h/L$: $\frac{\partial \mathsf{\eta }}{\partial Y}\frac{\partial \Theta }{\partial \mathsf{\eta }}=Bi\cdot \Theta$;
• at the heater surface: $\frac{{\lambda }_{hs}}{{\lambda }_{air}}\frac{\partial {\Theta }_{hs}}{\partial \overline{n}}=\frac{\partial {\Theta }_{air}}{\partial \overline{n}}-N_{rad}{Q}_{rad}$;
• at the inner surfaces of the solid material and air, parallel to the axis Oξ:

  $$\Psi =0, \frac{\partial \Psi }{\partial \mathsf{\eta }}=0, {\Theta }_w={\Theta }_{air}, {\lambda }_{w,air}\frac{\partial \mathsf{\eta }}{\partial Y}\frac{\partial {\Theta }_w}{\partial \mathsf{\eta }}=\frac{\partial \mathsf{\eta }}{\partial Y}\frac{\partial {\Theta }_{air}}{\partial \mathsf{\eta }}-N_{rad}{Q}_{rad};$$
• at the inner surfaces of the solid material and air, parallel to the axis Oη:

  $$\Psi =0,\frac{\partial \Psi }{\partial \xi }=0, {\Theta }_{w(g)}={\Theta }_{air}, {\lambda }_{w(g),air}\frac{\partial \xi }{\partial X}\frac{\partial {\Theta }_{w(g)}}{\partial \xi }=\frac{\partial \xi }{\partial X}\frac{\partial {\Theta }_{air}}{\partial \xi }-N_{rad}{Q}_{rad}.$$

It should be noted that the finite-difference technique is utilized to work out the system of control equations in Equations (10)–(17). A second-order-accurate central-differencing discretization scheme and second-order-accurate upwind difference scheme are employed to discretize various terms. The method of successive over-relaxation is applied to solve Equation (10). A Samarskii locally one-dimensional algorithm is employed to work out Equations (11)–(17), and then a Thomas algorithm is used.

To validate the original computational code (C++ programming language), benchmarking with earlier numerical works [24,25,26] of natural convection in a differentially heated cavity is conducted. The average Nusselt numbers at the hot wall for a broad range of Ra (10⁷–10⁹) are shown in

Table 2. Comparison of average Nusselt numbers (Pr = 0.7).

Ra	[24]	[25]	[26]	Present Data
107	16.79	16.523	–	17.13
108	30.506	30.225	28.78	33.06
109	57.35	–	62.0	60.54

.

The turbulent natural convection within a chamber with the right vertical massive wall is also considered. The results of the work done by Yedder and Bilgen [27] are compared (Figure 3) with the obtained results. It should be noted that the agreement is quite good. The computational domain of the size L × L is discretized $N_{\xi }\times N_{\eta }=120\times 120$ based on a sensitive grid resolution analysis.

3. Results

An analysis of the hydrodynamic and heat transfer for thermogravitational convection and surface radiation in a large-scale chamber was carried out. Calculations were conducted for the following key parameters values: ζ = 0.87, Ra = 10⁹, Pr = 0.71, h/L = 0.1, N_rad = 420.71, τ = 10,000, Po_v = 0.01. The impact of surface emissivity $\left({\tilde{\mathsf{\epsilon }}}_h=0.1-0.9\right)$ and height $\left(l/L=0.1-0.4\right)$ of the heater on the distribution of local (isotherms and streamlines) and integral (average convective and radiative Nusselt numbers) parameters was studied.

Figure 4 illustrates the results for isotherms and streamlines at different values of height of the heat source. One large convective cell was formed in the central part of the enclosure and occupied more than 80% of the airspace. It is interesting to note that a rise of l/L results in a development of minor convective cells, which were placed in the top and bottom right corners of the enclosure. At l/L = 0.4, two low-intensity vortices merged into one convective cell, which determined a clockwise circulation. These flow patterns were observed due to both the cooling of the large-scale enclosure from the outside and the geometric features of the analyzed area. The fluid motion intensifies with an increment of the heater height l/L. This is confirmed by an increase in the maximum value of Ψ inside the core of the main vortex ${|\Psi |}_{\max }^{l/L=0.1}=0.0165<{|\Psi |}_{\max }^{l/L=0.2}=0.0175<{|\Psi |}_{\max }^{l/L=0.4}=0.0201$. It can be seen that along the heater, the thermal boundary layer was formed. A thermal plume was generated at the top of the heat source, causing a stronger temperature gradient at the top of the enclosure and a much lower temperature gradient at the bottom of the enclosure. The plume displacement to the glass wall is observed. The growth of the height of the heater strongly affects the temperature distribution inside the large-scale enclosure. Moreover, an increase in l/L allows for enhancing the average temperature within the solution area.

The constant heat generation inside the power source is described by the Pomerantsev heat number (Po_v). An increment of Po_v leads to an increase in the heat-flow density. However, in the present work, the Pomerantsev heat number has a constant value. The distribution of Nu_total depending on the height of the heater is plotted in Figure 5. The heat exchange between the heater and the air rises for the higher values of l/L as the average total Nu increases with the heater height. For example, at τ = 10,000, the average total Nu is enhanced up to 2.03 times when the value of heater height changes from 0.1 to 0.4.

Table 3. Changes of integral characteristics with l/L at τ = 10,000.

Heater Height l/L	${\overline{\Theta }}_{\mathit{e}\mathit{n}\mathit{c}}$	${|\Psi |}_{\mathbf{max}}$	${\overline{\Theta }}_{\mathit{h}\mathit{s}}$
0.1	0.452	0.0165	0.519
0.2	0.464	0.0175	0.553
0.3	0.477	0.0186	0.580
0.4	0.492	0.0201	0.602

presents the variations of mean temperature inside enclosure ${\overline{\Theta }}_{enc}$, the maximum absolute magnitude of the stream function ${|\Psi |}_{\max }$, and the mean temperature inside the heater ${\overline{\Theta }}_{hs}$. Decreasing the heater height results in a reduction of both the average temperature inside the chamber and the mean temperature inside the heat-generating source.

To clarify the mechanism that governs the radiative heat transfer and hydrodynamics, the temperature field and flow patterns are shown in Figure 6. In this work, when performing the calculations of convective–radiative heat transfer, surfaces are supposed to be diffuse and gray. The phenomenon of thermal transmission between these walls (due to multiple reflection and absorption) is complicated in comparison with a similar process for absolutely black bodies. If, in the general case, the emissivity of the body depends on the wavelength, angle, and temperature of the body, then for the case of diffuse and gray surfaces, the spectral emissivity and absorbance depend only on temperature. A growth of the surface emissivity of the heater has no significant effect on the flow structure, except for a weak decrease in the power of convective circulation ${|\Psi |}_{\max }^{{\tilde{\mathsf{\epsilon }}}_h=0.9}=0.0172< {|\Psi |}_{\max }^{{\tilde{\mathsf{\epsilon }}}_h=0.1}=0.0178$. The air above the heater rises towards the top and glass walls. The maximum values of air velocity were observed directly above the heat-generating source.

The impact of surface emissivity of the heater on both the mean temperature inside the enclosure and the mean temperature inside the heater is also presented in

Table 4. Variations of various parameters for l/L = 0.2, τ = 10,000.

Heater Surface Emissivity	${\overline{\Theta }}_{\mathit{e}\mathit{n}\mathit{c}}$	${\overline{\Theta }}_{\mathit{h}\mathit{s}}$
0.1	0.463	0.566
0.5	0.464	0.553
0.9	0.465	0.543

.

In analyzing the data in the Table 4, we can conclude that under the conditions of convective–radiative energy exchange in an enclosure with a local energy unit under convective energy exchange with the ambient, it is possible to reduce the average temperature inside the energy source by increasing the values of heater surface emissivity ${\tilde{\mathsf{\epsilon }}}_h$.

In Figure 7, variations of Nu at the thermally-generating element surface are shown. This figure illustrates that enhancing the surface emissivity of the heater leads to different effects for the average radiative and convective Nu. It is clearly shown that the energy transfer due to radiation is significantly reduced with a decrease in ${\tilde{\mathsf{\epsilon }}}_h$. At τ = 10,000, Nu_conv decreases up to 43% at the changing of the heater surface emissivity from 0.1 to 0.9. It should be noted that the decrease in the values of the dimensionless integral heat transfer coefficient with time is associated with the formulation of initial and boundary conditions.

4. Conclusions

The combined influence of thermogravitational convection and surface thermal radiation on hydrodynamic and energy transmission within a rectangular large-scale chamber containing a local heat-generating element was studied numerically in this work. The obtained outcomes were found to be in excellent agreement with other available experimental and numerical data for turbulent free convection in enclosures. The system of governing dimensionless equations was solved in stream function–vorticity form by the finite-difference method. Based on the enclosure configuration, the effect of wall emissivity and the height of the thermally-producing unit on the distribution of integral (average radiative and convective Nusselt numbers) and local (isotherms and streamlines) parameters was investigated. A large convective cell was formed in the central part of the enclosure, occupying more than 80% of the airspace. At l/L = 0.4, a minor vortex (clockwise circulation) was generated near the glass wall of the enclosure. It was ascertained that raising the wall emissivity of the heater resulted in a slight reduction of the convective energy transference rate within an enclosure. The influence of the heater surface emissivity on both the mean temperature inside the enclosure and the mean temperature inside the heater was also investigated.