**The Optimization of the Thermal Performances of an Earth to Air Heat Exchanger for an Air Conditioning System: A Numerical Study**

#### **Adriana Greco <sup>1</sup> and Claudia Masselli 2,\***


Received: 2 November 2020; Accepted: 2 December 2020; Published: 4 December 2020

**Abstract:** The aim of this paper is to research the parameters that optimize the thermal performances of a horizontal single-duct Earth to Air Heat eXchanger (EAHX). In this analysis, the EAHX is intended to be installed in the city of Naples (Italy). The study is conducted by varying the most crucial parameters influencing the heat exchange between the air flowing in the duct and the ground. The effect of the geometrical characteristics of the duct (pipe length, diameter, burial depth), and the thermal and flow parameter of humid air (inlet temperature and velocity) has been studied in order to optimize the operation of this geothermal system. The results reveal that the thermal performance increases with length until the saturation distance is reached. Moreover, if the pipe is designed with smaller diameters and slower air flows, if other conditions remain equal, the outlet temperatures come closer to the ground temperature. The combination that optimizes the performance of the system, carried out by forcing the EAHX with the design conditions for cooling and heating, is: D = 0.1 m s<sup>−</sup>1; v = 1.5 m s<sup>−</sup>1; L = 50 m. This solution could also be extended to horizontal multi-tube EAHX systems.

**Keywords:** geothermal energy; earth-to-air; horizontal pipe; 2D model; air conditioning; renewable energy sources; thermal performances; energy efficiency; parametric study

#### **1. Introduction**

#### *1.1. General Concepts and Context*

Heating Ventilation and Air Conditioning (HVAC) systems contribute 10–20% to worldwide energy consumption [1]. A relevant amount (20–40%) of the energy consumption attributed to HVAC is required for building air conditioning. Thus, in addition to the increase in the use of renewable energy sources, energy saving solutions are also being adopted [2–6]. The basic imperative recommended by energy policies is to consider energetically improved solutions that often could not be satisfied only through the vapor compression technology, but which have limits linked to the use of refrigerants with high Global Warming Potential [7–10]. A viable path is a Not-In-Kind cooling technology [11–13] where the refrigerants are solid-state materials showing caloric effects, i.e., their Global Warming Potential is zero [14]. Although they are a promising solution, the bottleneck is the limited range of use of these technologies due to the limited caloric effects of the materials. This can constitute some difficulties in the development of environmental conditioning systems that operate on a wide range. Another valid solution is the utilization of renewable energy sources that must increasingly become a shared responsibility. Currently, only 14% of the global energy demand is satisfied by means of renewable energy [15], and many renewable energy sources are employable for these purposes. Among them, geothermal energy has been revealed to be very promising for employment as an answer to the energy demand attributable to HVAC systems. The recently emerging geothermal hybrid solutions with solar collectors of PV panels are also very promising [16–18] for building conditioning. Geothermal energy is very promising because of the properties of the soil to be, from a certain depth onward, be a constant temperature throughout the year. Specifically, during winter the temperature is higher than the temperature of the external air, and lower during summer [19]. For this reason, the soil assumes the double role of heat sink/source despite it functioning in cooling/heating operation mode.

There are three geothermal systems that exploit this property and could be used for air conditioning [20]: earth homes (the idea of using a partially buried building to minimize the heating and cooling loads), ground-source heat pumps (heat pumps where the water is the secondary fluid that flows underground to heat transfer with the soil), and Earth to Air Heat eXchangers (EAHX). The latter are geothermal systems formed by a certain number of pipes buried to a depth that allows the exploitation of the undisturbed temperature properties of the ground. The external air flows within the pipes and transfers heat from/to the soil.

An EAHX could be arranged through horizontal or vertical banks of pipe. These could be placed in parallel with each other or following other configurations. The specificity of the design depends on the project and the space available for the installation. However, they must be designed to optimize the performance. Although valid research exists [21] on numerical optimization of EAHX with vertical ducts, in most cases, EAHX is formed by horizontal ducts [22].

A well-designed EAHX should satisfy the requirement to bring the temperature of the external air as close as possible to the undisturbed temperature of the ground. Furthermore, efficiency and economic aspects should also be optimized [23].

#### *1.2. State of the Art*

There are several parameters that influence the performances of the EAHX: they are both geometric factors of the system and the thermo-physical parameters of the soil and air.

The geometric configuration of the tube of a heat exchanger is a key factor both because its performance depends on it, and because an adequate sizing optimizes both the excavation costs and the area occupied by the plant. Selamat et al. [24] researched, through numerical modelling, the best configuration for an EAHX during winter in the city of Dublin (Ireland). A sensitivity analysis was carried out by varying the key parameters influencing the thermal performance of the EAHX. The analysis clearly shows that an increase in the length of the tube (from 30 to 70 m), a decrease in the diameter of the tube (from 150 to 100 mm) and a decrease in the air velocity (from 15 to 5 m s−1) provide an improvement of the system heating capacity. The study was also focused on the thermal performance of 4 parallel tubes: 1.5 m distance, 30 m length and 125 mm diameter, with an air velocity of 8 m s−1. The simulations clearly show that with the multiple-duct configuration, due to heat loss from adjacent pipes, there was a 0.6 ◦C as decrement of the air temperature.

Mathur et al. [25] introduced an experimental study where an earth to air spiral heat exchanger was designed with a 60-m pipe length, a diameter of 0.1 m, and a 1-m distance among the spirals, and it was buried at a depth of 3 m. The results obtained were compared to straight tube EAHX-system results: they noted that the coefficient of performance of the latter configuration is 5.94 in the cooling phase and 1.92 in the heating phase, whereas for a spiral configuration these values are 6.24 and 2.11, respectively. As Benrachi et al. reported [26], the use of the spiral configuration can be an interesting alternative to straight tube exchangers if the space available for excavation is limited. This advantage in terms of occupied area tends to disappear when dealing with the installation of several pipes in parallel.

Many studies focus on parametric research into the influence of pipe diameter on system performance. Sodha et al. [27] made a comparison between single-tube and multiple-tube exchangers. It has been observed that, with the same air flow, the heating and cooling potential increases as the number of tubes of smaller diameter increase, because of the consequent augmentation of the heat

exchange area. With the increase in pipe diameter, the absolute value of temperature variation reduces, reducing the heating/cooling capacity of the system.

Wu et al. [28] estimated the daily cooling capacity of EAHX with a 0.4 m and 0.6 m pipe diameter. The results obtained (respectively 43.2 kWh and 74.6 kWh) underline the advantage that can be gained from using a smaller diameter. Mihalakakou et al. [29,30] evidenced that an increase in the diameter of the pipe corresponds to a decrease in the convective heat transfer coefficients which leads to a higher temperature at the outlet of the EAHX in summer and lower temperature in winter. Typical diameter values range from 0.1 m to 0.3 m, up to 1 m for commercial applications. Goswami and Biseli [31] suggested that the ideal performance for single-tube EAHX is achieved with a diameter of 0.3 m; for systems with parallel pipes, the ideal diameter varies between 0.2 m and 0.25 m. Sehli et al. [32], considering an EAHX placed in the Algerian arid climate and operating in cooling mode, estimated the variation of air temperature at the pipe outlet as a function of the tube form factor (the ratio between the length and the diameter of the pipe). They present the effect of both the form factor and the Reynolds number on the outlet air temperature, for an ambient temperature of 45 ◦C and a tube depth set at 4 m. As the form factor increases, the outlet air temperature decreases because a longer tube provides a longer path to the air, favoring heat exchange. As the Reynolds number increases, the outlet air temperature also increases, since the air remains in the tube for a shorter time with respect to the cases with smaller Reynolds numbers. The study therefore concludes that the optimal value of the form factor σ was 250.

Hanby et al. [33] reported that, for a fixed length of the duct, at fixed air flow, there is always an optimal diameter of the pipe that minimizes the energy consumption of the system. Niu et al. [34] found that, for an EAHX working in cooling mode, the advantage of a smaller diameter is not only linked to a faster lowering of the temperature in the duct, but also to a lower presence of humidity in air passing through the heat exchanger. With small pipe diameters, the air temperature drops faster than when using larger diameters. For this reason, the air temperature recorded at the outlet of the tube with a pipe diameter of 0.7 m is 17 ◦C, whereas for a diameter of 0.3 m the temperature is 13 ◦C.

The length of the tube is the main geometric factor affecting the performance of the earth to air heat exchanger. It is implicit that it must be chosen appropriately to reduce the initial costs of construction and installation of the system without compromising its performance. As the length of the tube increases, the temperature difference between the inlet and outlet ends of the air increases. Derbel and Kanoun [35] found out that the energy load of the heat exchanger increases with the length of the buried pipe; this is mainly due to the growth in the crossing time of the air flow. Furthermore, Benhammou and Draoui [36] asserted that the efficiency of heat transmission does not increase from a certain tube length onward: this is called the saturation length, and increases with increasing air flow rate. Lee et al. [37] noted that, with respect to their experimental test study focused on four locations of USA with different climates, there were no further advantages in using pipes longer than 70 m. Ahmed et al. [38] established that the pipe length is the greatest influencing factor on thermal performances of an EAHX. They chose four different pipe lengths (7.5 m, 15.0 m, 30.0 m and 60.0 m), achieving better performance with maximum length, while the system was operating in cooling mode.

The study of the temperature profile in the soil and its variation with depth is essential for the optimization of an EAHX system. The temperature distribution in the soil is mainly influenced by its physical properties and its structure. The soil temperature varies with time (daily and seasonal temperature) and with depth. The amplitude of this variation decreases with increasing soil depth until a certain value reaching the so-called Undisturbed Ground Temperature (UGT). This value is typically reached at a depth of 2–4 m, but in some latitudes, it is reached at greater depths. Hanby et al. [33] showed that there is an ideal diameter at each depth that allows the system to save energy. They also show that the energy saving increases by increasing the tube depth (from 2 to 4 m). Popiel et al. [39] measured the soil temperature at different depths and found that short-term variations occurred down to a depth of 1 m. Sanusi et al. [40] investigated on an EAHX buried at different depths: 0.5 m, 1.0 m

and 1.5 m in a humid and warm climate like that of Malaysia. The maximum soil temperatures at 0.5 m, 1.0 m and 1.5 m are respectively: 28.3 ◦C, 28.5 ◦C and 28.6 ◦C in the wet season and 30.3 ◦C, 30.1 ◦C and 30.1 ◦C in the dry season. The optimum depth for burying the EAHX considering both economic and thermal performance terms is 1 m. Mihalakakou et al. [41] analysed the influence of the depth of tube of the EAHX on the cooling performances considering the values of 1.2 m, 2.0 m, and 3.0 m and they observed that the best results are obtained for the deepest buried tube. Ahmed et al. [38] buried the duct at 0.6 m, 2.0 m, 4.0 m, and 8.0 m, finding that the best cooling performances are associated with a 8-m burieb pipe. Wu et al. [28] evaluated the thermal performances of the EAHX at two depths (1.6 m and 3.2 m) detecting that for air cooling in summer, the temperature at the EAHX outlet varies between 27.2 ◦C and 31.7 ◦C at 1.6 m and between 25.7 ◦C and 30.7 ◦C at 3.2 m. Generally, as the depth of the pipe increases, the potential for heating and cooling increases, but beyond a certain depth there is no noticeable increase in performance. Badescu [42] noticed that for depths greater than 4 m, the performance of the system remains unchanged, so a depth value of 2 m is indicated as a good compromise between the excavation costs (which increase with increasing depths) and the annual variation of the temperature (which decreases with increasing depths). Hermes et al. [43] conducted a research on the soil of Rio Grande, Brazil and underlined how problematic the installation of a heat exchanger can be in a coastal city, where the aquifer can be very close to the ground surface. They estimated that below the depth of 2 m the temperature remains constant both in the warm months and in the cold ones.

Lee et al. [37] studied the effect of the air velocity inside the tube on the outlet air temperature in four different locations (Key West, Peoria, Phoenix and Spokane) in both heating and cooling mode. They found that as the air flow velocity increases (from 2 to 14 m s<sup>−</sup>1) the air temperature at the outlet of the pipe in cooling mode also increases, as the air spends less time in the pipe and less time in contact with the ground. Likewise, the range and rate of this increase was a function of the different locations depending on the soil conditions.

Niu et al. [34] analyzed the effect of air velocity on the thermal performances of an EAHX in cooling mode. The study was conducted by varying the air velocity from 0.5 to 2.5 m s<sup>−</sup>1. The smaller the air velocity was, the faster the air temperature decreasing rate was, and the lower the outlet air temperature was. Bansal et al. [44] conducted an experimental study considering air velocity as 2.0, 3.0, 4.0, and 5.0 m s<sup>−</sup>1, in the case of operation of the EAHX in heating mode. They observed that 2 m s−<sup>1</sup> is the velocity at which the most noticeable temperature increase was registered. In the case of cooling [45], the most pronounced decrease in temperature was always recorded for an air speed of 2ms−1. For the cooling operation, Benhammou and Draoui [36] observed the increase in the outlet temperature of the exchanger as the velocity of the inlet air increased. As a consequence of increasing the speed from 1 to 3 m s<sup>−</sup>1, the outlet temperature increases by 5.6 ◦C. This phenomenon is attributed to the reduction of heat transfer from the air to the ground. While the air velocity rose from 1 to 3 m s−1, the average daily efficiency decreased by 31.6%, together with a significant increase in the COP which can be attributed to the increasing of the pressure losses.

#### *1.3. Aim of the Paper*

Beyond the generally accepted results on qualitative trends following the variation of the key parameters for the EAHX, the above introduced state of the art also presents peculiarities that in some cases may seem conflicting. Most simulations focus on only very limited configurations, most of them based on simple one-dimensional models that can hardly predict the real behavior of an EAHX. In addition, most research ignored the latent heat transfer related to the water vapor condensation during humid air cooling in the summer season.

The aim of this paper is to develop an extensive analysis of the parameters that would optimize the thermal performances of a horizontal single-duct system to give clear indications for the design of a real application for an EAHX. The behavior of an EAHX is influenced by the geographical and climatic conditions where the system is installed. There are national or global directives that suggest the climatic design parameters in projecting a ground to air heat exchanger system. The ASHRAE identifies [46] the reference climate design parameters such as intensity of the solar irradiation, temperature and relative humidity of the external air, both in winter and in summer operation modes for many cities. In this paper, the EAHX is intended to be installed in the city of Naples (Italy), which has a typical Mediterranean climate. The study is conducted by means of a numerical two-dimensional model based on finite element method and experimentally validated [47]. The step of the investigation concerning the optimization of geometrical factors of the duct and air velocity is conducted by imposing the design parameters for cooling and heating operation modes prescribed by ASHRAE. Subsequently, once the optimal parameters have been identified, the effect of the external air temperature on the thermal performances of the geothermal system is analyzed.

#### **2. Methods**

#### *2.1. The Design and the Assumptions Made*

The numerical tool employed in the analysis is a two-dimensional numerical model, based on a finite element method, of a horizontal single-tube earth to air heat exchanger surrounded by the ground domain. A schematic of the model design is visible in Figure 1. The horizontal disposition was adopted since a vertical EAHX is generally connected with much more expansive installation and maintenance costs, whereas the reason for the single-tube choice lies in the desire to preserve as much genericity as possible, to allow the optimization. To this purpose, the parameters of the tube such as burial depth (*z*), diameter (*D*), length (*L*), and inlet air velocity (*v*), are held as variables. Once the optimal combination of design parameters has been identified, the number of parallel tubes can be increased to supply larger demands of volumetric flow rate.

**Figure 1.** A schematic of the design of the numerical model of the horizontal Earth to Air Heat eXchanger (EAHX).

The mathematical model is forced on the earth to air heat exchanger considering the following assumptions:


Humid air is the fluid crossing the pipe, whose thermodynamic properties (such dry bulb temperature; relative and specific humidity) are punctually evaluated (timely and spatially) in the numerical simulation of the model. Furthermore, the tool allows us to identify and quantify the amount of water condensation of the humid flow.

#### *2.2. The Thermal Resistance Approach in the EAHX System*

According to the design of Figure 1 and following the approach of thermal resistance [48], the model could be schematized as the series of three resistances:

$$R\_{c\eta} = R\_{\text{ground}} + R\_{\text{pipe}} + R\_{c\prime} \tag{1}$$

Both the conductive thermal resistances of the surrounding ground and the pipe are evaluated as:

$$R\_{graund} = \frac{\ln\left(\frac{r\_{graud}}{r\_{pinv,ext}}\right)}{2\pi k\_{graund}L},\tag{2}$$

$$R\_{pipc} = \frac{\ln\left(\frac{r\_{pipec,int}}{r\_{pipec,int}}\right)}{2\pi k\_{ground}L},\tag{3}$$

The thermal resistance due to the convective heat exchange between the air flowing in the tube and the pipe is:

$$R\_{\rm c} = \frac{1}{2\pi r\_{\rm pipe, int} LLI'} \tag{4}$$

Considering a PVC pipe (k = 0.16 W m−<sup>1</sup> K<sup>−</sup>1; ρ = 1380 kg m<sup>−</sup>3; c = 900 J kg−<sup>1</sup> K<sup>−</sup>1), the orders of magnitude for *Rc* and *Rground* are 10−<sup>2</sup> <sup>÷</sup> 10−<sup>3</sup> K W−<sup>1</sup> whereas 10−<sup>4</sup> K W−<sup>1</sup> is the dimension proper of *Rpipe*.

For this reason, in our model, the thickness of the tube is not considered since the thermal resistance associated with the conduction in the thickness of the tube is negligible with respect to the two others. As a matter of fact, for modeling needs with a software CFD, such a small thickness, compared to the orders of magnitude of the dimensions of the other domains of the model (ground and air-occupied area in the tube), would have required a much smaller thickening of the FEM grid. Thus, the dimensions of the elements of the mesh should have assumed values small enough to make it really difficult to converge the model in finite or acceptable times.

The problem is very common in the literature and this approximation (calculation of the temperature field using CFD basing of FEM, therefore not through a resistive approach) is usually accepted for approaches similar to ours [45,49,50].

#### *2.3. The Mathematical System*

The fluid domain is described by the following differential equations:

• the mass conservation of the humid air:

$$\frac{\partial \rho}{\partial t} + \nabla \left(\rho \overrightarrow{\boldsymbol{\upsilon}}\right) = \dot{S}\_{\text{mass}} \tag{5}$$

where . *Sm* is a negative term that represents the mass of the condensed water;

• the conservation momentum is guaranteed by the Navier-Stokes equations for turbulent air flow:

$$
\rho \frac{\partial \overrightarrow{\boldsymbol{\upsilon}}}{\partial t} + \rho \left(\overrightarrow{\boldsymbol{\upsilon}} \cdot \nabla\right) \overrightarrow{\boldsymbol{\upsilon}} = \nabla \cdot \left[ -p \overrightarrow{\boldsymbol{I}} + (\mu + \mu\_T) \left[ \nabla \overrightarrow{\boldsymbol{\upsilon}} + \left( \nabla \overrightarrow{\boldsymbol{\upsilon}} \right)^T \right] \right] \tag{6}
$$

where μ*<sup>T</sup>* is the turbulent viscosity defined as:

$$
\mu\_T = \rho \mathbb{C}\_{\mu} \frac{\mathbb{K}}{\mathbb{g}'} \tag{7}
$$

with *C*μ, that is one of the constants of the K-εˆ model for turbulent flow [47].

• the energy equation for the air flow:

$$\frac{\partial(\rho E)}{\partial t} + \nabla \cdot \left[ \overrightarrow{\boldsymbol{v}} \left( \rho \boldsymbol{E} + \boldsymbol{p} \right) \right] = \nabla \cdot \left| \boldsymbol{k}\_{eff} \nabla \mathbf{T} - \sum\_{\mathbf{j}} \mathbf{h}\_{\mathbf{j}} \overrightarrow{\mathbf{J}}\_{\mathbf{j}} + \left( \boldsymbol{\pi}\_{\text{eff}} \overrightarrow{\boldsymbol{v}} \right) \right|, \tag{8}$$

where *keff* is the effective conductivity defined as the sum of the conventional thermal conductivity of the fluid (*k <sup>f</sup>*) and the thermal conductivity of the turbulent flow (*kT*) and thus modeled as:

$$k\_{eff} = k\_f + k\_T \tag{9}$$

• Using the K-εˆ model for turbulent flow, the turbulence kinetic energy equation is:

$$\frac{\partial(\rho K)}{\partial t} + \rho \overrightarrow{\boldsymbol{\nu}} \cdot \nabla K = \nabla \cdot \left( \left( \mu + \frac{\mu\_T}{\sigma\_K} \right) \nabla K \right) + \boldsymbol{P}\_K - \rho \boldsymbol{\varepsilon}, \tag{10}$$

where *PK* can be evaluated as:

$$P\_K = \mu\_T \left( \nabla \overrightarrow{\boldsymbol{\upsilon}} : \left( \nabla \overrightarrow{\boldsymbol{\upsilon}} + \left( \nabla \overrightarrow{\boldsymbol{\upsilon}} \right)^T \right) - \frac{2}{3} \left( \nabla \cdot \overrightarrow{\boldsymbol{\upsilon}} \right)^2 \right) - \frac{2}{3} \rho K \nabla \cdot \overrightarrow{\boldsymbol{\upsilon}} \tag{11}$$

The specific dissipation rate equation is:

$$\frac{\partial(\rho\varepsilon)}{\partial t} + \rho\overrightarrow{\boldsymbol{v}} \cdot \nabla\varepsilon = \boldsymbol{\nabla} \cdot \left( \left( \mu + \frac{\mu\_T}{\sigma\_\ell} \right) \nabla\varepsilon \right) + \boldsymbol{\mathcal{C}}\_{\ell 1} \frac{\boldsymbol{\xi}}{K} \boldsymbol{P}\_K - \boldsymbol{\mathcal{C}}\_{\ell 2} \rho \frac{\boldsymbol{\xi}^2}{K} \tag{12}$$

Table 1 reports the experimental constants of the K-εˆ model.


**Table 1.** The K-εˆ model constants.

Where the fluid flows in laminar motion, the turbulent flow parameters *kT* and μ*<sup>T</sup>* become zero and Equations (2) and (4) are incorporated respectively in:

σε 1.3

$$
\rho \frac{\partial \overrightarrow{\boldsymbol{\upsilon}}}{\partial t} + \rho \left(\overrightarrow{\boldsymbol{\upsilon}} \cdot \nabla\right) \overrightarrow{\boldsymbol{\upsilon}} = \nabla \cdot \left[ -p \overrightarrow{\boldsymbol{I}} + \mu \left[ \nabla \overrightarrow{\boldsymbol{\upsilon}} + \left(\nabla \overrightarrow{\boldsymbol{\upsilon}}\right)^{T} \right] \right] \tag{13}
$$

$$\frac{\partial(\rho E)}{\partial t} + \nabla \cdot \left[ \overrightarrow{\boldsymbol{v}} \left( \rho \boldsymbol{E} + \boldsymbol{p} \right) \right] = \nabla \cdot \left| k \boldsymbol{\nabla} \mathbf{T} - \sum\_{\mathbf{j}} \mathbf{h}\_{\mathbf{j}} \overrightarrow{\mathbf{J}}\_{\mathbf{j}} + \left( \boldsymbol{\tau} \cdot \overrightarrow{\boldsymbol{v}} \right) \right|, \tag{14}$$

The differential equation of heat transfer in ground domain is the energy equation for solid medium:

$$\frac{\partial (\rho\_{\rm soil} c\_{\rm soil} T\_{\rm soil})}{\partial t} = \nabla \cdot (k\_{\rm soil} \nabla T\_{\rm soil}) , \tag{15}$$

The soil humidity is considered balancing water and solid properties according to the porosity (ψ) with the following equation:

$$z\_{\rm soil} = \psi z\_{\rm liquid} + (1 - \psi) z\_{\rm solid}.\tag{16}$$

The above introduced mathematical model must be coupled with a number of boundary conditions forced on the domain. Specifically, the adopted ones are reported below:


$$\begin{array}{c} T\_{\text{ground}}(D, t) = T\_m - A\\ \cdot \exp\left[ -D \text{epth} \cdot \sqrt{\frac{\pi}{365 \cdot \alpha}} \right] \cdot \cos\left[ \frac{2\pi}{365} \cdot \left( t - t\_{\text{min}} - \frac{\text{Depth}}{2} \cdot \sqrt{\frac{265}{\pi \cdot \alpha}} \right) \right] \end{array} \tag{17}$$

• the *sun-air temperature* model is a 1st type boundary condition imposed at the top of the ground domain, as visible in Figure 1. This temperature takes into account the influence of both the incident solar radiation on the ground surface and the convective heat exchange with the external air, according to the following equation:

$$T\_{\rm sr}(\mathbf{x},0,t) = T\_{\rm air,ext}(t) + \frac{aG(t)}{\mathcal{U}\_c},\tag{18}$$

• at the inlet of the pipe, the temperature and relative humidity of the external air are imposed, whereas the inlet velocity of the air is a variable to be optimized through the investigation introduced in this paper (Figure 1).

The model is solved through the finite element method applied to the domain meshed in free triangular elements; the dimension of the mesh has been chosen according to a grid independence study that has been the object of a previous investigation [47].

The model has been experimentally validated in three case studies of horizontal ducts placed in three different countries (Algeria, Morocco, Egypt) and the maximum relative error between the experimental and the numerical data is 2.55%. Further details are reported in [47,52].

#### *2.4. Parameters of the Optimization and Operative Conditions*

Using the tool introduced in the previous section, a parametric analysis is carried out on the earth to-air heat exchanger, verifying the main design parameters such as: velocity of the air flow, length, diameter, and burial depth of the pipe affect the thermal performance of the EAHX itself, both in summer (cooling) and winter (heating) operation modes. In particular, we consider the influence of the three parameters on:


$$\varepsilon = \frac{T\_{out} - T\_{in}}{T\_{\text{ground}} - T\_{in}} \, \text{} \tag{19}$$

For this purpose, the Italian location of Naples (Lat. 40◦ 51 22.72" N; Long. 14◦ 14 47.08" E) has been considered. According to Köppen climate classification [53], which divides the areas of the globe into five main climatic zones mainly based on temperature and vegetation criteria, Naples belongs to a *Csa* climatic zone (hot-summer Mediterranean climate). The thermodynamical properties of the soil of the city of Naples are the ones typical of an Italian locality (k = 1.63 W m−<sup>1</sup> K−1; ρ = 1700 kg m−3; c = 1600 J kg−<sup>1</sup> K<sup>−</sup>1) with 37% as porosity. Through the Kusuda relation [51] and using the data coming from the ASHRAE climate database, in a previous investigation [54], the undisturbed temperature of the ground for the city of Naples was calculated (17.0 ◦C). Furthermore, always referring to ASHRAE climate data [46], the reference climate design parameters such as intensity of the solar irradiation, temperature and relative humidity of the external air, both in winter and in summer operation modes, were identified and listed in Table 2.


**Table 2.** The reference climate design parameters for the city of Naples.

The investigation has been conducted for different values of length, diameter, and burial depth of the tube. The speed of the air flow entering the pipe was also varied. Specifically, below are reported the values under which the EAHX has been tested:

$$z = \lfloor 1.5; 3.5; 5.5; 7.5; 9.5 \rfloor \text{ m},\tag{20}$$

$$L = \left[ 20; 50; 60; 80; 100; 120; 140; 160 \right] \text{m},\tag{21}$$

$$D = [0.1; 0.2; 0.3; 0.4] \text{ m},\tag{22}$$

$$v = \begin{bmatrix} 0.5; 1.0; 1.5; 2.0; 2.5 \end{bmatrix} \text{m s}^{-1} \text{ } \tag{23}$$

The regime of motion that is going to develop into the pipe, in dependence with the choice of different values of the parameters *D* and *v* is well ascertainable through the estimation of the corresponding Reynolds number, whose mathematical formulation is given by:

$$\text{Re} = \frac{vD}{v} \,\text{s} \tag{24}$$

The establishment of the fully turbulent regime inside a duct is guaranteed by Reynolds numbers higher than 104. Values lower than 2300 ensure the laminar one. In the range of numbers delimited by these two extremes, a transition zone takes place in which the fluid gradually evolves from laminar to turbulent, or it is in a metastable state. Table 3 reports the proper Reynolds number, evaluated at 20 ◦C, for each diameter and inlet air velocity couple.


**Table 3.** Reynolds numbers evaluated at 20 ◦C for each diameter and inlet velocity couple. In bold the couples for which the motion in not fully turbulent.

For almost all the couples of (*D*,*v*) tested, the motion is fully turbulent since the coupled Re are greater than 104, except the following ones (bold marked in Table 3):


which fall in the transition zone between laminar and fully developed turbulent flow. In any case, there are no cases where the flow is laminar. For this reason, according to a study in the literature [55,56], the effect of natural convection was not considered since it becomes relevant and comparable with a forced one, for Reynolds up to 600.

#### **3. Results**

Numerous series of numerical tests were carried out to develop maps of performances of the modelled earth to air heat exchanger with the final purpose of optimizing the thermal performances with respect to one working parameter at time. Specifically, to accurately evaluate the incidence of the different parameters on the thermal performances of the EAHX, a single parameter (burial depth, pipe length, diameter and air velocity) was varied at a time while the other ones remain constant.

#### *3.1. E*ff*ect of the Burial Depth*

The first parameter to be analyzed was the burial depth of the tube. The burial depth *z* was varied and the air temperature at the outlet of the pipe was evaluated while flowing, according to different velocities, in tubes with different diameters and lengths. The resulting data in all the tests show that if the pipe is buried from a depth of 1.5 m onwards, the outlet air temperature is not affected by the variation of z, whereas the diameter, length, and air velocity are. The reasons for such behavior can be found with the help of Figure 2.

**Figure 2.** Thermal penetration depth generated by the earth surface temperature.

Figure 2 shows a zoomed-in image of the investigation domain. From the figure it is clearly visible that the boundary condition (14) forced on the top of the domain, and representing both the influence of the solar radiation incident on the surface and the convective heat exchange, generates a thermal disturbance that has an influence up to 1.0 m of depth. For this reason, the tubes placed at 1.5 m and 9.5 m are not disturbed by the presence of the surface of the ground. Anyhow, as we detected and revised in previous studies published in the literature [22,38,57], the undisturbed ground temperature is generally observed at a depth of 2–4 m, but at some latitudes, depending on radiation and soil properties, it could be more than 4 m. Correspondingly, in many studies [58–60] it has been noticed that that diurnal variation (observation period of 24 h) of earth surface temperature does not penetrate more than 0.5 m, whereas, for annual variation (observation period of 365 days) it is no more than 4.0 m. The trends detected in this analysis are in accordance with the studies [58–61] because our numerical tests are run with fixed design conditions of solar radiation and external air temperature and relative humidity; thus, the simulation period is the time needed to the domain to reach thermal steady-state (about 3000–5000 s seconds), which is much less than diurnal variation.

#### *3.2. E*ff*ect of the Pipe Diameter*

The effect of the pipe diameter on the thermal performances of the EAHX is analysed in this subsection. In Figure 3 one can appreciate the trends of the air temperature at the pipe outlet at different tube diameters while the air flow enters the tube with: 0.5 m s−<sup>1</sup> (Figure 3a), 1.0 m s−<sup>1</sup> (Figure 3b), 1.5 m s−<sup>1</sup> (Figure 3c), 2.0 m s−<sup>1</sup> (Figure 3d), and 2.5 m s−<sup>1</sup> (Figure 3e). In the upper part of the figures, we show the temperature profile during summer, and, in the lower part, winter. The temperature profiles are reported for each diameter as a function of the tube length. For all the figures an analogous tendency can be noted: the air flow temperature increases/decreases through the tube length (in winter/summer). This increment/decrement is faster for the initial length of EAHX and thereafter it becomes moderate.

**Figure 3.** Outlet air temperature vs. length parametrized for diameter while the inlet air velocity is: (**a**) 0.5 m s<sup>−</sup>1; (**b**) 1.0 m s<sup>−</sup>1; (**c**) 1.5 m s<sup>−</sup>1; (**d**) 2.0 m s<sup>−</sup>1; (**e**) 2.5 m s<sup>−</sup>1.

Figure 3 reveals that, for a fixed length, the smaller the diameter is, the closer the outlet temperature is to the ground temperature (which thermodynamically constitutes an upper limit). As a consequence, at fixed length and fluid velocity, the air temperature during the summer always reaches lower values for smaller tube diameters, and the opposite trend is seen during the winter. The result is that decreasing the tube diameter at fixed air velocity increases the convective heat transfer coefficient, enhancing the heat transfer between the air flow and the tube wall. The difference is more marked corresponding to smaller tube length. As an example, at a fluid velocity of 0.5 m s−<sup>1</sup> the outlet temperature at a tube length of 60 m and diameter of 0.1 m is 17% lower/higher than that corresponding to a diameter of 0.4 m; at a length of 160 m the absolute value of the difference reduces to 3.7. This trend can be explained because, by increasing the tube length, the increment/decrement of air temperature becomes moderate because of the driving force of the heat exchange process, which is the difference between the temperature of the air and that of the undisturbed ground. This decreases, making the convective heat exchange less effective. Hence the influence of the decrease in the heat transfer coefficient is less sensitive. Observing Figure 3, one can note that the temperature profiles varying the tube diameter

are closer for low air velocities: at a speed of 0.5 m s−1, the temperature profiles corresponding to the different diameters tend to overlap at a tube length greater than 100 m. As an example, at a fluid velocity of 2.5 m s−<sup>1</sup> the outlet temperature at a tube length of 60 m and diameter of 0.1 m is 29% lower/higher than that corresponding to a diameter of 0.4 m; at a length of 160 m the absolute value of the difference reduces to 17%. Comparing these values with those found for the speed of 0.5 m s−<sup>1</sup> it is evident that the difference is more marked. Therefore, at higher air flow velocities the influence of the tube diameter is more marked. The results shown in Figure 3 are in agreement with the inherent literature. Mihalakakou et al. in their investigation [29], found that a reduction of the pipe diameter from 0.5 m to 0.25 m resulted in an improvement of the outlet air temperature drop by 1.5–2.5 ◦C during cooling mode operation. Moreover, they detected [30] that the reduction the diameter from 0.3 m to 0.2 m, during heating mode operation, reflected in an increase of outlet air temperature by 0.9–1.8 ◦C. They concluded that augmenting pipe diameter lowers convective heat transfer coefficients [62,63].

The effect of the diameter also reflects on the efficiency of the EAHX. Figure 4 shows the efficiency vs. length parametrized for diameter for different air velocities at inlet of the pipe in cooling mode. Similar trends have been observed in heating mode. In the cases plotted in the Figure 4c–e, the air flow always follows turbulent motion; here one can see that, for a fixed length, the lower the diameter, the higher the efficiency. For a tube length of L = 100 m, the efficiency is 98.6% if D = 0.1 m, with respect to 60.1% if D = 0.4 m; i.e., at this length the reduction of the diameter from 0.4 m to 0.1 m carries an increment of the efficiency of +32%, whereas the medium increment is +38.4%. The situation in Figure 4a,b is different. In this case, for D = 0.1 m, v = 0.5 m s−<sup>1</sup> and 1.0 m s−<sup>1</sup> an for D = 0.2 m and v = 0.5 m s<sup>−</sup>1, the flowing of the fluid according to a metastable state (transition between laminar and turbulent motion) results in a degradation of the efficiency for large lengths with respect to the closer turbulent points.

**Figure 4.** Efficiency of the EAHX vs. pipe length parametrized for diameter for different air velocities at inlet of the pipe: (**a**) 0.5 m s<sup>−</sup>1; (**b**) 1.0 m s<sup>−</sup>1; (**c**) 1.5 m s<sup>−</sup>1; (**d**) 2.0 m s<sup>−</sup>1; (**e**) 2.5 m s<sup>−</sup>1.

From this section an important conclusion can be drawn: to increase the volumetric flow rate, it is better to project the system with a certain number of parallel tubes with smaller diameters rather than working with only one bigger duct.

In any case, one must ensure to set *v* and *D* opportunely in order that the turbulent motion is guaranteed.

#### *3.3. E*ff*ect of the Pipe Length*

The length of the tube is a crucial parameter that influences the thermal performances of an earth to air heat exchanger. In this subsection, all the aspects inherent to this variable are analyzed.

Figure 5 reports the inlet-outlet temperature span of the EAHX operating in cooling as a function of the air velocity for different tube lengths.

**Figure 5.** Inlet-outlet temperature span of the EAHX operating in cooling mode vs. air velocity, parametrized along the length for: (**a**) D = 0.1 m; (**b**) D = 0.2 m; (**c**) D = 0.3 m; (**d**) D = 0.4 m.

Specifically, in Figure 5a–d, this effect is studied, respectively, for the following diameters: 0.1 m; 0.2 m; 0.3 m; 0.4 m. At a diameter of 0.1 m, increasing the length, one can observe how the curves approach each other and from a certain length onwards they almost overlap. At higher diameters, the curves approach but they did not overlap. There is a limit in which a further increase in length no longer leads to an improvement in thermal performances of the EAHX: this limit is called saturation length. The saturation length constitutes a knee point and it is defined as the length of the tube at which more than 90% of the global increase or decrease of the air temperature has been obtained. More generally, in correspondence with the saturation length, the thermal performances of the EAHX, such as outlet temperature, temperature span and efficiency, reach 90% of their upper limits. The aspect was plainly detected in many studies [22,41,64] where it was detected that, after a limit, further increments in pipe length do not lead to further improvements of thermal performances. The physical explanation lies in the fact that the air temperature can at most reach that of the ground; therefore, the convective heat exchange is more noticeable the greater the ΔT is. As the air flows in the pipe (and therefore as the length increases) its temperature gets closer and closer to that of the ground. This is the reason why, when the air temperature approaches the ground temperature, the heating and cooling powers associated with the heat exchanged by convection become very small. As a consequence, during the projecting of the earth to air heat exchanger the pipe lengths of the tubes

should be optimized to balance the trade-off between the investment cost for the construction and the expected thermal and energy performances.

To better understand this theme, in Figure 6 the saturation lengths with respect to inlet air velocity, parametrized for the diameter are plotted while the EAHX works in cooling mode. As a general trend, the Figures clearly show that for smaller diameters the saturation of the thermal performances is reached in correspondence of smaller lengths. Moreover, the figure shows that the incidence of the regime of fluid motion evidently emerges: if the fluid flows in the transition zone between laminar and turbulent motions, on equal diameter, the saturation lengths are higher than if the turbulent flow is fully developed. Specifically, if we collocate in turbulent motion, we can observe that, on equal velocity, the saturation length for pipes with D = 0.4 m is four times (+200%) the one established at D = 0.1 m. This is due to the above-mentioned influence of the diameter of the pipe on the thermal performances of the EAHX. Therefore, with a tube diameter of 0.1 m, a tube length lower than 80 m can be used. In contrast, with a tube diameter of 0.4 m, the tube length should always exceed 150 m.

**Figure 6.** Saturation lengths vs. inlet air velocity parametrized for the diameter while the EAHX works in cooling mode.

The pipe form factor (i.e., the ratio between the tube length and diameter) takes into account both the effect of tube length and diameter on the thermal performances of the EAHX. Figure 7 presents the effect of the pipe form factor on the air temperature variation in the EAHX during summer for different Reynolds numbers. For each Reynolds number, as the form factor increases the temperature variation also increases. This increase is faster for low form factors and then (for form factor greater than 500) becomes moderate: by increasing the form factor the length of the tube increases and the diameter decreases. This leads to a longer path over which heat transfer between the air flow and the ground can take place together with an increase of the heat transfer coefficient. At fixed form factor, the Reynolds number increases with the air velocity. For L/D less than 500, by increasing the Reynolds number the temperature difference decreases because air spends less time in the tube and thus in contact with soil. For L/D greater than 1000, all the curves merge together and therefore the influence of the Reynolds number on the thermal performance of the EAHX is negligible. This behavior can be found only for Reynolds numbers that follow in a turbulent regime. In all the explored range of L/D, the worst values of the temperature difference are those corresponding to Reynolds number that follow in the transition regime. Therefore, to optimize the thermal performance of an EAHX, a form factor L/D greater than 500 is recommended with a Reynolds number that follows in the turbulent regime. This result agrees with the results of the investigation by Sehli et al. [32].

**Figure 7.** Temperature variation in the EAHX during summer as a function of the tube form factor for different Reynolds numbers.

#### *3.4. E*ff*ect of the Air Flow Velocity*

Another design parameter playing a key role in the performance of the earth to air heat exchanger is the velocity of the air flow. As shown in Figure 8, the outlet temperatures are reported to be parametrized for velocity variation when the diameters are: (a) 0.1 m; (b) 0.2 m; (c) 0.3 m; (d) 0.4 m.

**Figure 8.** Outlet air temperature vs length parametrized for air velocity while the diameter is: (**a**) 0.1 m; (**b**) 0.2 m; (**c**) 0.3 m; (**d**) 0.4 m.

Figure 8 reveals that the outlet temperature is closer to the ground temperature (both for winter and summer) the air flow is slower. The reason of such behavior is that air flows with lower velocities result in a longer time at which the air is in contact with the tube wall: consequently, the heat transfer increases with reducing air velocity. For a fixed diameter, comparing the data carried out with the smaller and the larger velocities, a medium decrement of +2.5 ◦C in the temperature difference between the air outlet temperature and the ground temperature, is observed in correspondence with the smaller velocity. We can also note that the outlet temperature of the air approaches the undisturbed ground temperature more quickly if the velocity is smaller. For a fixed diameter, if the motion is turbulent, the most promising temperature drop is observed for v = 0.5 m s−1, as Figure 8c,d show. In Figure 8a,b the situation is slightly different where the fluid is not fully developed in turbulent: here is clearly noticeable that the best performances are detected: for v = 1.5 m s−1, if D = 0.1 m (Figure 8a); v = 1.0 m s−1, if D = 0.2 m (Figure 8b). In other words, there is a lower limit to the increment of heat transfer with velocity reduction: the latter mentioned points (D = 0.1 m, v = 1.5 m s−1; D = 0.2 m, v = 1.0 m s−1) represent these limits (one for each diameter) since, for a fixed diameter, they are the smaller (according to velocity growing) where the fluid is considered turbulent. The same behavior is seen for both working modalities of the EAHX.

The trend is also in agreement with a comparison to the literature [29,30,61,64,65] where a decrement of the temperature span was shown to lead to an increment of the flow speed. Furthermore, Niu et al. found the highest temperature drop to be at 0.5 m s−<sup>1</sup> because a slow flow is more in contact with the pipe walls, leading to a more efficient heat exchange [34].

Figure 9 reports the efficiency of the EAHX operating in cooling mode vs. length, parametrized for air velocity while: (a) D = 0.1 m; (b) D = 0.2 m; (c) D = 0.3 m; (d) D = 0.4 m. In Figure 9a,b it is clearly noticeable that reduced efficiency (9.5% medium reduction) is registered for the D and v couples that do not ensure the fully turbulent motion. In fully turbulent motion, the effect of the air velocity on the thermal performances of the EAHX reflects in an improvement of the efficiency for slower speed for a fixed diameter. The larger the diameter, the more remarkable the improvement is: for D = 0.4 m, a medium (with respect to the length) increment of +18% is registered if the velocity falls from 2.5 m s−<sup>1</sup> to 0.5 m s<sup>−</sup>1, whereas for smaller diameters the increment is progressively reduced.

The trends in Figure 10 are in agreement with the data of Figures 8 and 9 since, at fixed diameter and velocity (i.e., at fixed air flow rate), the heating power grows with the length since as the air passes through the tube its temperature get closer to the temperature of the ground and therefore inlet and outlet temperature span increases. At fixed D and L, . *Qh* augments with the air flow velocity, since the air flow increases. Similarly, . *mair* also grows with the diameter, at fixed v and L, i.e., the heating power follows the same trend.

**Figure 9.** Efficiency of the EAHX operating in cooling mode vs length, parametrized for air velocity with: (**a**) D = 0.1 m; (**b**) D = 0.2 m; (**c**) D = 0.3 m; (**d**) D = 0.4 m.

In Figure 10, the heating powers evaluated while the EAHX works in heating mode are reported as a function of the pipe length whereas the diameter is: (a) 0.1 m; (b) 0.2 m; (c) 0.3 m; (d) 0.4 m. The heating power has been evaluated according to the following relation:

$$
\dot{Q}\_h = \dot{m}\_{air} c\_p (T\_{out} - T\_{in})\_\prime \tag{25}
$$

**Figure 10.** Power of the EAHX operating in heating mode vs length, parametrized for air velocity with: (**a**) D = 0.1 m; (**b**) D = 0.2 m; (**c**) D = 0.3 m; (**d**) D = 0.4 m.

#### *3.5. E*ff*ect of the Air Temperature*

The evaluation of the effects of the burial depth, the pipe lengths, the pipe diameter and the air velocity at pipe inlet, illustrated in the previous sections, allows us to identify, with respect to our case-study, the combinations of such parameters that enhance the thermal performances of the EAHX. These results can be generalized to other cities with a Mediterranean climate. In order to optimize a horizontal single tube EAHX system placed in the city of Naples, its thermal performances must be formed by tubes with D = 0.1 m where the air flows at an inlet velocity of v = 1.5 m s<sup>−</sup>1, buried at z = 1.5 m. For the above chosen combination *(D,v,z*) the length that would satisfy the trade-off between thermal performances and installation and excavation costs is 50 m.

As next step of the investigation introduced in this paper, we studied the effect of the air temperature on the EAHX formed by one horizontal single tube and characterized by (*D*,*v*,*z*) selected above. Furthermore, we decided to leave a degree of freedom on *L* (not fixing it) to better study the effects of the thermal transient on the length.

Specifically, fixing the relative humidity rate to 50%, the EAHX was tested while the temperature at the inlet of the air assumes the following values:

$$T\_{\rm air,iu} = \left[1.9; 5; 10; 15; 20; 25; 31.9; 35\right] \,^\circ \text{C},\tag{26}$$

Figure 11 plots the outlet temperature (Figure 10a) and the specific humidity (Figure 10b) of the air flowing with 1.5 m s−<sup>1</sup> in the tube with D = 0.1 m as the diameter for different inlet temperatures.

**Figure 11.** (**a**) Outlet temperature and (**b**) specific humidity, parameterized with respect to inlet temperature of the air flowing with 1.5 m s−<sup>1</sup> in the tube with D = 0.1 m as diameter.

From Figure 11a one can observe that the highest temperature gradients are in the part of the tube defined as initial with respect to the fluid flow. In this area, for a fixed length, the gradients are larger for greater inlets ΔT between the air and the ground. Consequently, in the initial part of the tube the heat transfer is much more pronounced. Such trends correspond with the work of Niu et al. [34].

For L = 50 m, the maximum deviation detected between the air temperature and the ground temperature is −5.23% during heating operation mode (for air entering the tube with 2.3 ◦C) and +7.05% during cooling operation mode (for air entering the tube with 35 ◦C). This data confirms the choice of L = 50 m as optimal length of the tube both for cooling and heating operation modes, if the EAHX works with optimized (*D,v,z*). Therefore, if the tube length exceeds the saturation value, the inlet air temperature has a negligible influence on the thermal performance of the EAHX.

Figure 11b allows us to detect the length where the air achieves the saturation point (Φ= 100%) and consequently the water begins the condensation process. As long as the air does not saturate, the humidity ratio remains constant, then it begins to condense and the amplitude of the negative gradient is a function of the temperature, since higher temperature corresponds with larger specific

humidity. These circumstances occur only during cooling operation mode and while the inlet air temperature is 31.9 ◦C and 35 ◦C. The condensation lengths in the tube are: 32.5 m if Tair,in = 31.8 ◦C; 19.9 m if Tair,in = 35 ◦C.

The amounts of condensed water are reported in Figure 12 with respect to the length. According to the considerations made for Figure 11, for a fixed tube length, the larger the inlet temperature is, the higher the specific humidity, and the greater the amount of condensed water. Furthermore, from the knee length of the EAHX onward, amount of the condensed water at the end of the condensation process remains constant.

**Figure 12.** Amount of condensed water with length parameterized for inlet temperature of the air.

#### **4. Discussion and Conclusions**

The analysis introduced in this paper aims to explore the thermal performances of a horizontal single-tube earth to air heat exchanger to be placed in the city of Naples (Italy), by varying the most crucial parameters influencing the system. Both geometrical (pipe length, diameter, burial depth) and physical (velocity and temperature of the air flow at pipe inlet) properties are varied in an appropriate range, with the final aim to find the combinations that would optimize the operation of this geothermal system. The following conclusion can be drawn:


Nevertheless, if the tube length exceeds the saturation value, the influence of this parameter is negligible.

• If the EAHX works in cooling mode, the air entering at higher temperatures could during the flowing process, generate water condensation along the tube. For fixed length and relative humidity, the entity of such a phenomenon grows with the inlet temperature.

In conclusion, the investigation has identified the optimal combination of parameters that optimize the thermal performances of a single-tube horizontal EAHX system; consequently, an optimal air flow rate also emerges. In any case, for applications needing higher air flow rates, the optimal parameters of this optimization could be applied for each duct composing a horizontal multi-tube geothermal system. These results can be generalized to all cities with a Mediterranean climate.

**Author Contributions:** Conceptualization, A.G. and C.M.; methodology, A.G. and C.M.; formal analysis, A.G. and C.M.; investigation, A.G. and C.M.; data curation, A.G. and C.M.; writing—original draft preparation, A.G. and C.M.; writing—review and editing, A.G. and C.M.; visualization, A.G. and C.M. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **Nomenclature**

#### **Roman symbols**


#### **Greek symbols**


#### *Energies* **2020**, *13*, 6414


#### **Subscripts**


#### **References**


**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).
