Numerical Modeling of the Thermal Behavior of Subsea Hydro-Pneumatic Energy Storage Accumulators Using Air and CO₂

Abstract

This paper numerically models the thermal performance of offshore hydro-pneumatic energy storage (HPES) systems composed of a subsea accumulator pre-charged with a compressed gas. A time-marching numerical approach combining the first law of thermodynamics with heat transfer equations is used to investigate the influence of replacing air within an HPES system with carbon dioxide (CO₂). The latter is able to experience a phase change (gas–liquid–gas) during the storage cycle in typical subsea temperatures when limiting the peak operating pressure below the critical point. The influences of integrating a piston and an inner liner within the accumulator to mitigate issues related to gas dissolution in seawater and corrosion are explored. It is found that the energy storage capacity of subsea HPES accumulators increases substantially when CO₂ is used as the compressible fluid in lieu of air, irrespective of the accumulator set up. It is also noted that the length-to-diameter ratio of the accumulator has a considerable influence on the round-trip thermal efficiency for both air- and CO₂-based accumulators. Another factor influencing the round-trip thermal efficiency is the presence of the inner liner. Moreover, the CO₂-based HPES system yields a lower round-trip thermal efficiency over that of air.

1. Introduction

The United Nations established 17 Sustainable Development Goals (SDGs) in 2015, with the 13th highlighting the need to take action to combat climate change and its impact by reducing greenhouse gases (GHG) [1]. One of the main greenhouse gases is CO₂, which is a bi-product of the combustion of fossil fuels [2]. A reduction in fossil fuels is possible by shifting the electricity energy market towards renewable energy sources (RES). Although RES technologies such as wind and solar power are now mature enough to effectively combat climate change, the intermittency in supply presents severe engineering challenges when large RES farms are interconnected to electricity grids. Such challenges of RES could be mitigated through the integration of energy storage systems (ESS). Some of the main issues of ESS are their initial capital required and the cost of operation. Yet, grid-scale ESS are increasing in their development and deployment due to technology advancements and policy actions.

By the end of 2020, the global installed capacity of energy storage had reached 191.1 GW [3]. ESS are usually classified as follows: mechanical, thermal, electrochemical and chemical storage systems. The greatest share of stationary energy storage is covered by pumped storage hydro (PSH), accounting for over 99% of the global ESS installed capacity. Due to geographical restrictions, it is not always possible to install PSH systems. An alternative to this is compressed air energy storage (CAES). Compressed air energy storage systems have been around since the 1940s, but their potential was significantly studied in the 1960s and 1970s [4].

One of the ESS benefiting from recent technological advancements is the Li-ion battery. Main advantages over other ESS include a high round-trip efficiency and a high specific energy capacity. However, they still suffer from a number of drawbacks related to lifetime, dependence on rare materials, recyclability and fire safety [5].

As the name implies, CAES utilizes compressed air to store energy produced either from renewable resources or fossil fuels to be converted into electrical energy during periods of high energy demand [6]. The air is stored in a pressure containment system which could either be natural underground caverns, a salt dome or pressure vessel. The pressure vessels could be placed on land or offshore, either on a platform or subsea [7]. Deploying CAES systems offshore, co-located with RES generation plants such as offshore wind farms, avoids the additional space required on land to store energy. It also offers the opportunity to decouple the power transmission system from RES generation, enabling cost reductions in the electric cables used to transit the electricity to shore.

In various CAES configurations, energy from RES is utilized to compress ambient air through a compressor into the pressure containment. Then, air is expanded back to atmospheric conditions through a turbine coupled to an electrical generator connected to the electrical grid. Different CAES systems exist, including diabatic, adiabatic and isothermal. Diabatic CAES systems allow thermal energy to be dissipated to the surrounding during compression, while fuel is used to reheat the expanding air. Adiabatic CAES systems utilize thermal energy storage (TES) to store thermal energy generated during compression for later use during expansion. Isothermal CAES attempts to obtain isothermal expansion and compression [8]. An alternative form of CAES technology is hydro-pneumatic energy storage (HPES) systems, which in reality combine PSH and CAES by replacing the air compressor and turbine with hydraulic machinery. The pressure vessel is initially pre-charged with compressed air and the gas is compressed and expanded via a liquid piston [9,10,11].

The use of hydraulic pumps and turbines instead of compressors and expanders enables a higher power density to be achieved. This reduces the footprint of the machinery required, which is an important attribute to consider when installing storage systems offshore. Furthermore, HPES can achieve quasi-isothermal conditions when the compression ratio is limited [12].

D. Buhagiar et al. [13,14] presented an offshore-based HPES system known as the floating liquid-piston accumulator using seawater under compression (FLASC). The concept uses air as the compressible fluid and exploits the excellent heat source/sink properties of the surrounding seawater in the marine environment to approach isothermal conditions during the compression and expansion stages of the HPES cycle [13]. Tests on a scaled prototype at sea have demonstrated that the compressed air responded very quickly to the surrounding water temperature, proving that isothermal conditions can be attained with thermal efficiencies exceeding 94% [13,14].

Recent research has been exploring the possibility of increasing the storage density of HPES systems by using CO₂ as a compressible fluid storing the energy. Abuheiba et al. [15] studied the influence of replacing air with a mixture of CO₂ and nitrogen (N₂) within an onshore pressure vessel HPES system through numerical modeling. Within this study, mineral oil was used as the liquid piston since both gases are immiscible to it [10,15]. From their study, it was concluded that making use of 88% CO₂ and 12% N₂ by volume and an initial pressure of 73 bar led to increasing the energy density by 56.8% over that of pure N₂. Moreover, in their study it was noted that compressing the mixture from 30 bar to 128 bar led to 5.26 MJ/m³ to be stored, while the maximum energy stored of N₂ was only of 5.09 MJ/m³ [15]. Cutajar et al. [16] investigated the potential increase in the storage density of subsea HPES systems using CO₂ instead of air by assuming ideal isothermal conditions and following the pressure-enthalpy chart of the gas. It was concluded that replacing air with CO₂ can theoretically more than double the storage density of air-based HPES systems; however, this depends on the dryness fraction of the CO₂ at the end of the charging cycle, as well as on the sea temperature and depth [16].

While the works of Abuheiba et al. [15] and Cutajar et al. [16] have independently confirmed the potential of CO₂ in substantially increasing the energy storage density of HPES systems, it is still unclear whether high thermal efficiencies similar to those reported by Buhagiar et al. [13,14] can still be achieved in the marine environment when replacing air with CO₂. This paper aims to address this knowledge gap using time-resolved numerical modeling for the thermal behavior of a subsea accumulator undergoing an entire storage cycle, including both charging and discharging of hydraulic energy. Simulations are carried out for both air-based and CO₂-based accumulators so that a direct comparison of the HPES systems using the two different compressible fluids could be made. This paper is organized as follows: Section 2 delves into the background theory used to model the systems for both gases, Section 3 discusses the code developed and its flow, Section 4 introduces the numerical simulations and system parameters, Section 5 presents and discusses the results of the two gases for the numerical simulations and finally Section 6 highlights the conclusions derived from this study.

2. Background Theory

This section describes the subsea HPES accumulator considered in the present study. It also delves into the theory behind the mathematical modeling of the accumulator utilizing both air and CO₂.

2.1. Modeled Accumulator

The subsea HPES accumulator, presented in Figure 1, was assumed to consist of a long steel pipeline laid on the seabed. The pipeline was initially pre-charged with air or CO₂. Inside the steel pipeline, an inner high-density polyethylene (HDPE) liner was installed to prevent steel corrosion induced by seawater or CO₂. A piston was utilized to prevent any gas from dissolving into the seawater. The piston was assumed to be similar to a pipeline inspection gauge (PIG), which is well-known technology in the pipeline industry. The PIG also aided in avoiding water from being entrapped in the pipeline when undulations were present on the seabed.

The gas experienced four main stages: charging, discharging and two hold stages, with the typical pressure variation against time being presented in Figure 2. During the charging stage, a hydraulic pump was used to inject seawater into the accumulator during periods of excess energy supply, thereby compressing the gas and storing energy. During the hold stages, both the first and second, no work was completed on or by the system as heat is lost or gained from the surrounding walls of the system and the seawater. As the system discharged, the gas within the pressure chamber was allowed to expand to the initial pressure, causing the seawater to drive the recovery turbine and generate electricity. The large surface area of the pipeline enabled the heat exchange between the gas and the surrounding seawater in the subsea environment. In this way, the seawater acts as a heat sink during the charging stage and heat source during the discharging stage, to maintain quasi-isothermal conditions.

Prior to carrying out numerical modeling, the following assumptions were considered: (1) the pipeline length-to-diameter ratio is relatively large such that the heat transfer between the gas to the closed end of the pipeline and the gas to the PIG are negligible; (2) the contact point of the pipeline with the seabed is negligibly small such that the outer wall of the steel pipeline is completely in contact with the surrounding seawater; (3) the inner diameter of the outer pipeline and the outer diameter of the inner liner are in direct contact with one another; (4) frictional losses in the seawater flowing in and out of the accumulator are negligible; and (5) the mass of the PIG is ignored.

2.2. Model for Air-Based Accumulator

The thermal behavior of the proposed HPES accumulator was observed with respect to time, while modeling the system using a discrete time-marching approach. The change in internal energy ($\Delta U$) between each time step of the system could be expressed in terms of the work completed ($W$) by/on the compressed gas and the heat emitted/absorbed ($Q$) by the compressed gas. The latter relation is expressed by Equation 1, also referred to as the first law of thermodynamics.

$$\Delta U=W-Q$$

(1)

The work completed by/on the gas over each incremental time step was calculated from Equation (2). The latter equation was dependent on the change in volume ($\Delta V$) of the compression chamber resulting from the pumped seawater or expanding gas and the pressure of the chamber ($p_{ch}$). The change in volume, calculated from Equation (3), was dependent on the volume flow rate ($\dot{V}$). As shown from Equation (4), the volume flow rate was dependent on the pump/turbine pressure ($p_{p/t}$) which was calculated from Equation (5). Thus, the pump pressure was dependent on the chamber pressure ($p_{ch}$) of the previous time step and the frictional force imposed by the PIG ($F_f$). A positive frictional force was considered during charging and a negative force during discharging.

$$W=p_{ch\left(t-i\right)}\times \Delta V$$

(2)

$$\Delta V=\dot{V_t}\times \Delta t$$

(3)

$$\dot{V}=\frac{E_{c-d}}{t_{c-d}\times p_{p/t \left(c-d\right)}}$$

(4)

$$p_{p/t \left(c-d\right)}=\frac{\left(p_{ch \left(t-i\right)}\times A_{pist}\right)\pm F_f}{A_{pist}}$$

(5)

In the present study, the real hydraulic performance characteristics of the hydraulic pump and the recovery turbine were not modeled. A black box model approach was taken in which the power transferred by the compressible fluid to the pump and that recovered by the turbine were simply taken to the product of the pressure across the machines and the volume flow rate. X. Zhu et al. [17] presented the dependence of the frictional force between the PIG disc and the steel pipeline as a result of the percentage interreference ($\delta$) and the thickness-to-diameter ratio ($\xi$). In their work, the influence of the two mentioned parameters on the PIG disc frictional force was studied based on the differential pressure between the two surfaces of the disc. The data presented for the two parameters by X. Zhu et al. [17] were digitized to derive a simplified relation (Equation (6)). This relation calculated the frictional force experienced by the PIG disc while sliding in the steel pipeline or inner liner based on $\delta$ and $\xi$, while taking into account the static and kinetic friction. Since the digitized relations from X. Zhu et. al. [17] only yield friction values related to the diameter of the pipeline utilized for their experimentation, an additional correction was introduced in Equation (6) to convert the friction results from their diameter to the diameters considered in this study.

$$F_f=\left(C_{\delta }\times C_{\xi }\right)\times {\mu }_{s/k}\times \left(\frac{d_i}{0.9868}\right)$$

(6)

The heat being exchanged between the gas and the surrounding seawater was modeled through a number of resistances in series, as presented in Figure 3. The resistances under study only considered convection and conduction. The nodes of the inner liner and outer steel pipeline were placed in the middle of the respective parts. Equations (7) to (12) represent the resistances indicated in Figure 3.

$$R_{conv in}=\frac{1}{h_{in}{A}_i}$$

(7)

$$R_{cond 1}=\frac{\ln \frac{r_x}{r_i}}{2\pi k_{inl}{L}_{ch}}$$

(8)

$$R_{cond 2}=\frac{\ln \frac{r_m}{r_x}}{2\pi k_{inl}{L}_{ch}}$$

(9)

$$R_{cond 3}=\frac{\ln \frac{r_y}{r_m}}{2\pi k_{osp}{L}_{ch}}$$

(10)

$$R_{cond 4}=\frac{\ln \frac{r_o}{r_y}}{2\pi k_{osp}{L}_{ch}}$$

(11)

$$R_{conv in}=\frac{1}{h_{out}{A}_o}$$

(12)

For the external free convection heat transfer coefficient over the horizontal steel pipe, the Churchill and Chu model was used to calculate the average Nusselt number, while assuming uniform heating throughout the cylindrical wall and uniform wall temperature. Its equation is presented in [18]. Thus, the result of the average Nusselt number was used to calculate the external convective heat transfer coefficient from Equation (13).

$$h_{out}=Nu\times \frac{k_{sw}}{d_o}$$

(13)

To compute the internal convective heat transfer coefficient within the chamber, the system was considered as a non-flow system during charging and discharging stages for both air and CO₂. Few researchers carried out experimentation to determine the internal heat transfer correlations for non-flow systems under compression. The available literature only considers vertical pipes with liquid pistons, as, for example, presented by Neu et al. [19,20]. Therefore, for this study, correlations for flow systems were considered. The Gnielinski correlation, presented in [21], was considered for turbulent flow, while for a laminar flow a Nusselt number of 4.36 was used for air. Then, the internal convective heat transfer coefficient ($h_{in}$) was calculated from Equation (14) while substituting the Nusselt number depending on the type of flow.

$$h_{in}=Nu\times \frac{k_g}{d_i}$$

(14)

During both hold stages, $h_{in}$ was set to a fixed value, taking into consideration typical values found in the literature. The value was set to 10 W/m²·K throughout the study [22,23,24].

The resistances between nodes 1–2, 2–3 and 3–4 acted in series to one another; therefore, they could be added to each other as presented in the groups of Figure 3. Moreover, the heat transfer between each node is presented by Equations (15) to (17), where the temperature of each node was that of the previous time step.

$$Q_{1-2}=\frac{T_1-T_2}{R_{1-2}}\times {\Delta }_t$$

(15)

$$Q_{2-3}=\frac{T_2-T_3}{R_{2-3}}\times {\Delta }_t$$

(16)

$$Q_{3-4}=\frac{T_3-T_4}{R_{1-2}}\times {\Delta }_t$$

(17)

Equation (1) was modified to Equation (18) to obtain the temperature of the gas at the current time step, for air and CO₂ in the gaseous phase. Since in both cases no internal energy (${\dot{E}}_{gen}$) was being generated, it reached zero, while the energy input (${\dot{E}}_{in}$) and output (${\dot{E}}_{out }$) from the nodes remained for both nodes 2 and 3. The energy between nodes 1 and 2 acted as an input (${\dot{E}}_{in}$) into node 2, while the energy transfer between 2 and 3 acted as an output (${\dot{E}}_{out }$), depicted in Figure 3. Similarly, the energy between nodes 2 and 3 acted as an input (${\dot{E}}_{in}$) into node 3, while the energy transfer between 3 and 4 acted as an output (${\dot{E}}_{out }$). The energy equation, presented by Equation (19), transformed to Equation (20) to estimate the temperature of the HDPE inner liner. Similarly, Equation (21) was used to determine the temperature of the outer steel pipeline.

$$\begin{array}{l}m{c}_v\left(T_t-T_{t-i}\right)=W-Q_{1-2}\\ T_{1 \left(t\right)}=\frac{\left(W-Q_{1-2}\right)}{m{c}_v}+T_{1 \left(t-i\right)} \end{array}$$

(18)

$${\dot{E}}_{in}-{\dot{E}}_{out}+{\dot{E}}_{gen}=\rho cV\left(\frac{dT}{dt}\right)$$

(19)

$$T_{2 \left(t\right)}=\left(\frac{Q_{1-2}-Q_{2-3}}{{\rho }_{lnl}\times V_{inl}\times c_{inl}}\right)+T_{2 \left(t-i\right)}$$

(20)

$$T_{3 t}=\left(\frac{Q_{2-3}-Q_{3-4}}{{\rho }_{osp}\times V_{osp}\times c_{osp}}\right)+T_{3 \left(t-i\right)}$$

(21)

Air was assumed to behave as a real gas during the HPES cycle; therefore, the Beattie–Bridgeman equation, presented in [25], was used to calculate the pressure of the gas at the current time step for the new temperature and volume.

2.3. Model for Carbon Dioxide-Based Accumulator

The theory described for the work conducted acting on air in Section 2.2 applies for CO₂, irrespective of its phase. The resistances of the heat transfer, presented in Figure 3, between the gas and the surrounding seawater apply for CO₂, as well. Similar to air, the Churchill and Chu correlation was utilized for the external convective heat transfer coefficient. The following describes briefly the models used for CO₂ accumulator simulations only.

The internal convective heat transfer coefficient ($h_{in}$) for CO₂ depends on whether it is in a complete gaseous state or a liquid–vapor state, while the latter state is further divided into boiling or condensation. When a fluid changes its phase from liquid to vapor, it is referred to as evaporating or boiling, while when the vapor is transforming to liquid it is said to be condensing. Similar to air, the $h_{in}$ for CO₂ in its gaseous state was determined through the Gnielinski correlation, as it exhibited turbulent flow, while a Nusselt number of 4.36 was used for laminar. Carbon dioxide was expected to experience condensation during the charging and evaporation during the discharging under subsea conditions when the pressure was retained below the critical point. The determination of $h_{in}$ during condensation was carried out through the Thome et al. [26] model. For evaporation, $h_{in}$ was determined by initially calculating the Nusselt number from the Fang et al. study [27].

During both hold stages, the internal convective heat transfer coefficient was set to 10 W/m²·K when CO₂ was in its gaseous state, while for the two-phase region it was set to 1000 W/m²·K [22,23,24].

The temperature of CO₂ in the two-phase region was estimated using a different approach. The specific internal energy of CO₂ in the two-phase region could be described by Equation (22). Equation (23) provides a modified version of the first law of thermodynamics, where the specific internal energy of the current time step was in terms of the specific work completed, the heat transfer and internal energy of the previous time step. The specific enthalpy of the current time step for Equation (22) was calculated through Equation (24).

$$u_{\left(t\right)}=h_{\left(t\right)}-p_{ch\left(t\right)}{v}_{spec\left(t\right)}$$

(22)

$$u_{eqn\left(t\right)}=u_{\left(t-i\right)}+\frac{W_{\left(t\right)}}{m}-\frac{Q_{\left(t\right)}}{m}$$

(23)

$$h_{\left(t\right)}=h_{L\left(t\right)}+x\left(h_{V\left(t\right)}-h_{L\left(t\right)}\right)$$

(24)

In theory, if the temperature used to find the specific enthalpies at the liquid and vapor state corresponds to the actual temperature of the current time step, there should be no difference between the specific internal energy calculated from Equation (22) to that of Equation (23). On the other hand, if a difference is present, the parameters of Equation (22) are altered to a new temperature by increasing or decreasing it by a specified increment, $T_{inc}$, accordingly. This process was repeated until the difference between the results of the two mentioned equations, calculated by Equation (25), reached the acceptable minimal value. The temperature satisfying the condition was set as the current time step temperature. The temperatures of nodes 2 and 3 were calculated through Equations (20) and (21).

$$u_{dif\left(t\right)}=\left|\frac{u_{eqn\left(t\right)}-u_{\left(t\right)}}{u_{eqn\left(t\right)}}\right|\times 100$$

(25)

During the liquid–vapor state, the pressure was obtained from a Python^® library defined as a CoolProp (CP) library using Equation (26) according to the final temperature, which satisfied the desired percentage difference of Equation (25) [28]. Since for the liquid–vapor state the library was used, it was ideal to make use of the library for the gaseous state, as well.

$$p_{ch\left(t\right)}=CP.PropsSI\left({}^{\prime }{P}^{\prime },{}_{}{}^{\prime }T{}^{\prime }, T,{}^{\prime }{Q}^{\prime }, 1, R744\right)$$

(26)

2.4. Round-Trip Thermal Efficiency

The round-trip thermal efficiency of the system is defined as the ratio between the hydraulic energy extracted to that injected into the accumulator across an entire cycle. The energy stored during the charging stage was calculated by Equation (27), while the energy extracted during the discharging stage was equated by Equation (28). The thermal efficiency ($TE$) of the system was equated through Equation (29).

$$E_{st}={{\displaystyle \int }}_{V_1}^{V_2}pdV={\displaystyle \sum }_{i=1}^{N_{st}}\left(p_{ch\left(c-d\right)}-p_h\right)\times \left({\dot{V}}_{in\left(t\right)}\times \Delta t\right)$$

(27)

$$E_{dis}={\displaystyle \sum }_{i=1}^{N_{dis}}\left(p_{ch\left(c-d\right)}-p_h\right)\times \left({\dot{V}}_{out\left(i\right)}\times \Delta t\right)$$

(28)

$$TE=E_{dis}/E_{ch}$$

(29)

3. Implementation of Numerical Models

The mathematical equations described in Section 2 above were coded in Python^® V3.8 language using the Spyder Integrated Development Environment (IDE). Two separate codes were developed to simulate the behavior of the air-based and CO₂-based HPES systems. Both codes developed followed the flow presented in Figure 4. It was mainly divided into the four stages of an HPES cycle, being the charging stage, first hold stage, discharging stage and finally the second hold stage. At each stage, the code saved a selected set of results as an excel file for later analysis. The air code will be referred to as CAIR (code of air), while that of CO₂ will be referred to as C²O₂ (code of CO₂) in future sections.

4. Methodology

The modeled baseline HPES accumulator consisted of a 36-inch pipeline constructed from X70 grade steel. The external and internal diameters were taken to be equal to 0.914 m and 0.884 m, respectively. The internal volume was set to 1782.72 m³, corresponding to a total pipeline length of 2905 m. These dimensions were determined through the in-house code Smart PVB developed by Cutajar et al. [29]. The accumulator was designed for a storage capacity of 1 MWhr system when operating with air at a peak working pressure of 60 bar. Cutajar et al. [29] had determined that in the case of air-based HPES accumulators, a pressure ratio of 2.5 resulted in the lowest steel requirements per unit of storage capacity. Thus, the pre-charge pressure in this study for the air-based system was set to 24 bar. The simulations for the CO₂-based accumulator were limited to a peak working pressure of 60 bar and a dryness fraction ($x$) of 0.2 during charging, while returning to the initial conditions during discharging. This means that the charging stage for CO₂ was terminated when the pressure reached 60 bar or the dryness fraction reached 0.2, whichever value was reached first.

Jim Schmitz [30] suggested a value of 41 for the ratio of the outer diameter to the wall thickness of the HDPE liner in the case of pipelines used for transporting CO₂. The outer diameter of the HDPE inner liner was assumed to be equal to the inner diameter of the outer steel pipeline. Hence, the outer and inner diameters of the inner liner were set to 0.884 m and 0.841 m, respectively.

The default power values of the hydraulic pump and energy recovery turbine were both set to 250 kW, implying a charging and discharging time of approximately 4 h for the baseline 1 MWhr accumulator. The duration of both hold stages (post-charging and post-discharging) was also set to 4 h. The system was set to operate at a sea depth of 30 m, while having the surrounding seawater temperature at 283 K [29]. The properties of the seawater were considered for a fixed temperature of 283 K and a salinity of 35 g/kg.

The PIG disc dimensions were determined for a percentage interference ratio of 4.5 and a thickness-to-diameter ratio of 0.022, while utilizing the equations presented in the work of X. Zhu et al. [17]. The disc thicknesses with and without the inner liner are presented in

Table 1. Pipeline inspection gauge disc diameter, thickness, static and kinetic friction for the accumulators with and without the inner liner.

Description	Values	Units
Without Inner Liner
Diameter ($D_s$)	0.924	m
Thickness ($th$)	0.0195	m
Static Friction (${\mu }_s$) [17]	0.60	/
Kinetic Friction (${\mu }_k$) [17]	0.50	/
With Inner Liner
Diameter ($D_s$)	0.879	m
Thickness ($th$)	0.0185	m
Static Friction (${\mu }_s$) [31]	0.20	/
Kinetic Friction (${\mu }_k$) [31]	0.12	/

. This table also provides the static and kinetic friction between the PIG disc and the inner pipeline surface. The mentioned parameters presented above and in Table 1 are the default parameters for the baseline accumulator, for both air and CO₂.

A series of numerical experiments were initially carried out to test the codes for ideal isothermal and adiabatic conditions when the system did not include the inner liner nor the PIG. Moreover, numerical experiments were carried out to investigate the influence of the following design parameters on the thermal performance: the accumulator pipeline length-to-diameter ratio, the use of a PIG which introduces frictional losses and the use of the inner liner which increases the thermal resistances for the heat exchange between the compressed fluid and the surrounding seawater. The simulations carried out are summarized in

Table 2. Duration of the different simulations considered for this study.

Simulation Number	Simulation Description for Air and CO2	Duration (Minutes)
		Air	CO2
1.	Isothermal conditions	0.35	4.46
2.	Adiabatic conditions	0.70	1.74
3.	Baseline accumulator	0.81	28.21
4.	Influence of aspect ratio (36, 48, 60 and 84 inches)	0.81, 0.74, 0.71 and 0.64	28.21, 55.29, 71.29 and 31.88
5.	Influence of pipeline inspection gauge	0.77	39.34
6.	Influence of inner liner	0.75	53.15
7.	Combined effect of both inner liner and pipeline inspection gauge	0.75	83.74

. A verification process was carried out for both the air-based and the CO2-based HPES codes to determine the time step size. A time step of 0.4 s yielded negligible numerical errors, while retaining the computational time within acceptable limits. The dimensions of the different diameter pipelines, being of 48, 60 and 84 inches while maintaining the same volume of the 36-inch pipeline, were also determined through Smart PVB, as indicated in

Table 3. Different accumulator pipeline dimensions considered, apart from the baseline one.

Pipeline Diameter (inch)	Outer Diameter (m)	Inner Diameter (m)	Pipeline Length (m)
48	1.22	1.19	1603
60	1.52	1.48	1036
84	2.13	2.08	525

. Additionally, Table 3 indicates the length of the pipeline for the different diameters. The different design parameters are discussed in Section 5.

The simulations considered for this study were carried out on two desktop computers, one powered by an AMD Ryzen 7 5800H (3.20 GHZ) and 16 GB of RAM, and the other powered by an Intel Core i7-7800X and 48 GB of ram. Both desktop computers ran on Microsoft Windows 10, while having an SSD as the hard disc type. For this study, different conditions were simulated to be discussed and analyzed in Section 5. The computational time for each numerical simulation is shown in Table 2. From the computational time, it is evident that C²O₂ took longer to yield results than CAIR.

5. Results and Discussion

This section compares the predictions for the thermal performance of air-based and CO₂-based accumulators for the conditions listed in Table 2.

5.1. High Convective Heat Transfer Coefficient and No-Heat Transfer Study

For both air- and CO₂- based HPES systems, the internal and external convective heat transfer coefficients are initially set to relatively high values to test whether the codes are able to reproduce isothermal conditions (simulation 1, Table 2). The internal and external values for the air-based and CO₂-based accumulators are both set to 10,000 W/m²·K and 100,000 W/m²·K, respectively. A round-trip thermal efficiency of 99.98% is predicted for the air-based system while storing 1.0073 MWhr and recovering 1.0071 MWhr. Similarly, the CO₂-based HPES system yields an efficiency of 99.33% while storing 1.4461 MWhr and recovering 1.4365 MWhr. For the air system, the polytropic index is calculated through Equation (30) for both the charging and discharging stages, yielding to an index relatively close to 1, and thus confirming that the codes are able to predict isothermal conditions when a high heat exchange is present.

$$n_{avg}=\frac{{{\displaystyle \sum }}_{i=1}^N{n}_i}{N}$$

(30)

The results obtained for the CO₂-based HPES system at high convective heat transfer coefficient are presented through a pressure-enthalpy plot (Figure 5). It is evident that the charging and discharging stages follow closely the isotherm, corresponding to that of the surrounding seawater temperature. This confirms the ability to model the isothermal behavior, even during the phase change process involving CO₂.

In the second test (simulation 2, Table 2), the heat transfer across the accumulator boundary is set to zero in order to check whether the codes could represent adiabatic behavior. Adiabatically compressing air within the accumulator leads to a thermal efficiency of 100% while storing 0.808 MWhr. Similarly, the thermal efficiency for the CO₂ is equal to 100% while storing 0.96 MWhr.

Compressing and expanding air adiabatically should yield a polytropic index of 1.4 [32]. Calculating the polytropic index for the adiabatic system through Equation (30) results in a polytropic index equal to 1.42. Hence, in switching the heat transfer off, the code yields to an adiabatic behavior. Compressing and expanding the CO₂ without heat transfer results in the gas to exhibit temperatures which deviate from that of the surrounding seawater. The pressure–enthalpy plot for the CO₂ process is shown in Figure 6. The cycle of the HPES deviates away from the liquid–vapor dome, that is, the CO₂ remains in the gaseous state. Moreover, Figure 6 indicates that CO₂ follows the same path for the charging and discharging stages, highlighting that no energy is lost or gained.

5.2. Influence of the Accumulator Pipeline Length-to-Diameter Ratio

The pipe diameter is varied from 36 to 84 inches, as discussed in Section 4 and Table 2 (simulation 4), while maintaining the original volume fixed at 1783 m³. Varying the pipeline diameter causes the length-to-diameter ratio to vary as indicated in

Table 4. Pipeline diameter and the different length-to-diameter ratios along with the round-trip thermal efficiency for air and CO₂.

Pipeline Diameter (inch)	Length-to-Diameter Ratio	Round-Trip Thermal Efficiency (%)
		Air	CO2
36	3286	95.1	86.7
48	1347	89.8	79.1
60	700	82.7	70.0
84	188	68.7	68.3

. From the prediction for the thermal efficiency, it is evident that decreasing the length-to-diameter ratios has a negative effect on the round-trip thermal efficiency for both gases. Moreover, the thermal efficiency of CO₂ approaches that of air at a pipeline diameter of 84 inches, corresponding to the lowest length-to-diameter ratio.

The length-to-diameter ratio significantly influences the thermal performance and heat transfer across the accumulator, since for the greater diameters, the final temperature of the gas of the charging stage is much greater (see

Table 5. Root mean square and final temperature of the charging stage along with the pressure drop experienced during the first hold stage for both air and CO₂.

Pipeline Diameter (inch)	Root Mean Square (K)		Final Temperature (K)		Pressure Drop during First Hold Stage (bar)
	Air	CO2	Air	CO2	Air	CO2
36	4.66	3.26	289.19	290.13	1.59	7.88
48	9.75	5.51	296.40	292.36	3.22	10.59
60	16.37	9.29	306.21	294.79	5.25	13.59
84	29.38	18.32	328.19	309.67	9.21	10.08

). Equation (31) is the root mean square (RMS) of the difference between the temperature of the compressed fluid and that of the surrounding seawater. An RMS value of zero would indicate purely isothermal conditions as no temperature difference between the compressed fluid and the seawater would be present. The RMS of the charging stage increases for smaller length-to-diameter ratios, as the temperature of air or CO₂ deviates away from that of the surrounding seawater during charging. Consequently, a greater pressure drop during the first hold stage is experienced, as indicated in Table 5. The greater pressure drop implies that a larger fraction of the energy utilized to compress the air or CO₂ is eventually lost to the surrounding seawater in the form of thermal energy. Similar behavior is noted during the discharging stage, with the pipeline having a smaller length-to-diameter ratio encountering larger temperature drops with both air and CO₂.

Due to the different diameters, the CO₂ code is not always terminated at the same limit, as presented in

Table 6. Limit reached and final pressure of the CO₂ after the charging stage along with the energy stored and recovered for both air and CO₂.

System Diameter (inch)	CO2		Energy Stored (MWhr)		Energy Recovered (MWhr)
	Limit Reached	Final Pressure (bar)	Air	CO2	Air	CO2
36	$x$ = 0.2	53.34	1.00	1.55	0.95	1.34
48	$x$ = 0.2	56.24	0.98	1.60	0.88	1.26
60	$x$ = 0.2	59.52	0.96	1.68	0.80	1.18
84	60 bar	60.00	0.91	1.13	0.62	0.77

. This also influences the final pressure. Table 6 indicates that irrespective of the diameter size of the pipeline, CO₂ stores and extracts more energy than air. This is a consequence of the phase change undergone through the storage cycle. In the case of air, the energy stored decreases as the diameter increases. In the case of CO₂, this increases up to the 60-inch pipeline but reduces considerably at the 84-inch diameter pipeline.

$${\tau }_{\mathrm{iso}}=\sqrt{\left(\frac{1}{N}\right){\displaystyle \sum }_{i=1}^N{(T_{1{(\mathrm{t})}_{\mathrm{i}}}-T_4)}_{}^2}$$

(31)

5.3. Influence of the Pipeline Inspection Gauge

When studying the influence of the PIG (simulation 5, Table 2), it became evident that the efficiency is minimally influenced by the friction of the PIG as it travels through the pipeline during the charging and discharging stages. As noted from

Table 7. Round-trip thermal efficiency, energy stored and recovered for air and CO₂ while considering the accumulator set-up with and without the pipeline inspection gauge.

	Round-Trip Thermal Efficiency (%)		Energy Stored (MWhr)		Energy Recovered (MWhr)
	Air	CO2	Air	CO2	Air	CO2
Without PIG	95.07	86.67	0.9957	1.5457	0.9466	1.3396
With PIG	94.35	85.99	0.9994	1.5509	0.9430	1.3337

, the presence of the PIG only results in a round-trip thermal efficiency drop of less than 1% for both air and CO₂. The presence of the PIG, however, requires a marginally longer charging duration over that when it is not present, due to a fraction of the pump power utilized to overcome the PIG friction. This in effect allows the system to dissipate more heat to the surrounding seawater as the accumulator is charged. On the other hand, during the discharging phase, a fraction of the stored energy is utilized to overcome the PIG friction, leading to less energy to be recovered as useful energy. Irrespective of the PIG being present or not, CO₂ stores and extracts more energy than air; however, the thermal efficiency is less, as shown in Table 7.

5.4. Influence of the Inner Liner

The introduction of the inner liner (simulation 6, Table 2) influences the thermal efficiency of the accumulator, as could be noted from

Table 8. Influence of the inner liner over thermal efficiency, and energy stored and recovered for air and CO₂.

	Round-Trip Thermal Efficiency (%)		Energy Stored (MWhr)		Energy Recovered (MWhr)
	Air	CO2	Air	CO2	Air	CO2
Without Inner Liner	95.07	86.67	0.9957	1.5457	0.9466	1.3396
With Inner Liner	93.18	79.10	0.9907	1.2633	0.9231	0.9992

. This is especially the case for the CO₂-based accumulator, where it experiences approximately a drop of 7% in thermal efficiency. The equivalent drop for the air-based system is less than 2%. The inner liner also influences the energy stored and recovered, as may be observed from Table 8.

Figure 7 presents the plots of the accumulator pressure versus the compressed gas volume for an air-based system during charging and discharging. The disparity due on the presence of the liner is negligible. An important minor observation from Figure 7a is that for the same volume, the pressure is greater for the system with the inner liner than without during the charging stage. On the other hand, during the discharging stage shown in Figure 7b at the same volume, a lower pressure is attained at a given volume when the inner liner is present. This further highlights the influence of the inner liner on the system behavior. The maximum temperature during the charging stage is estimated to reach 289.66 K when the inner liner is not present. The corresponding maximum temperature is found to reach 292.25 K with the inner liner. Moreover, the minimum temperatures attained following the complete discharge are found to be 277.26 K and 275.32 K without and with the liner, respectively.

Figure 8 shows the plots for the pressure versus specific enthalpy for the CO₂-based accumulator, with and without the inner liner. From Figure 8a, it is evident that the inner liner causes the CO₂ not to liquify and the code is terminated at 60 Bar rather than at a dryness fraction of 0.2. During charging, the maximum temperature reached in the presence of the inner liner is equal to 295.13 K, while that with no liner present is estimated to be 290.13 K. Figure 8b indicates that the starting point of the discharging stage is influenced by the inner liner; thus, the thermal energy gained during the charging stage in the presence of the inner liner is not completely dissipated. Moreover, the sudden pressure drop during discharging in the presence of the inner liner indicates that the model is predicting that CO₂ experiences a sudden temperature drop as it expands. Yet, it is still observed that at the end of the discharging stage, both accumulators return to a similar final temperature.

In the presence of the inner liner, the minimum temperature predicted to be reached by the model is 265.67 K, while that without the inner is equal to 277.28 K during discharging. This confirms that the liner does increase the thermal resistance substantially, with the CO₂-based system more influenced by the inner liner than that operating with air.

5.5. Combined Influence of the Inner Liner and Pipeline Inspection Gauge

Studying the combined influence of the inner liner and PIG (simulation 7, Table 2), it is evident that the thermal efficiency for the accumulators using air as the compressible fluid are more thermally efficient, with efficiency values greater than 90%, irrespective of the inner liner and PIG being present (see

Table 9. Reduction influence of the inner liner and pipeline inspection gauge over the thermal efficiency, energy stored and recovered for air and CO₂.

	Thermal Efficiency (%)		Energy Stored (MWhr)		Energy Recovered (MWhr)
	Air	CO2	Air	CO2	Air	CO2
Without Inner Liner and PIG	95.07	86.67	0.9957	1.5457	0.9466	1.3396
With Inner Liner and PIG	92.94	78.87	0.9920	1.2648	0.9219	0.9975

). On the other hand, a CO₂-based system experiences a significant round-trip thermal efficiency drop with the introduction of the inner liner and PIG. The energy stored for the air-based accumulator is relatively not influenced by the presence of the inner liner and PIG. On the contrary, a CO₂-based system experiences a reduction in the energy stored of around 0.28 MWhr. Similar behavior can be observed from Table 9 for the energy recovered for air and CO₂.

Throughout the charging stage, the accumulator utilizing CO₂ without the inner liner and PIG yields the least RMS (Figure 9). This indicates that the temperature of CO₂ deviates from the surrounding seawater by a smaller magnitude than that for the remaining accumulators for the whole stage. A similar trend is observed for the discharging stage. The greatest deviation of the gas from the surrounding seawater temperature occurs for the discharging stage when CO₂ is utilized with both the inner liner and PIG.

A trend is visible for the two gases as air yields a higher RMS for charging than discharging, as opposed to what is observed for CO₂. It is also evident for the charging stage that accumulators utilizing CO₂ yield lower RMS values compared to when the air is the compressible gas. Similar behavior is noted for the discharging stage except for the system utilizing CO₂ and including the inner liner and PIG.

Operating with CO₂ as the compressible gas results in the greatest amount of energy to be recovered per unit accumulator volume, irrespective of whether or not the inner liner and PIG are present. This may be noted from Figure 10a. Similar behavior is present for the energy recovered per unit mass of steel of the outer pipeline (see Figure 10b). However, when the inner liner and PIG are introduced, the increase in recoverable energy per unit volume of steel from an accumulator operating with CO₂ is only being predicted to be marginally larger than that for a system operating with air.

6. Conclusions

This study provided a preliminary understanding of the expected transient thermal behavior of an HPES system consisting of a subsea pipeline to act as the accumulator for storing energy in the form of a compressible fluid. The numerical modeling compared the thermal behavior of the subsea accumulator utilizing CO₂ instead of air, while introducing an inner polymeric liner and a PIG to mitigate issues related to corrosion and as dissolution of the gas being compressed into the seawater.

The simulations indicate that the accumulator’s length-to-diameter ratio has a significant influence on the round-trip thermal efficiency. This is observed with air and CO₂-based systems. It is also noted that the friction of the PIG has a negligible effect on the efficiency when using either compressible fluid. However, a notable drop in the round-trip thermal efficiency is predicted when the inner anti-corrosion polymeric liner is introduced. The drop is considerably larger for the HPES system operating with CO₂ than air. Finally, the transient numerical modeling indicates that the achievable storage capacity of the CO₂-based accumulator is around 0.55 MWhr higher than that of air-based systems (around 55.2% increase) when no inner liner and PIG are included. However, this difference reduces to 0.27 MWhr (around 27.5% increase) when the liner and PIG are included. The presented study provides motivation for exploring in further depth the potential of utilizing CO₂ in subsea accumulators for energy storage applications. It is essential for future work to focus on more comprehensive numerical modeling able to capture the transient effects associated with the phase change across the energy storage cycle. It is also imperative that the numerical models are validated through physical experiments of scaled accumulators in real marine conditions. The measurement data from such tests would be essential to validate the numerical modeling tools and asses the validity of the heat transfer correlations used for this study. The data will allow for a better understanding of the underlying physics affecting the thermal efficiency as well as the effective storage capacity that can be realistically achieved in upscaled systems.