Numerical Simulation of Air Flow and Temperature Distribution in Volumetric Solar Receiver Consisting of a Porous Medium

Ceramics are materials with very favorable thermophysical properties, they are mainly used for the design of volumetric receivers for their high resistance to strong temperature gradients. Their plastic deformation is low and their melting temperature is very high (up to 2000℃), their density is low with a strong endurance. For that, we are principally interested to study air flow and convective heat transfer in a solar receiver constituted of a ceramic porous medium. This medium is presented by a series of periodic (3D) structures, in which, one element is constructed according to the Kelvin model.

Certainly, the radiations generated from heliostats field are focalized to the solar receiver, in which the heat is transferred from the solid phase to the air as it passes through the porous medium that absorb solar radiation in depth. The domain of receivers presents a great interest for a green electricity production, in this way, many researchers invested their efforts to increase the efficiency of porous media in volumetric receivers. Most of the works are concentrated in the study of heat transfer phenomena and pressure drop parameters in open cell foams numerically and analytically for different open cell foams with randomized forms. The medium has tortuosity’s which allow the increase of the heat exchange, and the reduction of the pressure losses compared to packed bed of spheres Kumar et al. [1].

Different forms of ceramic foams geometries have been a subject of thermic and dynamic studies; we cite for example, the cubic model of Dul’nev Lacroix [2] and Krishnan [3], Weaire–Phelan’s unit cell [4], Kelvin’s tetrakaidecahedral model Boomsma [5], Richardson [6], and 3D Laguerre-Voronoi model studied by Nie et al. [7]. The Kelvin’s model is widely used because of its useful construction. 

Among the works done on the kelvin structure (Tetrakaidecahedra), Lucci et al. [8] presented a parametric CFD analysis of Kelvin cell structures to study foams in a controlled environment. In particular their results are about 20% higher than that predicted by the correlation of Groppi et al. [9]. Zhang [10] proposed a modified Kelvin model for high mechanical property open-cell metal foams and investigated its application in thermal simulations. Ferrari et al. [11] analyzed three-dimensional fluid dynamics in order to evaluate the convective heat transfer coefficient in regular and non-regular structures. They also evaluated the pressure drop to extract the best performance given by the different structures.

Other research was focused on numerical simulations based on detailed foam structure models to predict effective properties of porous foams associated with fluid flow and/or heat transfer. Kuwahara et al. [12] studied numerically the interfacial convective heat transfer coefficient using a single structural unit. Petrasch et al. [13] determined the permeability, Dupuit-Forchheimer coefficient, and interfacial heat transfer coefficient in 10 ppi porous ceramic through a 3D digital sample representation generated from X-ray tomographic scans. They solved the governing equations of fluid flow via the finite volume direct pore-level numerical simulation. Mendes et al. [14] provided a simplified modeling approach for estimating the effective thermal conductivity of open cell foam using a finite volume method. Wu et al. [15] simulated numerically the convective heat transfer in ceramic foams represented by an idealized packed tetrakaidecahedra structure. Xu et al. [16] used the numerical method to analyze the typical influences of the porosity, average particle diameter and inlet velocity on the temperature distributions.

This work has as a principal objective, to participate in the performance of the kelvin structure in order to use it in the volumetric solar receiver. The main equations (continuity, momentum and energy) are solved simultaneously to determine dynamic and heat transfer characteristics between the flowing fluid and ceramic foam which is a periodic packed tetrakaidecahedron regular structure. The SSTkw turbulent model is used to appropriate the predict simulation. According to Menter et al. [17, 18], Vijay Garg and Ameri [19] this model gives results with high quality for a variety of heat transfer cases. Instationary and three-dimensional simulations are realised with the FLUENT software package. Calculations are made for one and two rows of structures. Our validation is based on the comparison with data and results provided in Wu et al. [15].

In porous media, in order to use the Navier Stokes equations, it is necessary to have a saturated medium, the principal formulas used to resolve the problem are derived from the second law of Newton and the first principle of thermodynamics.

The conservation of mass equation for a control volume in differential form is defined as:

$\frac{\partial \rho}{\partial \mathrm{t}}+\frac{\partial\left(\rho \mathrm{U}_{\mathrm{j}}\right)}{\partial \mathrm{x}_{\mathrm{j}}}=0$     (1)

The conservation of momentum in a static reference frame is given by:

$\frac{\partial\left(\rho \mathrm{U}_{\mathrm{i}}\right)}{\partial \mathrm{t}}+\frac{\partial\left(\rho \mathrm{U}_{\mathrm{i}} \mathrm{U}_{\mathrm{j}}\right)}{\partial \mathrm{x}_{\mathrm{j}}}=-\frac{\partial \mathrm{P}}{\partial \mathrm{x}_{\mathrm{i}}}$$+\frac{\partial}{\partial x_{j}}\left(\mu\left(\frac{\partial U_{i}}{\partial x_{j}}+\frac{\partial U_{j}}{\partial x_{i}}\right)-\left(\frac{2}{3} \delta_{i j} \mu \frac{\partial U_{k}}{\partial x_{k}}\right)-\rho \overline{u_{i} u_{j}}\right)$   (2)

The energy equation is presented by:

$\frac{\partial(\rho \mathrm{E})}{\partial \mathrm{t}}+\frac{\partial\left(\rho \mathrm{U}_{\mathrm{j}} \mathrm{E}\right)}{\partial \mathrm{x}_{\mathrm{j}}}=\frac{\partial}{\partial \mathrm{x}_{\mathrm{j}}}\left(\lambda_{\mathrm{c}}\left(\frac{\partial \mathrm{T}}{\partial \mathrm{x}_{\mathrm{j}}}\right)+\mathrm{U}_{\mathrm{i}} \tau_{\mathrm{ij}}-\rho \mathrm{C}_{\mathrm{v}} \overline{\mathrm{u}_{i} \mathrm{T}^{\prime}}\right)$     (3)

$E=h-\frac{p}{\rho}+\frac{U_{i}^{2}}{2}, \frac{p}{\rho}=R T$     (4)

To carry out this study, we took into account the variations of the following parameters: Reynolds number effect at the entrance and the porosity of a part of the computational domain that is presented by two rows of the matrix structure (Figure 4). The fluid enters the conduit, with an ambient temperature T_f equal to 300 K, by interacting with the structures of the source of the domain which are at the temperature T_S equal to 330 K. During this interaction the fluid extracts a quantity of heat until reaching the local equilibrium temperature (the temperature of the fluid is equal to that of the solid T_s = T_f).

5.1 Grid dependence

Because of the foam shape and the complex geometry, the flow is very tortuous; which generates strong temperature and speed gradients, so here we have to present a high-quality mesh. A study of the mesh sensitivity was conducted by checking the temperature which is a very sensitive parameter.

5.png

Figure 5. Grid dependence: Dimensionless temperature profiles

Figure 5 shows the dimensionless profiles of the temperature (T_s-T_f) /(T_s-T_fin) [15] for two grids; (13520) and (346485). The superposition of these profiles shows a very small difference (3.8%) between the values. Therefore, we opted for the first mesh to reduce the computation time.

5.2 Temperature behavior

The first results present the temperature distribution along the direction of flow. The air begins to extract the thermal energy from the structures until it reaches the local equilibrium temperature, which has been represented by the middle plane of the thermal field (Figure 6) and (Figure 7).

6.png

Figure 6. The temperature field along the flow direction

7.png

Figure 7. Variation of static temperature with horizontal position Z, Re =129

8.png

Figure 8. Temperature profiles with different velocities

9.png

Figure 9. Temperature profiles with one and two rows

Figure 8 thus presents the behavior in term of the temperature of the air in the matrix with varying inlet velocities. The elevated flow rates allow better temperature diffusion, with a distinct strong energy extraction. When the contact between the air and the structures lasts longer, the energy transfer is better and the fluid reaches the equilibrium temperature at a distance very far from the outlet.

We note that low velocities allow the extraction of energy more than the high ones. The velocity inlet (Reynolds number) affects directly the temperature behavior in the matrix medium [15].

In order to optimize the solar receiver geometry as mentioned before, we expose the results of the second geometry with two rows and compare them with the one’s one Figure (1-a, 1-b). We note that the fluid reaches the temperature of the solid in a very advanced position for the case of two rows, unlike for one row (The half of the first distance); where the propagation of heat generates a thermal equilibrium in position far downstream (Figure 9).

5.3 Pressure and velocity behavior

The shape of the structure of the Kelvin matrix causes irregularities which appear in the latter. The particles of fluid meet regions of separation and re-attachment which are at the origin of the dissipation of the kinetic energy which will be amplified by strong gradients of speed and temperature. (Figure 10) shows streamlines colored by velocity in the symmetric plane along the z axis.

10.png

Figure 10. Streamlines colored by velocity

This shows that there is a separation in streamlines in the foam, which is observed in the work of Nie et al. [7]. The majority of current lines encounter an impact zone in the foam, and depending on the speed - temperature gradients producing, the phenomena generated, favor the transfer by convection and rapid diffusion of heat into the fluid medium.

The Figure (11.a) shows the distribution of the mean velocity along the flow direction. Although the Reynolds number is small, the flow always remains turbulent because of the continuous change of the directions of the particles and the penetration in the pores by causing accelerations according to the conservation of the mass, this phenomenon is observed by  Nie et al. [7]. The same remarks are also observed in other studies, for example [20, 21]. A larger cross section of the spacer can cause a greater separation effect, and therefore a greater reduction in heat transfer.

The effects of friction indicate an internal heat generation followed by a reduction of the density of the fluid, to make pass the same mass flow, the continuity forces the velocity to increase and automatically the static pressure to decrease (Figure 11.b).

11.png

Figure 11. Velocity (a) and Pressure field (b)

According to Nie et al. [7], the complex structure of pores causes a complex pressure distribution inside the foam and stagnation areas are observed in front of the struts.

12.png

Figure 12. Static pressure profiles: one and two rows

Figure 12 presents the static pressure profile with one and two rows, which decreases linearly. The presence of more structures causes a high pressure drop. However, the decreasing porosity generates an increasing pressure drop.

The Figure 13 presents a comparison between Ergun correlation calculating the pressure drop and our simulation results concerning the same parameter.

The pressure drop increases with increasing velocity.  Compared with Ergun’s model, the predicted value agrees well at a velocity of 3.4 m/s. The curve follows the same draft as the proposed model, larger in the range of 1<V< 3.4 and smaller for V > 3.4 It can be noted that the simulated data are in good agreement With an Ergun's [22] model because it is the most successful model for predicting the pressure drop of porous medium [23] especially for high porosities: The Kelvin structure porosity is equal to 0.92.

• model 1 is a correlation of Ergun. [22]:           

$\frac{\mathrm{dp}}{\mathrm{dz}}=\frac{(1-\varepsilon)^{2} \mu \mathrm{V}}{150 \mathrm{d}^{2} \varepsilon^{3}}+\frac{(1-\varepsilon) \rho \mathrm{V}^{2}}{1.75 \mathrm{d} \varepsilon^{3}}$      (5)

13.png

Figure 13. Pressure drop comparison for one row and velocity varieties, our simulations and theoretical prediction Ergun [22]