A sensitivity study of hydrogen mixing with cushion gases for effective storage in porous media

Abstract

This study investigates the impact of methane, carbon dioxide, nitrogen, and their mixtures on hydrogen storage and recovery in depleted oil and gas reservoirs, with a focus on gas mixing and diffusion behaviour. Using the MATLAB Reservoir Simulation Toolbox (MRST), this research evaluates the effects of geological parameters such as permeability and porosity on hydrogen recovery and purity. Our results indicate a significant decline in hydrogen purity after 250 days of withdrawal, with methane providing the highest recovery rates but also exhibiting the greatest mixing with hydrogen. Permeability and porosity distributions were found to influence recovery and mixing dynamics in distinct and contrasting ways. While hydrogen purity and recovery were more sensitive to permeability, the computed mixing zone showed a marked sensitivity to the formation's porosity. By applying various control strategies for well injection rates and bottom hole pressures during injection and withdrawal phases, the optimal conditions for maximizing hydrogen recovery while minimizing gas mixing were identified. This study provides valuable insights into the optimization of hydrogen storage in subsurface geological formations, contributing to the advancement of sustainable energy solutions.

---

1 Introduction

The increasing global demand for sustainable energy sources is driven by the urgent need to address environmental challenges, particularly those related to carbon emissions and energy security. Hydrogen, as an energy carrier is considered one of the most promising alternatives, as its combustion produces only water, offering a clean and sustainable solution for various energy needs.¹ Despite these advantages, the widespread adoption of hydrogen energy has been hindered by challenges related to hydrogen storage, primarily due to hydrogen's low molecular weight.²

Hydrogen can be stored in different forms—solid, liquid, or gas—through either surface or subsurface facilities. Gaseous hydrogen storage can be classified into two main types: compressed hydrogen storage and underground hydrogen storage (UHS). In compressed hydrogen storage, hydrogen is kept in its gaseous form under high pressure within specially designed high-pressure cylinders or large storage vessels capable of withstanding the elevated pressures. In contrast, UHS involves storing hydrogen gas in subsurface geological formations. Liquid storage involves cooling and converting hydrogen gas into a liquid state at extremely low temperatures. The liquefied hydrogen is then stored in insulated cryogenic tanks that are designed to maintain these low temperatures and safely contain the hydrogen. This method provides a higher energy density compared to gaseous storage. Solid-state storage, on the other hand, involves storing hydrogen within solid materials such as metal hydrides, and chemical hydrides, or by adsorption onto porous materials like carbon. This approach allows for hydrogen to be stored at lower pressures and temperatures compared to gaseous or liquid forms, potentially increasing energy density while reducing infrastructure demands.^3–7

Among these methods, UHS is particularly advantageous (compared to compressed gaseous form, liquid and solid forms) due to its superior safety, large storage capacity (giga to tera watt scale), stable temperature and pressure conditions, compatibility with existing infrastructure, minimal environmental impact, long-term seasonal storage potential, versatility in accommodating multiple gases, cost-effectiveness, and the benefit of ongoing technological advancements.^8–11 UHS leverages natural underground formations such as salt caverns, aquifers, and depleted gas reservoirs to store hydrogen safely and efficiently. Each UHS option has its unique advantages and challenges.¹² Salt caverns, for instance, offer high security due to the impermeability and mechanical strength of salt rocks, allowing rapid hydrogen injection and withdrawal. However, they are geographically limited and require significant initial investment due to complex construction and maintenance requirements.^3,4 Aquifer storage is another option, using the voids in underground aquifers to store hydrogen. Aquifers are widely distributed across the Earth's surface, making this method more accessible and cost-effective in many regions. However, aquifer storage faces challenges related to hydrogen leakage, water contamination, and ecosystem disruption. There is also a lack of standardized procedures and regulatory support for this method, contributing to negative perceptions and legal uncertainties.¹³ Depleted gas reservoirs offer a more economical and reliable option for hydrogen storage, as they already possess existing infrastructure from previous natural gas extraction operations. These reservoirs provide large storage capacities and geological stability, making them ideal for large-scale hydrogen storage. However, hydrogen storage in depleted reservoirs requires careful management of injection and production processes to prevent excessive mixing with cushion gases (base gas), and the residual oil and gas fluids, which can affect the purity of the recovered hydrogen.¹⁴

The cushion gas serves as the foundational gas typically injected before hydrogen. Its key functions include (i) acting as a pressure buffer within the storage reservoir, (ii) enhancing safety by reducing the risk of hydrogen loss via residual trapping, and (iii) facilitating the controlled production and withdrawal of hydrogen from the reservoir.^15,16 The choice of cushion gas type plays a crucial role in depleted gas reservoirs as the remaining base gas makes it the most preferred option currently for UHS. Nitrogen (N₂), methane (CH₄), and carbon dioxide (CO₂) are commonly used, based on their distinct impacts on storage capacity, pressure maintenance, and hydrogen purity. For instance, nitrogen offers excellent pressure support due to its low molecular weight and high compressibility, while methane and carbon dioxide provide varying degrees of hydrogen recovery efficiency depending on their solubility and diffusivity in the reservoir.^17,18

As a result, recent investigations have focused on different criteria for the evaluation of cushion gas mixing for depleted gas reservoirs.^18–22 For example, a CMG GEM reservoir simulator was used by Terstappen²⁰ to analyse hydrogen mixing with CH₄ cushion gas in Dutch sandstone gas fields. The study revealed that molecular diffusion had a negligible impact on mixing compared to mechanical dispersion at typical reservoir flow rates. Kanaani et al.¹⁸ used CMG GEM to investigate the role of different cushion gases (N₂, CH₄, and CO₂) in depleted oil reservoirs. The study focused on the type and cushion gas composition. Additionally, the research examined the behaviour of hydrogen under various hydrodynamic reservoir phenomena, such as gravity segregation and overriding, across different scenarios. Their findings demonstrated that the injection of cushion gas mitigates the gravity segregation mechanism and reduces its negative impacts, with CH₄ being the most effective cushion gas, yielding a hydrogen recovery of 89.7%. Bo et al.¹⁹ also utilized the CMG GEM to explore the impact of geological heterogeneity on gas mixing in braided-fluvial reservoirs, using a mixture of 10% H₂ and 90% CH₄ as cushion gas. Due to the lateral continuity of these reservoirs, their study revealed that macro-scale geological heterogeneity significantly intensifies gas mixing, reducing the hydrogen fraction in the produced stream. Moreover, the study found that reservoir structure, particularly dip angle, and closure area, plays a more critical role than facies heterogeneity in determining CH₄ breakthrough and gas mixing behaviour across the tested conditions. In a more recent study, Ghaedi et al.²² employed Eclipse 300 to investigate the mixing dynamics of hydrogen with cushion gases (CH₄ and N₂) and its impact on the recovery factor (RF). Their study evaluated the effects of reservoir permeability, well perforation length, intelligent completion, and fractures on RF. They found that initiating hydrogen storage at higher hydrocarbon gas recovery factors led to improved RF. Additionally, hydrogen storage in lower-permeability reservoirs resulted in higher RF, provided that well pressure limits were not a concern. Shorter perforation lengths at the top of the formation yielded optimal RFs, while intelligent completions enhanced hydrogen recovery by enabling purer hydrogen production. However, the presence of natural fractures notably reduced hydrogen recovery, particularly during the initial storage cycles. Kobeissi et al.²¹ recently utilized nuclear magnetic resonance (NMR) (an experimental approach) to investigate the dispersion of hydrogen in various potential cushion gases, including CH₄, N₂, and CO₂. This study builds on their previous work,²³ where they measured the dispersion coefficient of hydrogen with N₂ as the sole cushion gas in Berea sandstone. The findings of their recent research provided valuable experimental data on dispersion coefficients, which serve as critical inputs for reservoir simulation models.

The majority of recent studies^{18–20,22,24}via numerical simulation (for evaluating the effect of mixing hydrogen with cushion gases in depleted gas reservoirs) have utilized commercial simulators such as CMG-GEM and Eclipse 300, which are well-established tools for multiphase flow and compositional simulations. However, for this study, the MATLAB Reservoir Simulation Toolbox (MRST) was selected due to its open-source nature, its flexibility and modularity via add-on modules (e.g., for compositional simulation and CO₂ storage), flow diagnostics and enhanced visualisation capabilities^25,26 (also see Section 2.1). Also, the use of MRST in this study offers the flexibility to explore novel approaches peculiar to the objectives of the research for a comprehensive investigation.

Therefore, this study's simulations examined the effects of geological factors, operational strategies, and numerical configurations on hydrogen storage efficiency, focusing on hydrogen purity, recovery, and the extent of mixing with cushion gases as the primary parameters for efficiency assessment. Previous simulation studies lack the quantification of the degree of mixing, without due consideration to the impact of injection and bottom hole pressure (BHP) profiles on the purity and recovery rates of hydrogen. The impact of permeability and porosity distribution on the effectiveness of hydrogen storage and withdrawal has also been largely overlooked. We account for these complexities, while also minimizing numerical dispersion to obtain a realistic indication of hydrogen's mixing tendencies during the storage and withdrawal phases.

2 Methodology

In this study, numerical simulations were conducted to evaluate the effects of geological, numerical, and operational parameters on the performance of underground hydrogen storage. The simulations followed a default configuration involving one year of hydrogen injection and one year of production for comparison. The injection was carried out at a base rate of 0.85 × 10³ mega m³ per year, with different injection rates over time: 25% of the base rate during the first 0.3 years, the base rate from 0.3 to 0.65 years, and 175% of the base rate after 0.65 years. During the production phase, BHP was used to control production, with a maximum allowable pressure (maxP) of 50 bar. The BHP was set to 90% of maxP during the first 0.3 years, 80% from 0.3 to 0.65 years, and 60% after 0.65 years. Our reasoning here was to account for the decline in reservoir pressure during the production phase, thus necessitating a reduction in the bottom hole pressure (BHP) to ensure continuous hydrogen production at a reasonably high rate – the reference case. However, BHP control in practice could be for reasons other than this; hence other BHP control scenarios are explored as presented in Table 1. Tables 1 and 2 provide details of the simulation parameters.

Table 1 Independent and dependent variables involved in this study

Category	Independent variable	Range	Dependent variable
Numerical factors	Grid size	5 × 5 × 5	Mixing zone
		10 × 10 × 10	Computing time
		15 × 15 × 15
		20 × 20 × 20
Geological factors	Average absolute permeability	10–50 mD	H2 purity
	Average porosity	0.1–0.3	H2 recovery
	Permeability distribution	Lognormal, uniform, bimodal	Mixing zone
	Porosity distribution	Normal, uniform, bimodal	Injection well BHP
Operational factors	Injection method	Step-accelerating (case A)	H2 purity
		Intermittent (case B)	H2 recovery
		Constant rate (case C)	Mixing zone
		Step-decelerating (case D)	Component flow rate
	Production method	Step pressurisation (case E)	Phase flow rate
		Intermittent (case F)
		Constant pressure (case G)
		Step depressurisation (case H)
	Cushion gas	100% CO2
		100% CH4
		100% N2
		33% each

---

Table 2 Default configuration of gas–oil model

Category	Name	Configuration
Grid	Resolution	15 × 15 × 15
	Reservoir size	1000 m × 1000 m × 1000 m
Rock	Horizontal permeability (kh)	30 mD
	Vertical permeability (kv)	3 mD
	Porosity (φ)	0.2
	Permeability distribution	Lognormal
	Porosity distribution	Normal
	Porosity error	0.02
Black oil fluid	Residual oil saturation (Sor)	0.2
	Residual gas saturation (Sgr)	0.05
	Oil Corey exponent for (no)	2.5
	Gas Corey exponent (ng)	2.5
Operating strategy	Injection	Case A
	Production	Case H
	Cushion gas	CH4
Initial mixture components	Gas–oil	CH4: 0.94
		n-Undecane: 0.03
		n-Dodecane: 0.03
Two-phase contact	Gas–oil contact	500 m
Basic reservoir conditions	Temperature (TR)	50 °C
	Pressure (PR)	50 bar

---

2.1 Simulation tool and platforms

This study utilized the MATLAB Reservoir Simulation Toolbox (MRST) due to its robust capabilities for simulating fluid flow in porous media, which is essential for accurately modelling hydrogen storage in depleted gas reservoirs.^27,28 MRST is an open-source platform that can be customized with specific features tailored to the requirements of hydrogen storage simulations. Originally designed for simulating oil and gas recovery and carbon capture, MRST was extensively modified in this study to accommodate hydrogen storage scenarios. Additionally, MRST offers powerful visualization tools that allow for detailed analysis of simulation outcomes, enhancing understanding and facilitating informed decision-making.^25,26

More specifically, the simulations were conducted using MRST's compositional module, which is based on Peng–Robinson's equation of state.¹⁷ The cushion gases tested in the simulations included methane, nitrogen, carbon dioxide, and various mixtures of these gases. Additionally, n-undecane and n-dodecane were employed to simulate the residual oil phase within the reservoir.

2.2 Modelling process

The compositional module in MRST employs the Peng–Robinson (P–R) equation of state, with the binary interaction coefficients sourced from the Unisim® process simulator. A simulation grid with a resolution of 15 × 15 × 15 cells was constructed (after a grid sensitivity study – Section 3.1), upon which different permeability and porosity distributions were applied to represent varying geological factors. Additionally, the grid resolution was adjusted to account for changes in numerical factors, allowing for a comprehensive analysis of the impact of these variables on the simulation results. The distribution of porosity differs somewhat from that of permeability. We specifically tested normal, uniform, and bimodal porosity distributions for porosity; whereas lognormal, uniform, and bimodal distributions were utilized for permeability (Fig. 1). The choice of bimodal and uniform permeability distributions was made to explore the sensitivity of hydrogen recovery and mixing behaviour to permeability variations across contrasting scenarios. The uniform distribution represents an idealized, homogeneous reservoir, while the bimodal distribution reflects a simplified model for reservoirs with two dominant permeability regimes, potentially mimicking high-permeability layers or zones interspersed with lower-permeability regions. It is worth mentioning that irrespective of the distribution type employed, the same average permeability and porosity values were maintained for comparability. It is also necessary to define the relative permeability function for the oil, gas and water phases as shown in Table 2.^29,30 The cubic model used herein provides a neutral framework, which can be adapted to more complex scenarios that incorporate specific geological heterogeneity.

Fig. 1 Permeability and porosity distributions of the numerical model. Property distributions were selected to enable a controlled sensitivity analysis rather than to replicate a particular sedimentary environment.

Since the compositional module in MRST assumes that two phases are miscible, implementing three-phase splitting would require a three-phase flash calculation, which was not pursued here. Therefore, the model was configured for oil and gas (as the miscible phases), with the water phase treated as immiscible. In setting the boundary conditions, the boundary fluid flow rate at the top of the reservoir was set to 0, as described by eqn (1) to simulate hydrogen trapping. Additionally, gravity effects were included in MRST, as well as vertical equilibration of the fluids (see included code). The reservoir depth in the z-direction was set to 1000 m, which is reasonable for depleted gas reservoirs.

Boundary condition = Flux_out = 0 (1)

It was important to quantify the degree of mixing between hydrogen and the cushion gas during the injection phase, necessitating the definition of a mixing zone. The hydrogen concentration was used as the primary metric to determine whether a cell block was included in the mixing zone, with the hydrogen concentration needing to satisfy eqn (2)

0.1 < c_H2 < 0.9 (2)

The hydrogen recovery and purity are calculated based on eqn (3) and (4), respectively during the production phase. Where Q_H2 is the hydrogen flow rate from the well, and Q_g,i is the flow rate of each gas component in the mixture. N represents the total number of gas components in the mixture.

2.3 Hydrogen injection/production strategy cases

During the hydrogen injection phase, four simulation scenarios are tested to compare the effects of different injection strategies: intermittent injection (case B), constant rate injection (case C), and step-deceleration injection (case D), using step-acceleration mode (case A) as the baseline (Fig. 2). The objective is to analyze how these injection modes influence the mixing zone ratio, (as this directly impacts the quality of subsequent hydrogen production), as well as the BHP, which affects both the safety of a project and hydrogen recovery efficiency. Similarly, during the hydrogen extraction phase, where production is controlled by BHP, the simulation evaluates four scenarios: intermittent production (case F), constant pressure (case G), and step pressurization (case E), with step depressurization mode (case H) serving as the baseline (Fig. 2). The injection and production strategies were designed to ensure that the integrals of the curves were equivalent across all cases. For the injection rate curves, this means that the area under all curves is the same, ensuring that the total amount of hydrogen injected remains constant. In the production phase, the integrals of the BHP curves were also kept consistent. This approach facilitates accurate comparisons across the different scenarios.

Fig. 2 Injection strategies tested include; case A: step-accelerating injection, case B: intermittent injection, case C: constant rate injection, case D: step-decelerating injection. Production strategies tested include; case E: step pressurization, case F: intermittent production, case G: constant pressure production, and case H: step depressurization.

2.4 Sensitivity analysis

In this study, a sensitivity analysis was conducted on key variables of interest, and the independent and dependent variables involved are shown in Table 1. The adopted ranges are based on existing simulation- and experimental studies on hydrogen storage.

Table 2 highlights the simulation setup used for the base case. Except for the varied parameters in the sensitivity analyses (Table 1), all other simulation parameters are maintained at the values shown in Table 2.

3 Results & discussion

3.1 Numerical factors

It was necessary to find a balance between simulation accuracy and computational cost. By performing the simulations with multiple grid resolutions, the sensitivity of the results to the total number of cells in the model could be ascertained. It is worth mentioning that numerical dispersion in coarser grids amplifies non-physical mixing due to discretization errors, making the simulated mixing zone appear larger than it would with finer grids (Fig. 3). Additionally, coarser grids fail to capture fine details of fluid fronts and interfaces, resulting in a blurred transition zone/smoothed concentration gradient and an apparent increase in the mixing zone. However, higher-resolution grids can better capture the subtle heterogeneities in the reservoir, thereby providing more accurate simulation results. Fig. 3a shows that, as the grid resolution increases, the mixing zone gradually decreases. While the difference in the mixing zone is more prominent between N = 5 and N = 10, this difference is minimised between N = 15 and N = 20. However, increasing N from 15 to 20 shows an explosive increase in computation time, because of the higher resolution mesh that needs to handle more data points and complex calculations and the simulation configuration itself (Fig. 3b). It was realised that at higher mesh resolutions, some computational time steps require significantly more iterations to converge. This was addressed by adjusting the injection rate and the time step size. Based on the comparison in Fig. 3, the 15 × 15 × 15 grid was determined as the best choice in terms of stability, accuracy, and computational time. Table 3 illustrates the computer specifications utilized in generating the computational results shown in Fig. 3.

Fig. 3 The effect of grid resolution on the (a) mixing zone; and (b) computational time in the hydrogen production phase.

Table 3 Software and hardware configuration

Software and hardware	Configuration
CPU	Intel i7-9750H
GPU	NVIDIA GeForce 1660Ti
Memory	16 GB
MATLAB version	2019a
MRST version	2023b

---

3.2 Default configuration

3.2.1 Hydrogen spatial distribution in the injection phase. The hydrogen concentration contour plots are presented in Fig. 4. During the early stage (day 1 to day 100), hydrogen is primarily concentrated near the injection point, with concentrations rising from 0.35 to 0.8 in this area. In contrast, hydrogen concentrations in regions further from the injection port remain below 0.1, forming a clear concentration gradient. As time progresses, hydrogen spreads deeper into the reservoir, and the concentration distribution becomes more uniform. By the middle stage (day 100 to day 200), the high-concentration region expands significantly, covering a larger portion of the reservoir. The transition zone of hydrogen concentration also extends throughout the reservoir.

Fig. 4 Reservoir hydrogen concentration distribution on days 1, 50, 100, 200, 280, and 365 during the hydrogen injection phase.

By day 200, hydrogen concentrations near the injection point reached approximately 0.9, while concentrations in areas further away ranged from 0.2 to 0.3. In the later stage (day 200 to day 365), hydrogen concentration stabilizes, and the high-concentration region encompasses most of the reservoir. Hydrogen is distributed more evenly, and the concentration gradient decreases. By the end of the simulation, hydrogen concentrations near the injection port remained around 0.9, while concentrations in distant areas ranged from 0.4 to 0.5. The final mixing zone percentage, calculated from eqn (2), was determined to be 12.1%. The pore space and hydrodynamic properties allow hydrogen to accumulate rapidly near the injection site and eventually, over injection and diffusion periods, allowing the hydrogen phase to reach a dynamic equilibrium state within the reservoir.

3.2.2 Hydrogen spatial distribution in the production phase. The hydrogen production phase is initialized using the state variables of the last injection time step. During the initial stages of extraction, the pressure decreases rapidly in the vicinity of the production well (Prod), creating a low-pressure zone with pressure values around 48 bar (Fig. 5). The pressure in the zone away from the production well is relatively higher and uniform and is maintained at 50 bar. The low-pressure zone then gradually expands, and the pressure continues to drop in the region near the production well. In the middle stage, the low-pressure zone near the production wells further expands, and the pressure continues to decrease to 44 bar. The pressure in the high-pressure zone also begins to decrease gradually, indicating the expanding effect of hydrogen extraction on the reservoir pressure. Finally, the low-pressure zone covers almost the entire reservoir, and the pressure near the production wells reaches a minimum value of 36 bar (see Fig. 5).

Fig. 5 Changes in pressure distribution during the production stage. The colour bar represents the pressure in Pa.

In Fig. 6, it can be observed that the area of high hydrogen concentration diminishes during the production process. As hydrogen is extracted from the lower section of the reservoir, the high concentration of hydrogen near the production well is maintained. However, by the end of the production phase, a considerable amount of hydrogen remains unextracted within the reservoir. During the early and middle stages of hydrogen production, hydrogen purity remains high. After approximately 200 days, a more significant decrease in the concentration can be observed. This decline occurs because the initial flow composition through the production well consists almost entirely of pure hydrogen. As the area of high hydrogen concentration continues to shrink (as depicted in Fig. 6), the lower concentration regions gradually encroach upon the production well. These regions are part of the original mixing zone, and lead to a decrease in hydrogen purity. Afterwards, more cushion gas is extracted, leading to an increase in cushion gas concentration in the produced gas (as shown in Fig. 7a).

Fig. 6 Variation of hydrogen concentration during the production phase. The colour bar represents the hydrogen concentration.

Fig. 7 (a) Volume fraction of each component; (b) flow rate of each component; (c) flow rate of each phase after the hydrogen production phase.

3.2.3 Species concentration, fluxes, and phase flow rates. The default production strategy is case H (as shown in Fig. 2). As the BHP reduces, hydrogen is produced by a pressure gradient, (each BHP reduction step yields a higher-pressure gradient between the production well and the reservoir), and more hydrogen is produced (as shown in Fig. 7b); i.e., two decreases in BHP result in two hydrogen flow rate surges (Fig. 7b). After some time, the well approaches a steady state.

Fig. 7c shows that the total amount of gas in the early and middle stages of the production process is consistent with the hydrogen flux as shown in Fig. 7b. This is because the hydrogen produced in the early and middle stages is of high purity, while a considerable amount of methane is produced in the later stages; thus, the hydrogen flow rate decreases, and the total gas flow rate increases (Fig. 7c). A small amount of oil is also produced during the entire production process, maintaining the same trend as the hydrogen flow profile. The above (Fig. 4–7) are the test results using the default configuration (Table 2). Next, the reservoir geological conditions and operating parameters are varied to observe their impacts on the storage and withdrawal performance.

3.3 Sensitivity analyses

3.3.1 Geological factors.
3.3.1.1 Effect of permeability on hydrogen recovery, hydrogen purity, and well BHP. It can be observed that regardless of the permeability distribution, increasing the permeability can improve the hydrogen recovery, but the accompanying negative consequence is the considerable reduction in hydrogen purity in the later stages of production (Fig. 8). At lower permeabilities, the lognormal and uniform distributions show similar performance (Fig. 9b); however, at higher permeabilities, the lognormal distribution maintains a higher purity because of its long tail (as captured in Fig. 1), i.e., some reservoir regions with very high permeability provide more efficient flow channels, aiding in high-purity hydrogen extraction. At the same time, the many regions of relatively lower permeability in the lognormal distribution help maintain the purity. The purity and cumulative recovery values reported in Fig. 8 are those at each time point; whereas, the time average values calculated over the entire production period are subsequently reported in Fig. 9.

Fig. 8 Hydrogen recovery and purity under different mean permeabilities using different distributions (a) lognormal distribution, (b) uniform distribution, and (c) bimodal distribution in the hydrogen production phase.

Fig. 9 (a) Final hydrogen recovery and (b) mean hydrogen purity at different mean permeability values employing 3 key distributions in the hydrogen production phase.

Among the tested distributions, hydrogen purity was consistently highest under the lognormal distribution, followed by the uniform distribution, with the bimodal distribution resulting in the lowest purity. At lower permeability levels (10 mD and 20 mD), hydrogen purity was relatively higher than at higher permeability values. At low permeability levels, the flow of both gases (hydrogen and cushion gas) is significantly restricted, reducing the fluid mobility and the extent of molecular diffusion and convective mixing between the two gases.

In our study, hydrogen recovery increases with increasing permeability, which is in contrast with the findings of Ghaedi et al.,²² where very high permeability leads to easier lateral diffusion of hydrogen, making recovery more difficult. This discrepancy may be attributable to the permeability range in their study (50 mD and 750 mD), whereas the present study focuses on a much narrower and lower range of 10 mD to 50 mD. This may suggest that within the lower permeability range, increasing permeability enhances recovery by facilitating gas flow. However, at higher permeability levels, hydrogen diffuses more easily in the lateral direction, ultimately leading to a decrease in recovery.

The BHP is another crucial parameter in underground gas storage, affecting many aspects such as pressure management, reservoir protection, fluid dynamics, economic benefits and safety. By effectively controlling and managing the BHP, the production process can be optimized to ensure safe and efficient storage. Under all distribution types, a lower permeability results in a significant increase in the BHP (Fig. 10). Low permeability results in restricted gas flow and rapid pressure build-up. The higher the permeability, the better the gas mobility and the slower the pressure increase. Importantly, pressure fluctuations may be more pronounced under bimodal distributions due to the complication of gas flow paths resulting from the mixing/fluid interactions in the high and low permeability regions. However, as the gas injection time increases, the BHP under all distribution types and permeabilities shows a gradually increasing trend and finally reaches dynamic equilibrium. The pressure increase under the lognormal distribution is higher than that of the unform distribution, and this may be again due to the high permeability regions associated with the lognormal distribution.

Fig. 10 Injection well BHP under different permeability distributions, (a) lognormal distribution, (b) uniform distribution, (c) bimodal distribution.

3.3.1.2 Effect of porosity on hydrogen recovery, hydrogen purity, and well BHP. Next, the change in porosity is studied. The hydrogen recovery rate under different porosity conditions changes minimally, and the curves almost overlap (Fig. 11). The hydrogen recovery rate gradually increases over time, eventually reaching nearly 36%. Similarly, the hydrogen purity under different porosity conditions does not change significantly, but this is limited to the first 200 days. In the later period (>200 days), the hydrogen purity is considerably affected by porosity changes. Lower porosity results in a reduced storage capacity per unit volume of the reservoir, which can extend the duration of high-purity hydrogen. This occurs because less mixing with cushion gas maintains higher purity for a longer period as shown in Fig. 11. Although increased porosity enhances the storage capacity of the reservoir, its influence on fluid flow resistance remains insignificant. Thus, the indicators of hydrogen storage efficiency are more sensitive to permeability than porosity – since the permeability remains constant under the different porosity conditions illustrated in Fig. 11 and 12. While Fig. 11 captures hydrogen recovery and purity variation with time, Fig. 12 represents the time-averaged values. Thus, the fluid diffusion and convection behaviour, remain similar across different porosity conditions.

Fig. 11 Hydrogen recovery and purity under different mean porosities using different distributions (a) lognormal distribution, (b) uniform distribution, and (c) bimodal distribution in the hydrogen production phase.

Fig. 12 (a) Hydrogen recovery and (b) purity under normal, uniform, and bimodal distributions in the production phase.

Next, we analyse the effect of porosity on the BHP (Fig. 13). Again, all porosity distribution patterns follow the same trend, with an increase in porosity only leading to a slight decrease in BHP. The stepwise increase in BHP corresponds to different phases of gas injection (i.e., the different injection rates employed at different periods – see Tables 1 and 2). BHP is primarily influenced by the gas injection rate and volume, with changes in porosity having minimal effect at a constant permeability. Additionally, reservoir heterogeneity due to porosity variations has a limited effect on the hydrogen flow paths; an increase in porosity does not significantly alter these flow paths. Although a higher porosity increases the reservoir's total pore volume, the fluid compressibility is more critical as determined by gas expansion within the pores.

Fig. 13 Effect of porosity on injection well BHP, (a) normal distribution, (b) uniform distribution, (c) bimodal distribution.

3.3.1.3 Effect on permeability and porosity on the mixing zone. Again, a lower permeability implies that the gases (hydrogen and cushion gas) will have limited mobility in the storage medium. This slows down the mixing of hydrogen and cushion gas because gas migration between different regions of the reservoir is restricted, as shown in Fig. 14a. At mean permeabilities of 10 mD to 30 mD, the size of the mixing zone for the bimodal distribution is higher than that of other distributions. However, at higher permeabilities (≥ 40 mD), the size of the mixing zone for the bimodal distribution begins to decrease compared to other distributions. This is because the hydrogen dispersion reaches the boundary of the high-permeability area, and part of the hydrogen enters the lower-permeability zone. During extraction, this portion of hydrogen responds more slowly because it is located in the lower part of the reservoir, while most activity occurs in the high-permeability zone. This explains the reversal in the size of the mixing zone displayed by the bimodally distributed reservoir.

Fig. 14 Effect of mean (a) permeability and (b) porosity on the mixing region using the 3 distributions in the hydrogen injection phase.

By examining Fig. 14b, it can be observed that as porosity increases, the percentage of the mixing zone gradually decreases. This suggests that higher porosity may reduce the complexity of gas mixing, leading to a decrease in the proportion of the mixing zone. Interestingly, the influence of porosity on the mixing zone appears to be more significant than that of permeability for the tested ranges. With reduced porosity, the available void space for advective transport (bulk flow) is smaller. This could slow down the overall movement of gases, allowing diffusion to play a more dominant role in mixing. Additionally, in lower-porosity media, there is a more tortuous fluid motion, and dead-end pores can trap gases temporarily, further promoting mixing when these trapped gases slowly diffuse into adjacent flow paths. This increases the contact time between hydrogen and cushion gases, leading to more extensive mixing. Higher porosity facilitates advective flow (transporting the gas in bulk while limiting dispersion and mixing), which in turn is a function of the respective gas mobilities. Thus, in high-porosity scenarios, fluids may flow more rapidly along the main flow paths rather than diffusing broadly across the reservoir, leading to a reduction in the proportion of the mixing zone. Across the porosity conditions (Fig. 13b), the normal distribution mostly exhibits the highest percentage of the mixing zone, followed by the uniform distribution, while the bimodal distribution has the lowest percentage.

3.3.2 Operational parameters.
3.3.2.1 Effect of injection schedule on hydrogen recovery, purity, mixing zone, and BHP. By properly adjusting the operating parameters, the hydrogen storage process can be optimized, storage efficiency and safety improved, and costs and environmental impacts, reduced. This requires comprehensive consideration of reservoir characteristics, and equipment performance, to develop an optimized operating strategy. We first study the injection method.

It can be observed in Fig. 15a that while case D has a higher recovery rate, case B has the worst recovery ∼30%; however, the purity of hydrogen in case B supersedes all other methods. This is because the periodic injection and cessation of case B causes the flow of hydrogen in the reservoir to be discontinuous and slowed down. This in turn may result in a change in the flow paths of gas in the reservoir. In contrast, case D has the best recovery performance, but the worst purity (Fig. 15b). This is because, a gradual reduction in the injection causes a continuous redistribution of the reservoir fluids, potentially increasing the mixing area of different gas components, and thus reducing the hydrogen purity. The gradual reduction of the injection rate in case D provides a reasonable pressure management system so that hydrogen can be recovered more effectively.

Fig. 15 Effect of injection method on the (a) hydrogen recovery and (b) purity; (c) mixing zone; and (d) BHP during the injection and production phases.

As can be observed in Fig. 15c, the mixing zone of case B is the largest. This is because periodic injection and stoppage causes frequent fluctuations in the pressure in the reservoir, thus, promoting gas mixing and expanding the mixing zone. In addition, during each stop of injection, the flow rate of the gas in the reservoir slows down, increasing the role of diffusion, and enabling the mixing between different gas species. When the injection rate is restarted, the new hydrogen gas entering the reservoir potentially disturbs the existing gas distribution, causing further mixing of the gas and expanding the mixing zone. Case D has the highest recovery rate, but due to the enhanced diffusion effect during the gradual reduction of the injection rate, the mixing zone is larger, and the hydrogen purity is lower. Case C has a constant injection rate, which reduces pressure fluctuations and gas mixing and maintains a higher hydrogen purity and a smaller mixing zone. Case A performs moderately; the gradual increase in the injection rate helps to stabilize the flow path, and the mixing effect is moderate, resulting in an intermediate recovery rate and purity.

The injection method changes the BHP signature as shown in Fig. 15d. The BHP for case A gradually increases with the injection rate. Case B's periodic injection and stoppage led to large fluctuations in the BHP. The pressure rises rapidly at each injection and drops rapidly when the injection stops. This method leads to higher reservoir pressure fluctuations, which may also affect the flow paths and mixing tendency of hydrogen. Case C – constant injection rate leads to a steady increase in BHP in the early stage and this is maintained at a relatively stable level. This method can provide a stable pressure environment, which helps maintain the purity of hydrogen and a stable flow path. However, case C lacks practicality, as it is difficult to maintain a stable injection rate for a prolonged period. The gradual reduction of the injection rate in case D results in a rapid increase in BHP at the beginning and then a gradual decrease. The high injection rate at the early stages provides a large pressure driving force, but this decreases at the later stages, affecting the flow and recovery of hydrogen as earlier demonstrated in Fig. 15a.

3.3.2.2 Effect of production schedule on hydrogen recovery, purity, mixing zone, and BHP. Next, the effect of the production method on the hydrogen production properties is analysed. It can be observed from Fig. 16a that the hydrogen purity of case E (step pressurization) gradually decreases over time, from about 100% to about 92%. The hydrogen purity of case F (intermittent) decreases the fastest, from about 100% to about 84% (particularly in the latter stages). The hydrogen purity of case G decreases rather slowly, from about 100% to about 94%, throughout the production duration; whereas, the hydrogen purity of case H (step dpressurisation) decreases more rapidly, from about 100% to about 90%. Constant production pressure (case G) performs best in terms of hydrogen purity, as this provides a stable flow path, reduces gas mixing, and maintaining high hydrogen purity and a relatively stable recovery. The hydrogen and methane flow rate curves (Fig. 16c and d), generally follow the same observations with the rate-controlled injection wells in Fig. 15d. Although Case F performs best in terms of hydrogen recovery, the latter stages of production (>250 days) is associated with a significant recovery of the cushion gas.

Fig. 16 Effect of hydrogen production BHPs on (a) hydrogen purity, (b): hydrogen recovery, (c): hydrogen flow rate, and (d) cushion gas (methane) flow rate during the production phase.

3.3.2.3 Effect of cushion gas type on hydrogen recovery, purity, and flowrates. As can be observed in Fig. 17, nitrogen performs best in terms of hydrogen purity and is suitable as a cushion gas for hydrogen storage. Carbon dioxide is not as favourable for hydrogen recovery compared to methane (Fig. 17b). The 1 : 1 : 1 mixture performs moderately well, combining the positive attributes of the other cushion gases and may be practical in hydrogen storage projects, depending on gas availability. However, it should be noted that the performance of the cushion gas mixture (Fig. 17b) is closer to that of methane, indicating gravity separation effects in the mixture over time (given the density differences). Additionally, the results of the present study (Fig. 17b) are in good agreement with that of Kanaani et al.¹⁸ In their study, different production cycles were implemented, and the highest hydrogen recovery was obtained using methane, followed by nitrogen, and then carbon dioxide (with the lowest recovery). However, our study additionally places significant emphasis on hydrogen purity, as certain industrial applications require hydrogen to meet strict purity standards, and the separation and purification processes cannot be overlooked. Among the various methods for separating hydrogen, membrane separation technology is an effective approach relying on the selective permeability of membranes and the diffusion coefficients of different gases.^31,32

Fig. 17 Impact of cushion gas type on (a) hydrogen purity; (b) hydrogen recovery; (c) gas flow rate; and (d) oil flow rate during the production phase.

In general, hydrogen is more selectively separated from carbon dioxide and nitrogen but less selectively from methane. As a result, separating hydrogen from a methane-containing mixture may incur higher costs, as it necessitates a larger membrane area or more advanced materials to achieve the desired purity. While carbon dioxide can be easily separated from hydrogen, at higher pressures, CO₂ may dissolve into the membrane material, affecting both membrane performance and longevity, leading to increased long-term maintenance costs. Nitrogen, having larger molecules than hydrogen, is relatively easy to separate, often without requiring highly selective membrane materials, which keeps the separation costs lower. As shown in Fig. 17a, using methane as a cushion gas results in the lowest hydrogen purity, and this is a crucial consideration when using methane as a cushion gas. In contrast, nitrogen produces inherently purer hydrogen potentially lowering any separation costs. However, the source of nitrogen for this purpose remains a source of concern. It is also interesting to observe the difference in associated oil production as a result of the different cushion gases. CO₂ leads to a more significant residual oil production, and this may be related to its higher miscibility (relative to other cushion gases) with the residual oil phase at the reservoir's temperature and pressure.

3.3.2.4 Effect of cushion gas type on hydrogen mixing zone and reservoir pressure. As shown in Fig. 18a, the computed mixing zone for the four types of cushion gases is considerably similar, all hovering around 12%. However, methane shows the highest mixing percentage, while nitrogen exhibits the lowest. This is attributed to the fact that gases with lower molecular weights possess higher diffusion coefficients, allowing them to diffuse and mix with other gases more rapidly upon contact. Diffusion, thus serves as a key driver of gas mixing. This finding aligns with the work of Kobeissi et al.,²¹ where the diffusion coefficients were reported as ∼0.48 cm² min⁻¹ for nitrogen, ∼0.53 cm² min⁻¹ for methane, and ∼0.44 cm² min⁻¹ for CO₂, measured at a fixed interstitial velocity of 6 (u_m/cm min⁻¹)^1.2 and a fixed reservoir pressure of 50 bar. Similarly, Kanaani et al.,¹⁸ utilising a more complex geological model, noted that reducing the molecular weight of the cushion gas mixture improves hydrogen recovery, although this also promotes greater mixing and results in reduced purity as indicated herein. This similarity in our observations suggests that the fundamental attributes of the hydrogen–cushion gas interactions can be similar, irrespective of the geological complexity incorporated into the model. Nonetheless, it is worth mentioning that these observations are to be analysed in conjunction with the geological conditions of the reservoir (e.g., porosity, permeability) and the injection parameters (e.g., injection rate and pressure) for optimal decision-making in hydrogen storage projects. It is thus important that numerical simulations should be designed to be robustly sensitive to physical phenomena, particularly when grid resolution and computational resources are constrained – thus, allowing accurate distinction between the effects of different cushion gases.

Fig. 18 Effect of cushion gas on (a) the mixing zone during the hydrogen injection phase; and (b) the average reservoir pressure during hydrogen production (1 bar ≡ 100 000 Pa).

Fig. 18b shows the average reservoir pressure over time when different cushion gases are used during hydrogen production. During the entire process, (when nitrogen is used), the reservoir pressure drops the slowest, showing a strong pressure support effect and this could lengthen the production period at a high purity level. The performance of the mixed gas is between nitrogen and other gases; the pressure is relatively well maintained compared to CO₂ and CH₄. The reservoir pressure (when methane is the cushion gas) drops the fastest. Methane generally has a lower compressibility, so its pressure support effect is weaker particularly where long-term production is required. Kanaani et al.¹⁸ strongly support this observation (i.e., nitrogen demonstrates the best ability to maintain pressure, while carbon dioxide performs worst). Factors such as nitrogen's low molecular weight, good compressibility, and chemical inertness contribute to this outcome. Although methane is lighter, it is less effective than nitrogen because it diffuses and leaks more easily. Carbon dioxide, on the other hand, struggles to maintain sufficient pressure over time due to its higher molecular weight and greater solubility. Mixed gases, by combining the benefits of different pure gases, offer a more balanced approach and may be suitable for smaller-scale storage applications.

4 Study limitations

While we primarily aimed to isolate and analyse the hydrogen–cushion gas mixing and diffusion dynamics under the influence of key parameters (grid resolution, permeability, porosity, bottomhole pressures and injection flow rates), using a simple cubic model, this may not be fully representative/capture additional complexities of a full-scale geological model. Furthermore, the inclusion of geomechanical effects (which are important for designing injection operations), robust relative permeability and capillary pressure data (which are currently scarce in published literature) could be pursued to extend the model's reliability. Using Monte Carlo simulations for uncertainty quantification in our model (through multiple realizations) would constitute a more robust and generalizable approach, enabling more confident inferences from our preliminary results.

5 Conclusions

This study investigates the effects of different cushion gases on the mixing and diffusion characteristics of hydrogen in porous media, particularly in depleted oil and gas reservoirs. Using MRST, a comprehensive analysis of geological conditions, operating parameters, and numerical factors was performed to gain insights into hydrogen's behaviour in underground storage systems. While the computed mixing zone was primarily employed to evaluate the injection phase, the hydrogen purity and recovery were used principally to assess the production phase. The findings provide key insights for advancing underground hydrogen storage.

• Influence of cushion gas: the choice of cushion gas significantly impacts hydrogen recovery and purity. Nitrogen achieved the best balance between recovery and purity, outperforming methane and carbon dioxide. Methane, with its low molecular weight and lower solubility in the residual oil phase, provided enhanced hydrogen recovery. However, the extracted hydrogen had lower purity, potentially incurring higher separation costs. Mixtures of methane and carbon dioxide also demonstrated promising results, improving recovery while maintaining high hydrogen purity.

• Influence of geological and reservoir properties: the distribution of permeability and porosity played a crucial role in hydrogen mixing and recovery. A lognormal permeability distribution led to higher hydrogen purity. While a bimodal distribution improved recovery in certain cases, the presence of regions with very high permeability led to increased mixing effects, reducing hydrogen purity. In contrast, porosity had a lesser impact on recovery and purity but considerably influenced the size of the mixing zone.

• Effect of operating strategy: injection and production strategies are critical in optimizing storage and recovery performance. Constant-rate injection and constant-pressure production maintained steady hydrogen flow and high purity. Staged injection and production methods improved recovery but also widened the mixing zone, reducing hydrogen purity. Stepwise accelerated injection and depressurized production can be considered when constant injection and production are not feasible.

This study underscores the importance of selecting the appropriate cushion gas, understanding geological properties, and optimizing operational strategies to enhance the performance of underground hydrogen storage systems. These findings offer valuable guidance for future hydrogen storage projects for a potential transition to a hydrogen-based economy.

Data availability

The MATLAB code developed for this study can be found at: https://github.com/EEUoE/H2Storage MATLAB R2019a and MRST 2023b were used.

Author contributions

Junhan Lu: methodology, software, data curation, writing – original draft. Nasiru Salahu Muhammed: methodology, writing – original draft, writing – review draft. Jude A. Okolie: methodology, writing – original draft, writing – review draft, data curation. Emmanuel Epelle: conceptualization, methodology, software, writing – original draft, writing – review draft, project administration.

Conflicts of interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgements

The corresponding author gratefully acknowledges the financial support of the University of Edinburgh through a Chancellor's Fellowship.

References

1. P. Edwards, V. Kuznetsov and W. I. David, Hydrogen energy, Philos. Trans. R. Soc., A, 2007, 365(1853), 1043–1056.
2. A. Jaiswal and S. Sahu, Hydrogen Storage Challenge in the Hydrogen-Based Civilization, Hydrogen Fuel Cell Technology for Mobile Applications, 2023, pp. 157–181.
3. L. Zhang, C. Jia, F. Bai, W. Wang, S. An and K. Zhao, et al., A comprehensive review of the promising clean energy carrier: hydrogen production, transportation, storage, and utilization (HPTSU) technologies, Fuel, 2024, 355, 129455.
4. N. S. Muhammed, A. O. Gbadamosi, E. I. Epelle, A. A. Abdulrasheed, B. Haq and S. Patil, et al., Hydrogen production, transportation, utilization, and storage: recent advances towards sustainable energy, J. Energy Storage, 2023, 73(PD), 109207.
5. F. Qureshi, M. Yusuf, M. Arham Khan, H. Ibrahim, B. C. Ekeoma and H. Kamyab, et al., A State-of-The-Art Review on the Latest trends in Hydrogen Production, Storage, and Transportation Techniques, Fuel, 2023, 340, 127574.
6. S. Ahmad, A. Ullah, A. Samreen, M. Qasim, K. Nawaz and W. Ahmad, et al., Hydrogen production, storage, transportation and utilization for energy sector: a current status review, J. Energy Storage, 2024, 101(PA), 113733.
7. E. I. Epelle, K. S. Desongu, W. Obande, A. A. Adeleke, P. P. Ikubanni and J. A. Okolie, et al., A comprehensive review of hydrogen production and storage: a focus on the role of nanomaterials, Int. J. Hydrogen Energy, 2022, 47(47), 20398–20431.
8. M. Aslannezhad, M. Ali, A. Kalantariasl, M. Sayyafzadeh, Z. You and S. Iglauer, et al., A review of hydrogen/rock/brine interaction: implications for hydrogen geo-storage, Prog. Energy Combust. Sci., 2023, 95, 101066.
9. E. I. Epelle, W. Obande, G. A. Udourioh, I. C. Afolabi, K. S. Desongu and U. Orivri, et al., Perspectives and prospects of underground hydrogen storage and natural hydrogen, Sustainable Energy Fuels, 2022, 6(14), 3324–3343.
10. C. Sambo, A. Dudun, S. A. Samuel, P. Esenenjor, N. S. Muhammed and B. Haq, A review on worldwide underground hydrogen storage operating and potential fields, Int. J. Hydrogen Energy, 2022, 47(54), 22840–22880 , available from: https://www.sciencedirect.com/science/article/pii/S0360319922021899.
11. N. S. Muhammed, B. Haq, D. Al Shehri, A. Al-Ahmed, M. M. Rahman and E. Zaman, A review on underground hydrogen storage: insight into geological sites, influencing factors and future outlook, Energy Rep., 2022, 8, 461–499.
12. D. Zivar, S. Kumar and J. Foroozesh, Underground hydrogen storage: a comprehensive review, Int. J. Hydrogen Energy, 2021, 46(45), 23436–23462.
13. S. M. Jafari Raad, Y. Leonenko and H. Hassanzadeh, Hydrogen storage in saline aquifers: opportunities and challenges, Renewable Sustainable Energy Rev., 2022, 168, 112846.
14. M. S. A. Perera, A review of underground hydrogen storage in depleted gas reservoirs: insights into various rock-fluid interaction mechanisms and their impact on the process integrity, Fuel, 2023, 334(P1), 126677.
15. I. Izadi Amiri, D. Zivar, S. Ayatollahi and H. Mahani, The effect of gas solubility on the selection of cushion gas for underground hydrogen storage in aquifers, J. Energy Storage, 2024, 80, 110264.
16. M. Saeed and P. Jadhawar, Optimizing underground hydrogen storage in aquifers: the impact of cushion gas type, Int. J. Hydrogen Energy, 2023, 52, 1537–1549.
17. G. Wang, G. Pickup, K. Sorbie and E. Mackay, Numerical modelling of H₂ storage with cushion gas of CO₂ in subsurface porous media: filter effects of CO₂ solubility, Int. J. Hydrogen Energy, 2022, 47(67), 28956–28968.
18. M. Kanaani, B. Sedaee and M. Asadian-Pakfar, Role of Cushion Gas on Underground Hydrogen Storage in Depleted Oil Reservoirs, J. Energy Storage, 2022, 45, 103783.
19. Z. Bo, S. Hörning, J. R. Underschultz, A. Garnett and S. Hurter, Effects of geological heterogeneity on gas mixing during underground hydrogen storage (UHS) in braided-fluvial reservoirs, Fuel, 2024, 357(PC), 129949.
20. R. J. Terstappen, Analysis of mixing during hydrogen storage in gas reservoirs: a reservoir simulation study, Master of Science, Delft University of Technology, 2021.
21. S. Kobeissi, N. N. A. Ling, K. Yang, E. F. May and M. L. Johns, Dispersion of hydrogen in different potential cushion gases, Int. J. Hydrogen Energy, 2024, 60, 940–948.
22. M. Ghaedi, P. Østebø Andersen and R. Gholami, Mixing dynamics and recovery factor during hydrogen storage in depleted gas reservoirs, Gas Sci. Eng., 2024, 357(June), 129949.
23. K. Yang, S. Kobeissi, N. Ling, M. Li, L. Esteban and E. F. May, et al., Measurement of hydrogen dispersion in rock cores using benchtop NMR, Int. J. Hydrogen Energy, 2023, 48(45), 17251–17260.
24. E. R. Okoroafor, S. D. Saltzer and A. R. Kovscek, Toward underground hydrogen storage in porous media: reservoir engineering insights, Int. J. Hydrogen Energy, 2022, 47(79), 33781–33802.
25. S. Krogstad, K. Lie, O. Møyner, H. M. Nilsen, X. Raynaud and B. Skaflestad, MRST-AD – An Open-Source Framework for Rapid Prototyping and Evaluation of Reservoir Simulation Problems, in SPE Reservoir Simulation Symposium, Society of Petroleum Engineers, 2015.
26. K.-A. Lie and O. Moyner, Advanced Modeling with the MATLAB Reservoir Simulation Toolbox, ed. K.-A. Lie and O. Møyner, Cambridge University Press, 2021.
27. J. O. Helland, H. A. Friis, M. Assadi and S. Nagy, Machine learning for underground gas storage with cushion CO₂ using data from reservoir simulation, IOP Conf. Ser.: Mater. Sci. Eng., 2023, 1294(1), 012058.
28. E. Sarı and E. Çiftçi, Underground hydrogen storage in a depleted gas field for seasonal storage: a numerical case study of the Tekirdağ gas field, Fuel, 2024, 358(PB), 130310.
29. H. M. Goda and P. Behrenbruch, Using a Modified Brooks-Corey Model to Study Oil-Water Relative Permeability for Diverse Pore Structures, in SPE Asia Pacific Oil and Gas Conference and Exhibition, SPE, 2004.
30. F. Lomeland, E. Ebeltoft and W. H. Thomas, A new versatile relative permeability correlation, Int. Symp. Soc. Core Anal., Toronto, Canada, 2005, pp. 1–12.
31. M. Shahbaz, T. Al-Ansari, M. Aslam, Z. Khan, A. Inayat and M. Athar, et al., A state of the art review on biomass processing and conversion technologies to produce hydrogen and its recovery via membrane separation, Int. J. Hydrogen Energy, 2020, 45(30), 15166–15195.
32. N. Pal and M. Agarwal, Advances in materials process and separation mechanism of the membrane towards hydrogen separation, Int. J. Hydrogen Energy, 2021, 46(53), 27062–27087.

---