Thermodynamic Analysis of Compressed Air Energy Storage (CAES) Reservoirs in Abandoned Mines Using Different Sealing Layers

Abstract

Million cubic meters from abandoned mines worldwide could be used as subsurface reservoirs for large scale energy storage systems, such as adiabatic compressed air energy storage (A-CAES). In this paper, analytical and three-dimensional CFD numerical models have been conducted to analyze the thermodynamic performance of the A-CAES reservoirs in abandoned mines during air charging and discharging processes. Unlike other research works, in which the heat transfer coefficient is considered constant during the operation time, in the present investigation a correlation based on both unsteady Reynolds and Rayleigh numbers is employed for the heat transfer coefficient in this type of application. A tunnel with a 35 cm thick concrete lining, 200 m³ of useful volume and typical operating pressures from 5 to 8 MPa were considered. Fiber-reinforced plastic (FRP) and steel were employed as sealing layers in the simulations around the fluid. Finally, the model also considers a 2.5 m thick sandstone rock mass around the concrete lining. The results obtained show significant heat flux between the pressurized air and the sealing layer and between the sealing layer and concrete lining. However, no temperature fluctuation was observed in the rock mass. The air temperature fluctuations are reduced when steel sealing layer is employed. The thermal energy balance through the sealing layer for 30 cycles, considering air mass flow rates of 0.22 kg s⁻¹ (charge) and −0.45 kg s⁻¹ (discharge), reached 1056 and 907 kWh for FRP and steel, respectively. In general, good agreements between analytical and numerical simulations were obtained.

1. Introduction

Large scale energy storage systems are required to facilitate the penetration of variable renewable energies in the electricity grids [1,2,3,4]. Underground space from abandoned mines can be used as underground reservoirs for underground pumped storage hydropower (UPSH) and compressed air energy storage (CAES) systems [5,6,7,8,9,10,11]. Pumped storage hydropower (PSH) is the most mature large-scale energy storage technology, and the round trip efficiency is typically in the range of 70–80% [12,13]. Diabatic CAES (D-CAES) is an alternative to PSH that requires lower capital costs and the round trip energy efficiency is around 40–50% [14]. D-CAES systems use natural gas to heat the compressed air in the decompression period. However, Adiabatic CAES (A-CAES) allows the storage of the thermal energy during the air compression period, avoiding the consumption of natural gas, therefore increasing the round trip energy efficiency up to 70–75% [15,16,17].

During the operation of A-CAES plants, the air is stored in the underground reservoir at high pressures during the charge period. Then, in the discharge period, the pressurized air is released, heated (stored heat), and expanded to generate electricity as the pressure within the reservoir is reduced. During the storage period, between charging and discharging, the air pressure and temperature vary depending on the storage time. The amount of stored energy depends on the reservoir volume and the thermodynamic conditions. The air temperature variations depend on the thermal conductivity of the sealing layer, concrete lining and rock mass. Therefore, the air temperature and pressure fluctuations are essential to design the reservoir volume, to select the appropriate compressor and turbines and to ensure safe operating conditions. Kushnir et al. [18] developed a study to determine the temperature and pressure variations within compressed air energy storage caverns and showed that the heat transfer reduces the temperature and pressure variations during the compression and decompression periods, leading to a higher storage capacity. Thermodynamic models were performed to determine temperature and pressure variations within adiabatic caverns of compressed air energy storage plants [19]. The results showed that the storage volume is highly dependent on the air maximum to minimum pressure ratio. Safaei and Aziz [20] carried out a thermodynamic analysis of three compressed air energy storage systems. They concluded that A-CAES with physical heat storage is the most efficient option with an exergy efficiency of 69.5%.

A pilot-scale demonstration of A-CAES was built in Switzerland [21,22]. The cavern was 120 m long and 4.9 m in diameter with a usable volume of 1942 m³. A packed-bed of rocks for thermal energy storage was located inside the cavern. Charging temperatures of 550 °C and air pressures between 0 and 0.7 MPa were employed in the test and constant heat transfer coefficients of 5 and 10 W m⁻² K⁻¹ were considered in the simulations. The estimated round trip efficiency of the pilot plant reached 63–74%. Schmidt et al. [23] analyzed the effect of the cyclic loading on the geomechanical performance of CAES systems for lined and unlined tunnels in closed mines. The cycling loading operation was simulated for 10,000 cycles considering an operating pressure within the reservoir from 4.5 to 7.5 MPa. They observed moderate deformations and small thickness of plastic zones, while an increase of the initial volume of less than 0.5% was reached.

Jiang et al. [24] carried out experimental and numerical investigations of a lined rock cavern to analyze the thermodynamic process within the underground cavern. The volume of the cavern was 28.8 m³ with 5 m length and 2.9 m in diameter. A 2 cm thick RFP was employed as sealing layer around the compressed air. The air pressure reached 6 MPa at constant mass flow rate (0.12 kg s⁻¹). Zhou et al. [25] developed a numerical simulation for the coupled thermo-mechanical performance of a lined rock cavern for CAES systems. An air pressure range from 4.5 to 6.9 MPa was considered for 100 cycles with a thickness of concrete lining of 0 and 0.1 m. They concluded that significant air temperature fluctuation was observed in sealing layer and concrete lining and no fluctuation occurred in the rock mass. An analytical study was proposed to estimate the geomechanical performance caused by air temperature and pressure variations in CAES systems [26]. Khaledy et al. [27] studied the thermo-mechanical responses of CAES systems in rock salt caverns. They concluded that the stability of the salt cavern is affected by the internal operating conditions with the air temperature increasing the volume convergence. Li et al. [28] carried out a failure study for gas storage by thermo-mechanical modelling in salt caverns. The results showed an affected area by the air and temperature operating conditions up to 10 m in the rock mass from the cavern walls. Round trip efficiencies between 70.5 and 75.1% were estimated by Barbour et al. [29] in A-CAES systems for continuous operation with packed bed thermal energy storage. The thermodynamic and geomechanical performance was studied by Rutqvist et al. [30] in lined rock caverns. They found that 97% of the thermal energy injected in the charge period could be recovered in the discharge period. Deng et al. [31] designed a new CAES systems with constant gas pressure and temperature in closed coal mines. A power output of 18 MW and a generating time of 1.76 h were obtained for a mine tunnel with 10,000 m³ of volume at 500 m depth. The ADELE project was studied in Germany to install an A-CAES plant with a storage capacity of 360 MWh and output power of 90 MW [2].

In this paper, abandoned mines are proposed as underground reservoirs for large scale energy storage systems. A 200 m³ tunnel in an abandoned coal mine was investigated as compressed air reservoir for A-CAES plants, where the ambient air is stored at high pressure. The thermodynamic response of A-CAES reservoirs was analyzed considering three solids around the pressurized air: a 20 mm thick sealing layer, a 35 cm thick concrete lining, and a 2.5 m thick sandstone rock mass. One-dimensional analytical and three-dimensional CFD numerical models were conducted considering fiber-reinforced plastic (FRP) and steel as sealing layers and a typical operational pressure from 5 to 8 MPa. In the analytical model, air mass flow rates of 0.22 and −0.45 kg s⁻¹ have been considered for 30 cycles of compression and decompression, respectively. To reduce the computational time, the numerical simulation has been conducted for five cycles with air mass flow rates of 50 and −75 kg s⁻¹. The air temperature and pressure fluctuations were estimated in the simulations. The temperature and the heat flux were also analyzed on the contact surfaces of the solids considering both FRP and steel sealing layers and the heat convection from the air to the sealing layers was calculated. Finally, the energy balance through the sealing layer was obtained during air charging and discharging.

2. Materials and Methods

2.1. Problem Statement

Underground space in abandoned mines may be used as compressed air storage systems for CAES plants. The simplified schematic diagram of the CAES system is shown in Figure 1. The compressor and turbine facilities are installed above the ground, while the compressed air reservoir is underground. The ambient air is compressed during off-peak periods and stored at high pressure in the underground reservoir (charge period). Then, when the electricity is required, the compressed air is released, heated, and expanded during peak hours in conventional gas turbines, driving a generator for power production (discharge period). The charge and discharge processes are carried out with a typical operating pressure range from 5 to 8 MPa. To reduce air leakage, two different sealing layers, FRP, and steel, have been employed in the present study. The air temperature and pressure fluctuations are estimated for both FRP and steel sealing layers.

2.2. Analytical Model

A one-dimensional analytical model has been developed in MATLAB to study the thermodynamic performance during the operation time of the CAES system in a closed mine. Figure 2 shows the scheme of the 50 m long model. A lined tunnel with a usable volume of 200 m³ and a cross section of 8 m² has been considered. A 20 mm thick sealing layer, a 35 cm thick concrete lining and a 2.5 m thick sandstone rock mass have also been considered in the model around the fluid. The evolution of the temperature (T_a), density (ρ_a) and air pressure (Pa) over time has been analyzed considering different air mass flow rates (ṁ). The heat flux through the contact surfaces, i.e., air-sealing layer (${\dot{Q}}_1)$, sealing layer-concrete (${\dot{Q}}_2)$, and concrete-sandstone (${\dot{Q}}_3)$, has been estimated. Finally, the temperature on the sealing layer wall (T₁), concrete wall (T₂), and sandstone wall (T₃) have also been estimated during the operation time considering an external temperature (T₄) of 300 K.

The energy equation has been applied following Equation (1), where $\dot{Q}$ is the net heat transfer, $\dot{W}$ is the net work done, e is the specific energy, ρ is the density in kg m⁻³, V is the reservoir volume, v is the air velocity and S is the cross section of the tunnel. In the air domain, the first term of Equation (2) represents the heat convection from the air to the sealing layer, which depends on the film coefficient of heat transfer (h), in Wm⁻² K⁻¹, the sealing layer surface in contact with the fluid (A₁) and the temperature difference between the air and the concrete wall (T_a-T₁). The second term of Equation (2) corresponds to the energy variation within the reservoir where the air pressure increases, which depends on the air mass inside the reservoir (m), the specific heat at constant volume (C_v) and the air temperature (T_a). Finally, the third term represents the input and output of energy in the control volume due to the energy provided by the air jet in the form of heat (T₀) and kinetic energy. The air mass inside the reservoir at an instant of time, t, depends on the initial air mass (m₀) and the air mass flow rate (ṁ), and is obtained by applying Equation (3).

$$\dot{Q}-\dot{W}=\frac{\partial }{\partial t}{\displaystyle \int }\rho e dV+{\displaystyle \oint }\rho \widehat{e }\left(\widehat{v}dS\right)$$

(1)

$$-h{A}_1\left(T_a-T_1\right)=\frac{d}{dt}\left(m{C}_v{T}_a\right)-\dot{m}\left(C_v{T}_0+R{T}_0+\frac{v^2}{2}\right)$$

(2)

$$m=m_0+\dot{m} t$$

(3)

After some algebra, the air temperature (T_a) is calculated by applying Equation (4). The air temperature and density increase as the air pressure increases in the reservoir, and therefore the wall temperature increases by heat convection.

$$\frac{d{T}_a}{dt}=\frac{\dot{m}}{m}\left[\frac{C_p}{C_v}{T}_0+\frac{1}{2C_v}{\left(\frac{\dot{m}}{m}L\right)}^2-T_a\right]-\frac{h2\pi r_{1}L}{C_{v}m}\left(T_a-T_1\right)$$

(4)

The temperatures on the sealing layer wall (T₁), the concrete wall (T₂) and sandstone wall (T₃) are obtained applying the equations of non-stationary heat transfer in sealing layer, concrete and sandstone domains (solid regions) considering an external temperature (T₄) of 300 K and a temperature of the air mass flow inlet (T₀) of 310 K. The heat flux reaches the sealing layer by convection. A part of this heat is transmitted to the concrete by conduction while another part increases the temperature of the sealing layer. A similar mechanism takes place for the concrete lining. Equations (5)–(7) represent the heat transmission on the sealing layer, concrete and sandstone formation, respectively.

$$\frac{d{T}_1}{dt}=\frac{2U_h}{m_s{C}_{ps}}\left(T_a-T_1\right)-\frac{2U_s}{m_s{C}_{ps}}\left(T_1-T_2\right)-\frac{d{T}_2}{dt}$$

(5)

$$\frac{d{T}_2}{dt}=\frac{2U_s}{m_c{C}_{pc}}\left(T_1-T_2\right)-\frac{2U_c}{m_c{C}_{pc}}\left(T_2-T_3\right)$$

(6)

$$\frac{d{T}_3}{dt}=\frac{2U_c}{m_{ss}{C}_{pss}}\left(T_2-T_3\right)-\frac{2U_{ss}}{m_{ss}{C}_{pss}}\left(T_3-T_4\right)$$

(7)

where U_h, U_s, U_c, and U_ss are the convection transmittance, sealing layer transmittance, concrete transmittance, and sandstone transmittance, respectively. These thermal coefficients are evaluated according to

$$U_h=h2\pi r_{1}L$$

(8)

$$U_s=\frac{2\pi K_{s}L}{Ln\left(\frac{r_2}{r_1}\right)}$$

(9)

$$U_c=\frac{2\pi K_{c}L}{Ln\left(\frac{r_3}{r_2}\right)}$$

(10)

$$U_{ss}=\frac{2\pi K_{ss}L}{Ln\left(\frac{r_4}{r_3}\right)}$$

(11)

In these expressions, h is the film coefficient of heat transfer, in Wm⁻² K⁻¹; L is the tunnel length, in m; K_s, K_c and K_ss are the thermal conductivity of sealing layer, concrete and sandstone, in Wm⁻¹ K⁻¹; Cp_s, Cp_c, and Cp_ss are the specific heat at constant pressure for the sealing layer, reinforced concrete and sandstone, in J kg⁻¹ K⁻¹; and r₁, r₂, r₃, and r₄ are the (equivalent) radius of sealing layer, concrete, sandstone and external walls, in m (Figure 2b). Note that due to the geometrical characteristics of the tunnel, the heat transfer formulation in cylindrical coordinates has been used for Equations (8)–(11).

For the evaluation of the convection film coefficient, the correlation proposed by Woodfield et al. [32] and Heath et al. [33] for the measurement of heat transfer coefficients in high-pressure vessels during gas charging was considered. In particular, a mixed (natural and forced) convection is modelled as a combination of both unsteady Reynolds and Rayleigh numbers. For the present investigation, a slight correction has been introduced for the exponent of the Rayleigh number, resulting the Nusselt number in the following expression:

$$Nu=0.56 R{e}^{0.67}+0.104 R{a}^{0.34}$$

(12)

where Re is defined as the Reynolds number of the incoming mass flow rate ($Re=4\dot{m}/2\mu \pi r_1$) and Ra is computed from the instantaneous thermal properties of the air $(Ra=g\beta \left(T_a-T_1\right)C_p{\rho }^2{L}^3/\left(\mu K\right)$), being g the gravity acceleration, β is the volumetric thermal expansion coefficient, Cp the specific heat coefficient at constant pressure, ρ the density, L the tunnel length, μ the dynamic viscosity and K the thermal conductivity. The film coefficient is thus computed as

$$h=Nu\frac{K}{L}$$

(13)

Finally, a backward Euler explicit discretization has been employed over Equations (3)–(7) to resolve the coupled system of equations. A typical time-step of 0.01 s was applied to provide an accurate temporal description of the different temperatures. The model was resolved iteratively until a maximum prescribe pressure was reached inside the reservoir.

In addition, according to Kushnir et al. [18] and Zhou et al. [25], the air has been considered as a real gas through a compressibility factor (Z) that it is computed using the Berthelot gas state equation. Finally, considering a density value uniformly distributed in the volume, the pressure value of the CAES may be obtained at any instant by applying Equations (14) and (15).

$$P_a=\left(\frac{R}{V}\right)Zm{T}_a$$

(14)

$$Z=1-\frac{9}{128}\left(\frac{P_a}{P_c}\right)\left(\frac{T_c}{T}\right)\left(\frac{6T_c^2}{T^2}-1\right)$$

(15)

where T_c and P_c the air temperature and the air pressure at critical conditions, assumed to be 132.65 K and 3.76 MPa. The contact surfaces between air-sealing layer (A₁ = 2πr₁L), sealing layer-concrete (A₂~2πr₂L), concrete-sandstone (A₃~2πr₃L), and sandstone-exterior (A₄~2πr₄L), are 359, 365, 545, and 1610 m², respectively. Note that, in order to model a complete cycle of compression (storing energy) and expansion (releasing energy) of the CAES system, the energy balance for the air domain—Equation (4)—needs to be rewritten as a function of the discharge mass flow rate (ṁ_out < 0) and the air temperature within the reservoir, according to Equation (16):

$$\frac{d{T}_a}{dt}=\frac{{\dot{m}}_{out}}{m}\left[\left(\frac{C_p}{C_v}-1\right)T_a+\frac{1}{2C_v}{\left(\frac{{\dot{m}}_{out}}{m}L\right)}^2\right]-\frac{h2\pi r_{1}L}{C_{v}m}\left(T_a-T_1\right)$$

(16)

Air mass flow rates of 0.22 kg s⁻¹ and −0.45 kg s⁻¹ have been considered in the charge and discharge periods, respectively, for both RFP and steel sealing layers. The thermal properties and the volume of air, reinforced concrete, sandstone rock mass and sealing layers considered in the model are shown in

Table 1. Thermal properties and volume of air, concrete, sealing layers and sandstone rock mass.

Material	Specific Heat (J kg−1 K−1)	Thermal Conductivity (W m−1 K−1)	Volume (m3)	Density KN m−3
Air	1006	0.0242	200.53	0.0117
Reinforced concrete	1000	1.60	195.26	23
Sandstone	711	5.00	2605	25
Sealing layer (FRP)	384	0.40	7.24	8.82
Sealing layer (Steel)	500	45	7.24	76.5

[24,25]. Note significant differences between the thermal conductivity of FRP and steel.

2.3. CFD Numerical Model

2.3.1. Model Geometry, Mesh and Boundary Conditions

A 3D numerical model of a horseshoe-shaped tunnel located inside a closed mine was conducted to simulate the compression and decompression processes (air charging and discharging) that occur in the tunnel during the operation of the CAES plant. The computational domain is a tunnel with 8 m² of cross section and 50 m in length and includes both the fluid area and the solid areas around the fluid. The configuration of the simulated tunnel is illustrated in Figure 3. To improve the geomechanical performance, a 4 m² circular cross-section has been designed for the compressed air (blue area in Figure 3). The tunnel is reinforced with a 35 cm thick concrete lining and is finished with a dead-end. To avoid thermal leakage, a 20 mm thick sealing layer is considered between the air and concrete lining. Two types of sealing materials are considered: FRP and steel. The external part of the model corresponds to the 2.5 m thick sandstone rock mass existing in the mine. The total cross section of the model is 57.13 m² (8 × 8 m) and the model volume is 2587 m³. In addition, the useful volume of air reaches 200.53 m³. The entire geometry is meshed with 2,690,038 hexahedral and tetrahedral cells. Finer mesh was defined in the sealing layer, concrete and air zones, with higher grid density in those regions where the gradients of the flow characteristics are extremely important.

Four different type of materials (air, sealing layer, concrete, and sandstone) are used to simulate the heat transfer process between the air inside the mine and the surrounding media. In addition, as in the analytical model, FRP and steel are considered as sealing layers. The air is defined as real-gas to allow the simulation of the compression process. To reproduce the concrete layer existing between the sealing layer and the sandstone, a solid material zone is defined. The sandstone and sealing layer zones are also defined as solid materials. In addition, a series of thermal properties are imposed on the sandstone, concrete and sealing zones to simulate the conduction heat transfer.

The numerical simulation of the compression/expansion process was carried out with the commercial software ANSYS Fluent^® 16.2. This code was used to solve the full 3-D Unsteady Reynolds-Averaged Navier–Stokes (URANS) equations for compressible flow. The coupling between velocity and pressure was resolved by the SIMPLE algorithm. The spatial and temporal derivatives of the governing equations for the fluid flow were calculated by means of a second-order discretization, and the pressure interpolation PREssure STaggering Option (PRESTO) scheme has been used. Turbulent closure was established by the RNG k–ε model together with standard wall functions. In addition, the resolution of the energy equation for both fluid and solid volumes was activated. A constant mass flow rate was imposed as an inlet boundary condition for each model and an interface condition was established for solid/solid and solid/fluid surfaces. To reduce the computational time, air mass flow rates of 50 kg s⁻¹ and −75 kg s⁻¹ were considered in the numerical simulations. In Section 3.3, these mass flow rates are also employed in the analytical model to carry out the comparative analysis. In addition, a constant temperature of 300 K was considered in the sandstone surface. Moreover, a no-slip shear condition was imposed on the walls. The material properties used in the numerical model have been indicated in Table 1. A fixed time step of 5 × 10⁻² s was selected for the simulations to ensure their stable convergence. The convergence criteria is based on the residual values of the solution solved in the numerical domain. A typical threshold of 10⁻⁵ was set for continuity, momentum and energy equations.

2.3.2. Mesh Sensitivity Study

A mesh sensitivity analysis has been carried out to optimize the numerical model results. Four different mesh sizes between 1.44 × 10⁶ and 3.47 × 10⁶ cells were used to simulate the air compression process from the ambient conditions to 8 MPa using FRP as sealing layer. The thermal energy balance through the sealing layer was analyzed for the first cycle of compression.

The results for the four scenarios are indicated in

Table 2. Mesh sensitivity analysis.

Mesh Cells ×106	Thermal Energy (kWh)	Error (%)	Computational Time (h)
1.44	53.96	7.14	56
2.01	56.76	2.33	78
2.69	57.72	0.68	116
3.47	58.11	-	147

. Finally, a mesh model with 2.69 × 10⁶ cells was selected to simulate the air charge and discharge periods in the underground reservoir. When a model with 3.47 × 10⁶ cells is considered, an increase in the computational time of 26.72% is obtained regarding the selected mesh size.

3. Results and Discussion

3.1. Analytical Model Results

The analytical model has been simulated for 30 cycles. The air is compressed from atmospheric conditions to 8 MPa in the first cycle. In the following cycles, the air pressure varies from 5 to 8 MPa at a stable mass flow rate considering an injected air temperature of 310 K and daily compression and decompression cycles. During an operation cycle, the air is compressed for 8 h at a mass flow rate of 0.22 kg s⁻¹ and stored for 6 h. Then, the air is released for 4 h at a mass flow rate of −0.45 kg s⁻¹ and stored for 6 h.

The variations in temperature for the air, sealing layer and concrete lining are shown in Figure 4 using FRP as sealing layer for 30 continuous cycles, which is equivalent to one month of operation. The results are also presented more in detail for the first cycle (Figure 4a) and the 30th cycle (Figure 4c). The air temperature increases up to 322 K during the first compression cycle. However, due to mixing with the injected air and heat exchange between the compressed air and the lining, the air temperature within the reservoir decreases to 307 K from the fifth cycle (five days). During the discharge period the air temperature is also stable, reaching minimum values of 294 K. The temperature of the FRP and concrete during the charge period reach 305 and 303 K, respectively, decreasing to 296 and 300 K during the discharge phase.

As shown in Figure 5, the air pressure varies from 5 to 8 MPa over the entire 30 days. Because of air temperature, the pressure decreases to 7.8 MPa during the storage period in the first cycle (Figure 5a). This reduction decreases to 7.9 MPa from the fifth cycle. The air pressure increases slightly in the storage period after the decompression stage.

The surface heat flux by convection between the compressed air and the FRP and between the FRP and concrete is shown in Figure 6. Note the thickness of the sealing layer (20 mm). Due to the temperature effect, the surface heat flux increases in the first cycle to 150 and 140 W m⁻² for FRP and concrete lining, respectively. The heat flux decreases to −90 W m⁻² in the decompression stage. As observe in Figure 6c in more detail, from second cycle the heat flux through the FRP and concrete is stable, reaching maximum values of 50 and −100 W m⁻² in the charge and discharge periods, respectively. The surface heat flux is decreasing rapidly during the storage periods. Due to the thickness of the sealing layer and the long air charging and discharging time, the surface heat flux in the sealing layer and the concrete lining is very similar.

A comparative analysis for the air temperature, pressure and surface heat flux is shown in Figure 7 for the first cycle considering both FRP and steel as sealing layers. The results are also presented for cycle 30th in Figure 8. The air temperature fluctuation between the compression and decompression periods decreases when steel is employed as sealing layer. Because of the thermal conductivity, the temperature in the contact surface between the compressed air and sealing layer is lower when steel is employed as sealing layer (Figure 8a). The reduction in air pressure in the first cycle is greater during the storage period when FRP is used as sealing layer (Figure 7b). However, as shown in Figure 8b, this reduction is very similar from second cycle. Regarding the surface heat flux, the results are very similar for both sealing layers, reaching values 50 and −100 W m⁻² in the compression and decompression stages (Figure 8c). The thermal energy balance in the 30th cycle (day 30) through the sealing layer reaches 25.27 and 20.48 kWh for FRP and steel, respectively. Although the system loses slightly more thermal energy in the compression stage when steel is used, the lining contributes with more energy to the compressed air in the decompression stage when steel is used as sealing layer (Figure 8c).

3.2. Numerical Model Results

To validate the results of the analytical model, a 3D CFD numerical model has been performed and the results are shown in Figure 9. The numerical model has been simulated for 5 cycles, considering an injected air temperature of 310 K and stable air mass flow rates of 50 kg s⁻¹ and −75 kg s⁻¹ in the charge and discharge periods, respectively. The model geometry and thermal properties of air and solids are the same as those used in the analytic model. As presented in Figure 9a, the maximum air temperature decreases from 410 to 340 K between the first and fifth cycle. The air temperature fluctuation between air charging and discharging reaches 45 K in the fifth cycle. Likewise, the FRP temperature decreases from 385 to 322 K. However, the concrete lining temperature is more constant throughout the process. The surface heat flux for FRP and concrete is shown in Figure 9c. Due to convective effects, the heat flux increases in FRP during the first cycle. Then, from third cycle is stabilized in a maximum value of 1000 W m⁻² in the compression period and −1000 W m⁻² in the decompression period. The numerical model results using steel as sealing layer are shown for five cycles in Figure 10.

As in the analytical model, the air temperature fluctuation between air compression and decompression periods is reduced down to 41 K when steel is employed as sealing layer (Figure 10a). In general, the air temperature is much lower than the previous scenario. The surface heat flux is presented in Figure 10c, varying in the steel surface between 2000 W m⁻² and −2000 W m⁻² for the charge and discharge periods, respectively. Due to the high thermal conductivity (45 W m⁻¹ K⁻¹), the heat flux through the sealing layer is greater, and therefore the temperature in the steel surface is lower. The temperature in the steel coincides with the temperature on the concrete surface, reaching a maximum value of 320 K in the first cycle.

The distribution of the air temperature within de reservoir and the detail of the temperature in the sealing layer are shown in Figure 11 after the decompression period at 5 MPa in the first cycle, using FRP as sealing layer.

Due to the low thermal conductivity of FRP, the temperature reaches 390 K at the top of the reservoir. However, the temperature in the FRP decreases to 330 K at the bottom. A minimal increase in temperature is observed in the concrete lining. The distribution of the air temperature within de reservoir and the detail of the temperature in the sealing layer are shown in Figure 12 after the decompression period at 5 MPa in the first cycle, using steel as sealing layer. The temperature in the steel reaches 330 K at the top of reservoir. In this scenario, a temperature increase up to 320 K has been obtained for the concrete lining.

The distribution of the surface heat flux around the compressed air is shown in Figure 13 after the first compression cycle at 8 MPa for both FRP and steel sealing layers. The surface heat flux reaches a maximum value of 7500 W m⁻² in the steel surface, while in the FRP surface the maximum value is 3000 W m⁻².

3.3. Comparative Analysis

A comparative analysis between the analytical and numerical models have been carried out. Air mass flow rates of 50 kg s⁻¹ and −75 kg s⁻¹ were also considered in the analytical model to carry out the comparative study. As indicated previously, although both analytical and numerical models use the same model geometry, the air mass flow rates of the numerical model are higher to reduce the computational time. Figure 14 shows a comparative analysis for five cycles between analytical and CFD results using FRP as sealing layer. A comparison between FRP and steel is shown afterwards in Figure 15 for one cycle. The obtained heat transfer coefficient is also indicated in Figure 15b. In general, good agreements have been obtained between both analytical and numerical simulations.

To justify the accuracy of the simulations, the results obtained in the analytical model have been compared with experimental and numerical models. Figure 16 shows the comparison of results obtained for air temperature (Figure 16a) and pressure (Figure 16b) during the first cycle with an experimental model developed by Jiang et al. [24]. Dot points, corresponding to the experimental values, are compared with the results obtained in solid blue lines in Figure 16. In addition, Figure 17 shows the comparison of results obtained for air and wall temperatures (Figure 17a) and air pressure (Figure 17b) during the first cycle with a numerical model conducted by Zhou et al. [25]. The analytical model (solid lines) has been employed to reproduce the numerical model considered by Zhou et al. (dashed lines), reporting good accuracy as shown in Figure 17.

4. Conclusions

Analytical and numerical models were established to investigate the thermodynamic performance of an underground reservoir located in abandoned mines for A-CAES plants. Analyzing air temperature and pressure fluctuations during the operation time in A-CAES plants is essential to design the underground reservoir volume and the turbomachinery. Typical operating pressures from 5 to 8 MPa were considered in the air charging and discharging periods. A 20 mm thick sealing layer, 35 cm thick concrete lining and 2.5 m thick sandstone rock mass have been considered around the compressed air. In addition, FRP and steel have been employed as sealing layer. Unlike other research works, in which the heat transfer coefficient is considered constant during the operation time, in the present investigation a correlation based on both unsteady Reynolds and Rayleigh numbers is employed for the heat transfer coefficient.

The results obtained show greater variations in air temperature between the air compression and decompression when FRP is used as sealing layer. Thus, significant temperature variations in the sealing layer and only a 15 cm thickness of the concrete lining is affected by a temperature rise during operation when steel is employed as sealing layer. The volume of concrete affected is reduced when FRP is used as sealing layer around the fluid. Regarding the sandstone rock mass, no temperature fluctuation was observed in the simulations. In addition, the heat flux increases and the air temperature within the reservoir is lower when steel is employed as sealing layer.

Finally, the thermal energy balance in the 30th cycle through the sealing layer considering air mass flow rates of 0.22 and −0.45 kg s⁻¹ reached 25.27 and 20.48 kWh for FRP and steel, respectively. In general, good agreements were obtained between analytical and numerical simulations. Unlike a 1D analytical model, in the 3D CFD numerical model it is possible to analyze the distribution of the thermodynamic responses in the entire domain, assuming a significant advantage to design the reservoir.