# Synergies of storing hydrogen at the crest of CO2 or other gas storage

Sabrine Ben Rhouma, Salaheddine Chabab, Daniel Broseta

► To cite this version:

Sabrine Ben Rhouma, Salaheddine Chabab, Daniel Broseta. Synergies of storing hydrogen at the crest of CO2 or other gas storage. Greenhouse Gases: Science and Technology, 2024, 10.1002/ghg.2278. hal-04589273

HAL Id: hal-04589273

https://brgm.hal.science/hal-04589273

Submitted on 27 May 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_1.png" /></figure>

# Synergies of storing hydrogen at the crest of CO₂ or other gas storage

**Sabrine Ben Rhouma**, BRGM, 3 Ave Cl. Guillemin, Orleans, France

**Salaheddine Chabab**, Université de Pau et des Pays de l'AdourE2S UPPA, LaTEP, Pau, France and Université de Pau et des Pays de l'AdourE2S UPPA, CNRS, LFCR, Pau, France

**Daniel Broseta**, Université de Pau et des Pays de l'AdourE2S UPPA, CNRS, LFCR, Pau, France

**Abstract:** There are mutual benefits in storing H₂ in sedimentary reservoirs jointly with another gas serving as a cushion gas, such as the CO₂ of a carbon capture and storage (CCS) operation or the natural gas of seasonal storage or left in a depleted hydrocarbon reservoir. When H₂ occupies the crest of the reservoir, the presence of either gas is beneficial to the other. H₂ reinforces the sealing efficiency of the caprock due to its very favorable interfacial properties with respect to brine and rock-forming minerals. H₂ storage safety and capacity are also increased with cushion gases such as CO₂, which alleviate the buoyancy pressure at the top of the gas column. The potential drawback of this storage scheme is gas/gas mixing, which can, however, be strongly reduced if, by an appropriate choice of well completion and placement, H₂ is positioned in the upper zones of the reservoir, and its injection rate is kept below a critical value corresponding to the incipient fingering instability of the H₂/cushion gas mixing zone. This value, which depends on reservoir permeability, the dip angle of the mixing front, and how density varies with viscosity in the mixing front, turns out to be well above practical injection rates. Therefore, dispersive mixing is the only cause of front spreading, which is acceptable for not-too-heterogeneous reservoirs. The mutual benefits identified in this study are the strongest when the cushion is made up of dense CO₂, which suggests that the crest of offshore CO₂ storage reservoirs is a good candidate for H₂ storage. © 2024 The Authors. *Greenhouse Gases: Science and Technology* published by Society of Chemical Industry and John Wiley & Sons Ltd.

**Keywords:** caprock integrity; CO₂; geological storage; underground hydrogen storage; viscous fingering

## Introduction

With the increasing concern for reducing greenhouse gas emissions and balancing the seasonal production of renewable electricity, carbon capture and storage (CCS) and hydrogen energy are being considered as promising solutions. However, to fully utilize the potential of hydrogen, its large-scale storage becomes a necessity.¹ The

underground hydrogen storage (UHS) scheme examined in this work combines hydrogen storage with storage of another gas in the *same* sedimentary reservoir, for example, an aquifer or a depleted hydrocarbon reservoir. This other gas, which acts as a cushion gas (CG) for the stored H₂, could, for example, be the CO₂ captured from industrial gaseous effluents and permanently stored in a sedimentary reservoir as

Correspondence to: Sabrine Ben Rhouma, Bureau de Recherches Géologiques et Minières (BRGM), 3 Ave Cl. Guillemin, 45100 Orléans, France.

E-mail: s.benrhouma@brgm.fr

Received December 5, 2023; revised April 13, 2024; accepted April 26, 2024

Published online at Wiley Online Library (wileyonlinelibrary.com). DOI: 10.1002/ghg.2278

<figure><img src="image_2.png" /></figure>

This is an open access article under the terms of the Creative Commons Attribution-NonCommercial-NoDerivs License, which permits use and distribution in any medium, provided the original work is properly cited, the use is non-commercial and no modifications or adaptations are made.
<figure><img src="image_2.png" /></figure>

# Synergies of storing hydrogen at the crest of CO₂ or other gas storage

**Sabrine Ben Rhouma**, BRGM, 3 Ave Cl. Guillemin, Orleans, France

**Salaheddine Chabab**, Université de Pau et des Pays de l'AdourE2S UPPA, LaTEP, Pau, France and Université de Pau et des Pays de l'AdourE2S UPPA, CNRS, LFCR, Pau, France

**Daniel Broseta**, Université de Pau et des Pays de l'AdourE2S UPPA, CNRS, LFCR, Pau, France

**Abstract:** There are mutual benefits in storing H₂ in sedimentary reservoirs jointly with another gas serving as a cushion gas, such as the CO₂ of a carbon capture and storage (CCS) operation or the natural gas of seasonal storage or left in a depleted hydrocarbon reservoir. When H₂ occupies the crest of the reservoir, the presence of either gas is beneficial to the other. H₂ reinforces the sealing efficiency of the caprock due to its very favorable interfacial properties with respect to brine and rock-forming minerals. H₂ storage safety and capacity are also increased with cushion gases such as CO₂, which alleviate the buoyancy pressure at the top of the gas column. The potential drawback of this storage scheme is gas/gas mixing, which can, however, be strongly reduced if, by an appropriate choice of well completion and placement, H₂ is positioned in the upper zones of the reservoir, and its injection rate is kept below a critical value corresponding to the incipient fingering instability of the H₂/cushion gas mixing zone. This value, which depends on reservoir permeability, the dip angle of the mixing front, and how density varies with viscosity in the mixing front, turns out to be well above practical injection rates. Therefore, dispersive mixing is the only cause of front spreading, which is acceptable for not-too-heterogeneous reservoirs. The mutual benefits identified in this study are the strongest when the cushion is made up of dense CO₂, which suggests that the crest of offshore CO₂ storage reservoirs is a good candidate for H₂ storage. © 2024 The Authors. *Greenhouse Gases: Science and Technology* published by Society of Chemical Industry and John Wiley & Sons Ltd.

**Keywords:** caprock integrity; CO₂; geological storage; underground hydrogen storage; viscous fingering

## Introduction

With the increasing concern for reducing greenhouse gas emissions and balancing the seasonal production of renewable electricity, carbon capture and storage (CCS) and hydrogen energy are being considered as promising solutions. However, to fully utilize the potential of hydrogen, its large-scale storage becomes a necessity.¹ The

underground hydrogen storage (UHS) scheme examined in this work combines hydrogen storage with storage of another gas in the *same* sedimentary reservoir, for example, an aquifer or a depleted hydrocarbon reservoir. This other gas, which acts as a cushion gas (CG) for the stored H₂, could, for example, be the CO₂ captured from industrial gaseous effluents and permanently stored in a sedimentary reservoir as

Correspondence to: Sabrine Ben Rhouma, Bureau de Recherches Géologiques et Minières (BRGM), 3 Ave Cl. Guillemin, 45100 Orléans, France.

E-mail: s.benrhouma@brgm.fr

Received December 5, 2023; revised April 13, 2024; accepted April 26, 2024

Published online at Wiley Online Library (wileyonlinelibrary.com). DOI: 10.1002/ghg.2278

<figure><img src="image_3.png" /></figure>

This is an open access article under the terms of the Creative Commons Attribution-NonCommercial-NoDerivs License, which permits use and distribution in any medium, provided the original work is properly cited, the use is non-commercial and no modifications or adaptations are made.
part of a CCS operation, or those effluents themselves, or the natural gas of a seasonal gas storage or left in a depleted hydrocarbon reservoir, etc. The CG, usually the same gas and in a comparable amount as the working gas (WG; in this work, H₂), ensures acceptable H₂ withdrawal rates with minimal pressure drop in the reservoir and prevents brine production when this reservoir is (or is in contact with) an aquifer.²,³ The use of CGs different from the WG has already been examined, particularly in the context of natural gas storage⁴⁻⁷ and compressed air storage.⁸

Most of the gas storages or other gas-bearing reservoirs cited above are well characterized and possess most of the required infrastructure: their use to also store H₂ can be, therefore, an attractive option, provided mixing between the two gases is limited to avoid expensive treatment of the back-produced gas. In the case of a CO₂ cushion, the other reason for limiting its mixing with H₂ is to prevent (or limit) the growth of methanogenic and/or acetogenic microorganisms and the *in situ* production of methane and acetate, which occur in the mixing zone.⁹⁻¹²

A large amount of work is currently being devoted to the modeling, mostly by reservoir simulation, of H₂ storage in aquifers or depleted gas reservoirs with CGs such as CO₂, CH₄, or N₂ or their mixtures.¹³⁻¹⁸ The purpose is to help define the best storage scenarios, that is, scenarios that minimize mixing between H₂ and the CG (and thus the contamination of the back-produced H₂-rich gas) while maximizing recovery efficiency, that is, the quantity of H₂ back produced relative to that injected in consecutive H₂ injection/production sequences. At the same time, the pressure variations in the reservoir must lie within acceptable limits. Reservoir pressure should not exceed the limit for mechanical or capillary failure of the seal rock and too low a pressure means that the amount of H₂ stored in the reservoir is low as well. These processes, which control the performance and safety of the H₂ storage process, are sensitive to many parameters: the reservoir geological settings and petrophysical properties (particularly the permeability heterogeneities),¹⁸ the number, position, and completion of injection and production wells, and the CG type and composition.¹⁵,¹⁷ Most of this work relies on compositional reservoir simulation, in which the relevant fluid properties (density, phase composition, viscosity, interfacial tension [IFT], etc.) must be described accurately.¹⁶,¹⁷ Particular reservoirs or reservoir settings are considered in most studies, and it

may be difficult to disentangle the contributions of the various factors and draw definitive conclusions as to the physical mechanisms at play and their effects on the parameters of interest. These parameters strongly depend on the flow regimes, themselves controlled by reservoir heterogeneities and the various forces involved, namely the viscous, gravity, and capillary forces.¹⁴ The dispersive properties of the porous reservoir also have an important impact on the extent of mixing between the two gases.

This study examines the storage scheme depicted in Fig. 1, where H₂ is stored on top of a denser CG (e.g., CO₂), itself overlying an aquifer. The reservoir has top and lateral sealing boundaries, provided, for instance, by an anticline structure with a dome or crest at its top,¹⁹,²⁰ with one single vertical well allowing injection/withdrawal of the H₂ in the upper region of the reservoir. The gas/gas mixing front thus resembles the surface of a more or less deformed bubble that expands and shrinks around the well upon H₂ injection and withdrawal, respectively. This scheme is different from a horizontal (or slightly tilted) reservoir layer saturated with the CG beyond some distance from the vertical well, which injects and withdraws the WG.⁷,⁸

The issues addressed in this study are (i) mixing between the two gases and (ii) storage safety and capacity. On the one hand, H₂ is less dense, by one order of magnitude or more than most CGs, which stabilizes the gas/gas mixing front; but the viscosity contrast between the two fluids is such that viscous fingering and therefore enhanced mixing are expected above a certain injection rate. This rate or velocity limit depends on various factors such as the angle of inclination of the velocity with respect to the vertical and the viscosity and density contrasts of the two fluids,²¹,²² but it is not known whether this limit lies within practical values. This limit is determined and discussed for various H₂/CG couples and storage conditions in the next section. Below this limit, the only cause of mixing is dispersion between the two stored gases, also addressed in that section.

On the other hand, storage safety is ensured as long as the buoyancy pressure in the gas in contact with the seal rock does not exceed the capillary entry pressure of that gas into the seal rock.²³ The latter pressure is the minimum buoyancy pressure in the (nonwetting) gas phase (compared to brine pressure) required for the gas phase to intrude into the brine-saturated seal rock. It is reminded that water (or brine) is present as an irreducible phase that wets the pore walls throughout
part of a CCS operation, or those effluents themselves, or the natural gas of a seasonal gas storage or left in a depleted hydrocarbon reservoir, etc. The CG, usually the same gas and in a comparable amount as the working gas (WG; in this work, H₂), ensures acceptable H₂ withdrawal rates with minimal pressure drop in the reservoir and prevents brine production when this reservoir is (or is in contact with) an aquifer.²,³ The use of CGs different from the WG has already been examined, particularly in the context of natural gas storage⁴⁻⁷ and compressed air storage.⁸

Most of the gas storages or other gas-bearing reservoirs cited above are well characterized and possess most of the required infrastructure: their use to also store H₂ can be, therefore, an attractive option, provided mixing between the two gases is limited to avoid expensive treatment of the back-produced gas. In the case of a CO₂ cushion, the other reason for limiting its mixing with H₂ is to prevent (or limit) the growth of methanogenic and/or acetogenic microorganisms and the *in situ* production of methane and acetate, which occur in the mixing zone.⁹⁻¹²

A large amount of work is currently being devoted to the modeling, mostly by reservoir simulation, of H₂ storage in aquifers or depleted gas reservoirs with CGs such as CO₂, CH₄, or N₂ or their mixtures.¹³⁻¹⁸ The purpose is to help define the best storage scenarios, that is, scenarios that minimize mixing between H₂ and the CG (and thus the contamination of the back-produced H₂-rich gas) while maximizing recovery efficiency, that is, the quantity of H₂ back produced relative to that injected in consecutive H₂ injection/production sequences. At the same time, the pressure variations in the reservoir must lie within acceptable limits. Reservoir pressure should not exceed the limit for mechanical or capillary failure of the seal rock and too low a pressure means that the amount of H₂ stored in the reservoir is low as well. These processes, which control the performance and safety of the H₂ storage process, are sensitive to many parameters: the reservoir geological settings and petrophysical properties (particularly the permeability heterogeneities),¹⁸ the number, position, and completion of injection and production wells, and the CG type and composition.¹⁵,¹⁷ Most of this work relies on compositional reservoir simulation, in which the relevant fluid properties (density, phase composition, viscosity, interfacial tension [IFT], etc.) must be described accurately.¹⁶,¹⁷ Particular reservoirs or reservoir settings are considered in most studies, and it

may be difficult to disentangle the contributions of the various factors and draw definitive conclusions as to the physical mechanisms at play and their effects on the parameters of interest. These parameters strongly depend on the flow regimes, themselves controlled by reservoir heterogeneities and the various forces involved, namely the viscous, gravity, and capillary forces.¹⁴ The dispersive properties of the porous reservoir also have an important impact on the extent of mixing between the two gases.

This study examines the storage scheme depicted in Fig. 1, where H₂ is stored on top of a denser CG (e.g., CO₂), itself overlying an aquifer. The reservoir has top and lateral sealing boundaries, provided, for instance, by an anticline structure with a dome or crest at its top,¹⁹,²⁰ with one single vertical well allowing injection/withdrawal of the H₂ in the upper region of the reservoir. The gas/gas mixing front thus resembles the surface of a more or less deformed bubble that expands and shrinks around the well upon H₂ injection and withdrawal, respectively. This scheme is different from a horizontal (or slightly tilted) reservoir layer saturated with the CG beyond some distance from the vertical well, which injects and withdraws the WG.⁷,⁸

The issues addressed in this study are (i) mixing between the two gases and (ii) storage safety and capacity. On the one hand, H₂ is less dense, by one order of magnitude or more than most CGs, which stabilizes the gas/gas mixing front; but the viscosity contrast between the two fluids is such that viscous fingering and therefore enhanced mixing are expected above a certain injection rate. This rate or velocity limit depends on various factors such as the angle of inclination of the velocity with respect to the vertical and the viscosity and density contrasts of the two fluids,²¹,²² but it is not known whether this limit lies within practical values. This limit is determined and discussed for various H₂/CG couples and storage conditions in the next section. Below this limit, the only cause of mixing is dispersion between the two stored gases, also addressed in that section.

On the other hand, storage safety is ensured as long as the buoyancy pressure in the gas in contact with the seal rock does not exceed the capillary entry pressure of that gas into the seal rock.²³ The latter pressure is the minimum buoyancy pressure in the (nonwetting) gas phase (compared to brine pressure) required for the gas phase to intrude into the brine-saturated seal rock. It is reminded that water (or brine) is present as an irreducible phase that wets the pore walls throughout
<figure><img src="image_3.png" /><figcaption>**Figure 1.** Schematics of an H₂ injection and storage on top of a cushion gas (CG) in an aquifer. Left: reservoir/well configuration, with the CG overlying an aquifer, and the H₂ “bubble” being injected and pushing the CG, with velocity U at the gas/gas mixing front, whose dip or inclination with respect to the horizontal is characterized by the angle α. The vertical scale is different from the horizontal scale: the dome structure looks exaggerated. Right: water and gas hydrostatic pressure variations along the vertical; the buoyancy pressure, which depends on the densities and column heights of the two gases, is maximum at the crest of the reservoir.</figcaption></figure>

the gas column and connects the seal rock (caprock) to the underlying aquifer. The pressure profiles in both phases obey the law of hydrostatics (Fig. 1) – which also holds for low-to-moderate flow rates, an approximation known as the vertical equilibrium approximation. The capillary entry pressure depends on the pore structure of the caprock and on the interfacial properties of the seal rock/brine/gas system, which have been thoroughly investigated for the seal rock components and the gases of interest, including H₂, CO₂, CH₄, N₂ and some of their mixtures. These investigations, briefly reviewed in the third section, clearly demonstrate the very favorable interfacial properties of seal rock/brine/H₂ systems, and, therefore, adding H₂ at the top of a gas column renders the seal rock less prone to capillary failure while allowing some amount of H₂ to be safely stored in addition to the gas already present. In addition, the presence of a dense CG such as CO₂ at the bottom of the gas column alleviates the buoyancy pressure at its top as compared to a gas column made up of H₂ only, thus allowing a larger height (or quantity) of working H₂ to be stored safely. The calculation procedure and illustrative examples of those storage heights and quantities are given in the third section.

The last section summarizes the results and implications of those findings and discusses perspectives for H₂ storage in CO₂- or other gas-bearing reservoirs.

## Spreading of the gas-gas mixing zone

When H₂ is withdrawn from the reservoir, both gravity forces (i.e., the strong density contrast between H₂ and the CG) and viscous forces (i.e., the favorable viscosity contrast) stabilize the mixing front, irrespective of fluid velocity: the mixing front is unconditionally stable. The mixing zone spreads only under the effect of dispersion, which is presented and discussed in the next subsection. However, when H₂ is injected into the reservoir, the two forces have opposite effects: gravity stabilizes the front but not the viscous forces as the displacing fluid, H₂, is more mobile (less viscous) than the displaced CG. The mixing front is thus conditionally stable, that is, it is stabilized by gravity at a low enough velocity, but prone to viscous fingering above a certain injection velocity. This velocity, which corresponds to the incipient fingering instability of the front, is evaluated in the subsection ‘Rate conditions
<figure><img src="image_4.png" /><figcaption>**Figure 1.** Schematics of an H₂ injection and storage on top of a cushion gas (CG) in an aquifer. Left: reservoir/well configuration, with the CG overlying an aquifer, and the H₂ “bubble” being injected and pushing the CG, with velocity U at the gas/gas mixing front, whose dip or inclination with respect to the horizontal is characterized by the angle α. The vertical scale is different from the horizontal scale: the dome structure looks exaggerated. Right: water and gas hydrostatic pressure variations along the vertical; the buoyancy pressure, which depends on the densities and column heights of the two gases, is maximum at the crest of the reservoir.</figcaption></figure>

the gas column and connects the seal rock (caprock) to the underlying aquifer. The pressure profiles in both phases obey the law of hydrostatics (Fig. 1) – which also holds for low-to-moderate flow rates, an approximation known as the vertical equilibrium approximation. The capillary entry pressure depends on the pore structure of the caprock and on the interfacial properties of the seal rock/brine/gas system, which have been thoroughly investigated for the seal rock components and the gases of interest, including H₂, CO₂, CH₄, N₂ and some of their mixtures. These investigations, briefly reviewed in the third section, clearly demonstrate the very favorable interfacial properties of seal rock/brine/H₂ systems, and, therefore, adding H₂ at the top of a gas column renders the seal rock less prone to capillary failure while allowing some amount of H₂ to be safely stored in addition to the gas already present. In addition, the presence of a dense CG such as CO₂ at the bottom of the gas column alleviates the buoyancy pressure at its top as compared to a gas column made up of H₂ only, thus allowing a larger height (or quantity) of working H₂ to be stored safely. The calculation procedure and illustrative examples of those storage heights and quantities are given in the third section.

The last section summarizes the results and implications of those findings and discusses perspectives for H₂ storage in CO₂- or other gas-bearing reservoirs.

## Spreading of the gas-gas mixing zone

When H₂ is withdrawn from the reservoir, both gravity forces (i.e., the strong density contrast between H₂ and the CG) and viscous forces (i.e., the favorable viscosity contrast) stabilize the mixing front, irrespective of fluid velocity: the mixing front is unconditionally stable. The mixing zone spreads only under the effect of dispersion, which is presented and discussed in the next subsection. However, when H₂ is injected into the reservoir, the two forces have opposite effects: gravity stabilizes the front but not the viscous forces as the displacing fluid, H₂, is more mobile (less viscous) than the displaced CG. The mixing front is thus conditionally stable, that is, it is stabilized by gravity at a low enough velocity, but prone to viscous fingering above a certain injection velocity. This velocity, which corresponds to the incipient fingering instability of the front, is evaluated in the subsection ‘Rate conditions
for viscous fingering avoidance (front stability)'. Below this velocity, the mixing front is stable and subject to dispersive mixing only, just as in withdrawal sequences. This information can be used to optimize the injection/withdrawal strategy and minimize the negative effects of gas-gas mixing.

## Front spreading by dispersive mixing

Dispersive mixing is, strictly speaking, defined for two fluids that differ only by the presence or absence of tracer molecules; this is referred to as tracer dispersion. It is also used for two miscible phases that differ in density and viscosity as long as the composition profile in the mixing zone can be described by a convection-dispersion differential equation. This seems to be the case in practice whenever the front is stable. Early studies,24–26 mostly carried out with miscible liquid/liquid systems, indicate that longitudinal dispersion coefficients $K_L$ are similar or even reduced (relative to tracer dispersion coefficients) when the two fluids have different densities and viscosities.

These coefficients depend on the interstitial velocity $V$ of the fluids as follows:

$$ \frac{K_L}{D} = \frac{1}{\tau} + \left( \frac{\alpha_L V}{D} \right)^m, \quad (1) $$

where $D$ is the intermolecular diffusion coefficient between the two fluids, $\tau$ is the tortuosity of the porous medium, and $\alpha_L$ is the longitudinal dispersivity, a length scale of velocity variations closely related to permeability heterogeneities. The exponent $m$, experimentally found to be very slightly larger than 1, is often considered equal to 1. There are two contributions to dispersive mixing. One, $D/\tau$, is dominant at low-velocity $V$ and is related to molecular diffusion, $\tau > 1$ reflecting the fact that fluid motion is hindered by the porous substrate: $\tau \approx 1.4$ in granular media but several units in consolidated porous media such as Berea sandstones.27 The other term is the hydrodynamic (or mechanical) contribution, which reflects variations in the velocity field, themselves related to permeability heterogeneities. It is dominant for Peclet numbers $\alpha_L V/D$ exceeding a few units.

$D$, which is not larger than $1 \text{ cm}^2/\text{s}$ under near-ambient conditions,28,29 strongly decreases with increasing pressure $P$, $D \approx C/P$, where $C$ is a constant usually below unity ($D$ in $\text{cm}^2/\text{s}$, $P$ in bar).29 As to dispersivities $\alpha_L$, low values are measured in the

laboratory: the recent careful core flood experiments conducted by Yang and co-workers27 with $H_2$ and $N_2$ indicate values as low as $\alpha_L = 0.3$ mm in a (homogeneous) Berea sandstone. Field-scale dispersivities are much higher, reflecting permeability heterogeneities – values in the meter range are not uncommon.30 In practical situations of gas injection or withdrawal, the Peclet number is well above 1: dispersive mixing is dominated by the mechanical dispersion term in Eqn 1.

Dispersivities can and should be directly measured in the field, for example, by carrying out push-pull tracer tests around the well. Assuming that $H_2$ is withdrawn or injected with downhole velocities $V$ and $-V$, respectively, and that the front is horizontal and moves in the vertical direction, the extension $Z(t)$ of the $H_2$/CG mixing zone is expected to increase with the time $t$ spent during those injection and withdrawal sequences as follows: If there are no shut-in periods, the spreading length $Z(t)$ is roughly proportional to the square root of the cumulative quantity of $H_2$ circulated (i.e., the sum of the injected and withdrawn $H_2$, both counted positively) over the time $t$, which is proportional to $|V|t$. For continuous injection and withdrawal operations without shut-in periods, and with typical practical values $|V| \approx 1 \text{ cm/day}$ and $\alpha_L \approx 10^{-3}$, a few tens of centimeters, the extent $Z(t)$ of the $H_2$/CG mixing front following 1 year of injection/withdrawal sequences thus lies in the meter range.

## Rate conditions for viscous fingering avoidance (front stability)

We first recall and discuss the classical result for *immiscible* displacements in porous media,21,22 which can be applied for $H_2$ displacing brine in an aquifer.31 We then focus on the *miscible* displacements of interest. It is important here to emphasize that the fingering instability is not captured by conventional reservoir simulations unless some noise is introduced in the permeability field and very fine grid blocks are used.32 In other words, a stable front – that is, the absence of fingers – is *systematically* observed in conventional reservoir simulations of displacements in homogeneous porous media, irrespective of the flow rate.

In a downward displacement of a fluid by an *immiscible* less dense and less viscous fluid, fingering occurs at the front between the two fluids when the
for viscous fingering avoidance (front stability)'. Below this velocity, the mixing front is stable and subject to dispersive mixing only, just as in withdrawal sequences. This information can be used to optimize the injection/withdrawal strategy and minimize the negative effects of gas-gas mixing.

## Front spreading by dispersive mixing

Dispersive mixing is, strictly speaking, defined for two fluids that differ only by the presence or absence of tracer molecules; this is referred to as tracer dispersion. It is also used for two miscible phases that differ in density and viscosity as long as the composition profile in the mixing zone can be described by a convection-dispersion differential equation. This seems to be the case in practice whenever the front is stable. Early studies,24–26 mostly carried out with miscible liquid/liquid systems, indicate that longitudinal dispersion coefficients $K_L$ are similar or even reduced (relative to tracer dispersion coefficients) when the two fluids have different densities and viscosities.

These coefficients depend on the interstitial velocity $V$ of the fluids as follows:

$$ \frac{K_L}{D} = \frac{1}{\tau} + \left( \frac{\alpha_L V}{D} \right)^m, \quad (1) $$

where $D$ is the intermolecular diffusion coefficient between the two fluids, $\tau$ is the tortuosity of the porous medium, and $\alpha_L$ is the longitudinal dispersivity, a length scale of velocity variations closely related to permeability heterogeneities. The exponent $m$, experimentally found to be very slightly larger than 1, is often considered equal to 1. There are two contributions to dispersive mixing. One, $D/\tau$, is dominant at low-velocity $V$ and is related to molecular diffusion, $\tau > 1$ reflecting the fact that fluid motion is hindered by the porous substrate: $\tau \approx 1.4$ in granular media but several units in consolidated porous media such as Berea sandstones.27 The other term is the hydrodynamic (or mechanical) contribution, which reflects variations in the velocity field, themselves related to permeability heterogeneities. It is dominant for Peclet numbers $\alpha_L V/D$ exceeding a few units.

$D$, which is not larger than $1 \text{ cm}^2/\text{s}$ under near-ambient conditions,28,29 strongly decreases with increasing pressure $P$, $D \approx C/P$, where $C$ is a constant usually below unity ($D$ in $\text{cm}^2/\text{s}$, $P$ in bar).29 As to dispersivities $\alpha_L$, low values are measured in the

laboratory: the recent careful core flood experiments conducted by Yang and co-workers27 with $H_2$ and $N_2$ indicate values as low as $\alpha_L = 0.3$ mm in a (homogeneous) Berea sandstone. Field-scale dispersivities are much higher, reflecting permeability heterogeneities – values in the meter range are not uncommon.30 In practical situations of gas injection or withdrawal, the Peclet number is well above 1: dispersive mixing is dominated by the mechanical dispersion term in Eqn 1.

Dispersivities can and should be directly measured in the field, for example, by carrying out push-pull tracer tests around the well. Assuming that $H_2$ is withdrawn or injected with downhole velocities $V$ and $-V$, respectively, and that the front is horizontal and moves in the vertical direction, the extension $Z(t)$ of the $H_2$/CG mixing zone is expected to increase with the time $t$ spent during those injection and withdrawal sequences as follows: If there are no shut-in periods, the spreading length $Z(t)$ is roughly proportional to the square root of the cumulative quantity of $H_2$ circulated (i.e., the sum of the injected and withdrawn $H_2$, both counted positively) over the time $t$, which is proportional to $|V|t$. For continuous injection and withdrawal operations without shut-in periods, and with typical practical values $|V| \approx 1 \text{ cm/day}$ and $\alpha_L \approx 10^{-3}$, a few tens of centimeters, the extent $Z(t)$ of the $H_2$/CG mixing front following 1 year of injection/withdrawal sequences thus lies in the meter range.

## Rate conditions for viscous fingering avoidance (front stability)

We first recall and discuss the classical result for *immiscible* displacements in porous media,21,22 which can be applied for $H_2$ displacing brine in an aquifer.31 We then focus on the *miscible* displacements of interest. It is important here to emphasize that the fingering instability is not captured by conventional reservoir simulations unless some noise is introduced in the permeability field and very fine grid blocks are used.32 In other words, a stable front – that is, the absence of fingers – is *systematically* observed in conventional reservoir simulations of displacements in homogeneous porous media, irrespective of the flow rate.

In a downward displacement of a fluid by an *immiscible* less dense and less viscous fluid, fingering occurs at the front between the two fluids when the
Darcy (or superficial) velocity normal to the front exceeds the following critical value²¹

$$U_{im} = \frac{g\Delta\rho}{\Delta(1/\lambda)}\cos(\alpha), \quad (2)$$

where g is the acceleration due to gravity, α the dip angle of the mixing front (see Fig. 1), Δρ the density difference between the displaced and displacing phases, Δ(1/λ) the difference between the inverse mobilities on both sides of the front. The Darcy (or superficial) velocity U is equal to the interstitial (or front) velocity V multiplied by the porosity of the porous medium.

When gravity effects become negligible (Δρcos(α) decreases to 0), the sign of the denominator in Eqn 2 or, equivalently, the position with respect to 1 of the mobility ratio M<sub>r</sub> = λ<sub>displacing</sub>/λ<sub>displaced</sub>, controls whether there is viscous fingering (when this ratio is larger than 1) or a stable, piston-like front (the ratio is smaller than 1), irrespective of velocity. This is the classical Saffmann–Taylor result. α > 90° corresponds to a front with a fluid, here H<sub>2</sub>, pushing upwardly a denser and more viscous CG: the front is unconditionally unstable, that is, both gravity and viscous forces favor an unstable front, irrespective of front velocity. This situation, not depicted in Fig. 1, can be overlooked if the injection occurs near the crest of the reservoir.

For H<sub>2</sub> being injected at the top of an aquifer and thus downwardly displacing brine in a porous rock, the density and viscosity contrasts are such that Δρ ≈ ρ<sub>w</sub>, the brine density, and Δ(1/λ) = 1/λ<sub>w</sub> - 1/λ<sub>H2</sub> ≈ μ<sub>w</sub>/k<sub>w</sub>, where μ<sub>w</sub> is brine viscosity and k<sub>w</sub> = k<sub>rw</sub> the brine permeability on the brine-rich side of the front (k is the absolute permeability of the medium and k<sub>rw</sub> the brine relative permeability on the brine-rich side of the front). The above equation thus reduces to³¹

$$U_{im} = \frac{gk_{rw}\rho_{w}}{\mu_{w}}\cos(\alpha). \quad (3)$$

The stability criterion U < U<sub>im</sub> can be seen as a condition on the ratio of viscous forces μ<sub>w</sub>U/k<sub>rw</sub> to gravity or buoyancy forces gΔρcos(α). This ratio must be larger than unity for the displacement front to be stable, where unity corresponds to the stability limit U<sub>im</sub>. As an example, critical velocities of H<sub>2</sub> downwardly displacing water ranging from U<sub>im</sub> = a few to tens of centimeters per day are inferred from the permeabilities measured by Yekta et al.³³ in a k = 50 mdarcy sandstone under various (shallow and deep)

storage conditions. Smaller values are indeed expected for nonhorizontal fronts (0 < α < 90°).

The above classical result<sup>21,22</sup> is obtained by a linear stability analysis, in which the pressure field is calculated across a small protrusion at the front between the displacing (upper) and displaced (lower) phases. In the stable situation (no fingering), the pressure gradient is smaller in the upper phase than in the lower phase, whereas in the unstable situation (with fingering) the pressure gradient is larger. The critical velocity U<sub>im</sub> is determined by equating both gradients. If the velocity is below that critical value, the front is stable, and if it is above the front is unstable and fingering occurs. This analysis also provides the width or peak-to-peak distance of the incipient fingers, which increases with brine/water IFT and the proximity of U to the critical value U<sub>im</sub>.³¹

If the two phases are miscible but the gas/gas mixing front is sharp, Eqn 2 can be used, which reduces to

$$U_{im} = kg\frac{\Delta\rho}{\Delta\mu}\cos(\alpha), \quad (4)$$

where Δρ = ρ<sub>CG</sub> - ρ<sub>H2</sub> is the density difference and Δμ = μ<sub>CG</sub> - μ<sub>H2</sub> the viscosity difference between the two gases.

Indeed the front is not sharp: the fluid composition gradually changes over the mixing zone, from pure H<sub>2</sub> at the top to pure CG at the bottom. Applying the above reasoning to any porous layer of thickness dz, a more stringent stability condition is obtained, that is, U < U<sub>m</sub> with

$$U_m = kg\cos(\alpha)\left(\frac{d\rho}{d\mu}\right)_{min} = kg\cos(\alpha)\left(\frac{d\mu}{d\rho}\right)_{max}, \quad (5)$$

where (dρ/dμ)min is the minimum value of dρ/dμ as the fluid composition varies over the gas-gas mixing zone. This minimum is reached for a particular composition, that is, at a particular position in the mixing front, which is where the instability is triggered and fingers start to develop.

This condition was first derived in 1964 by Dumoré,³⁴ who applied it to miscible liquid/liquid (oil/solvent) displacements. Using simple models for the viscosity and density of oil/solvent mixtures, Dumoré observed that this minimum is reached when the fraction of the most viscous (and denser) liquid, that is, oil, is maximum, meaning that the onset of instability occurs in the bottom part of the miscible front; in addition, he observed rather low values of U<sub>m</sub> because contrasted liquid/liquid mixtures typically
Darcy (or superficial) velocity normal to the front exceeds the following critical value²¹

$$U_{im} = \frac{g\Delta\rho}{\Delta(1/\lambda)}\cos(\alpha), \quad (2)$$

where g is the acceleration due to gravity, α the dip angle of the mixing front (see Fig. 1), Δρ the density difference between the displaced and displacing phases, Δ(1/λ) the difference between the inverse mobilities on both sides of the front. The Darcy (or superficial) velocity U is equal to the interstitial (or front) velocity V multiplied by the porosity of the porous medium.

When gravity effects become negligible (Δρcos(α) decreases to 0), the sign of the denominator in Eqn 2 or, equivalently, the position with respect to 1 of the mobility ratio M<sub>r</sub> = λ<sub>displacing</sub>/λ<sub>displaced</sub>, controls whether there is viscous fingering (when this ratio is larger than 1) or a stable, piston-like front (the ratio is smaller than 1), irrespective of velocity. This is the classical Saffmann–Taylor result. α > 90° corresponds to a front with a fluid, here H<sub>2</sub>, pushing upwardly a denser and more viscous CG: the front is unconditionally unstable, that is, both gravity and viscous forces favor an unstable front, irrespective of front velocity. This situation, not depicted in Fig. 1, can be overlooked if the injection occurs near the crest of the reservoir.

For H<sub>2</sub> being injected at the top of an aquifer and thus downwardly displacing brine in a porous rock, the density and viscosity contrasts are such that Δρ ≈ ρ<sub>w</sub>, the brine density, and Δ(1/λ) = 1/λ<sub>w</sub> - 1/λ<sub>H₂</sub> ≈ μ<sub>w</sub>/k<sub>w</sub>, where μ<sub>w</sub> is brine viscosity and k<sub>w</sub> = k<sub>rw</sub> the brine permeability on the brine-rich side of the front (k is the absolute permeability of the medium and k<sub>rw</sub> the brine relative permeability on the brine-rich side of the front). The above equation thus reduces to³¹

$$U_{im} = \frac{gk_{rw}\rho_{w}}{\mu_{w}}\cos(\alpha). \quad (3)$$

The stability criterion U < U<sub>im</sub> can be seen as a condition on the ratio of viscous forces μ<sub>w</sub>U/k<sub>rw</sub> to gravity or buoyancy forces gΔρcos(α). This ratio must be larger than unity for the displacement front to be stable, where unity corresponds to the stability limit U<sub>im</sub>. As an example, critical velocities of H<sub>2</sub> downwardly displacing water ranging from U<sub>im</sub> = a few to tens of centimeters per day are inferred from the permeabilities measured by Yekta et al.³³ in a k = 50 mdarcy sandstone under various (shallow and deep)

storage conditions. Smaller values are indeed expected for nonhorizontal fronts (0 < α < 90°).

The above classical result<sup>21,22</sup> is obtained by a linear stability analysis, in which the pressure field is calculated across a small protrusion at the front between the displacing (upper) and displaced (lower) phases. In the stable situation (no fingering), the pressure gradient is smaller in the upper phase than in the lower phase, whereas in the unstable situation (with fingering) the pressure gradient is larger. The critical velocity U<sub>im</sub> is determined by equating both gradients. If the velocity is below that critical value, the front is stable, and if it is above the front is unstable and fingering occurs. This analysis also provides the width or peak-to-peak distance of the incipient fingers, which increases with brine/water IFT and the proximity of U to the critical value U<sub>im</sub>.³¹

If the two phases are miscible but the gas/gas mixing front is sharp, Eqn 2 can be used, which reduces to

$$U_{im} = kg\frac{\Delta\rho}{\Delta\mu}\cos(\alpha), \quad (4)$$

where Δρ = ρ<sub>CG</sub> - ρ<sub>H₂</sub> is the density difference and Δμ = μ<sub>CG</sub> - μ<sub>H₂</sub> the viscosity difference between the two gases.

Indeed the front is not sharp: the fluid composition gradually changes over the mixing zone, from pure H<sub>2</sub> at the top to pure CG at the bottom. Applying the above reasoning to any porous layer of thickness dz, a more stringent stability condition is obtained, that is, U < U<sub>m</sub> with

$$U_m = kg\cos(\alpha)\left(\frac{d\rho}{d\mu}\right)_{min} = kg\cos(\alpha)\left(\frac{d\mu}{d\rho}\right)_{max}, \quad (5)$$

where (dρ/dμ)min is the minimum value of dρ/dμ as the fluid composition varies over the gas-gas mixing zone. This minimum is reached for a particular composition, that is, at a particular position in the mixing front, which is where the instability is triggered and fingers start to develop.

This condition was first derived in 1964 by Dumoré,³⁴ who applied it to miscible liquid/liquid (oil/solvent) displacements. Using simple models for the viscosity and density of oil/solvent mixtures, Dumoré observed that this minimum is reached when the fraction of the most viscous (and denser) liquid, that is, oil, is maximum, meaning that the onset of instability occurs in the bottom part of the miscible front; in addition, he observed rather low values of U<sub>m</sub> because contrasted liquid/liquid mixtures typically
exhibit large viscosity variations while density changes only slightly. Later stages and the nonlinear evolution of the instability have been investigated by means of various theoretical and numerical approaches,³⁵ in which the density and viscosity profiles are shown to play an important role. The recent X-ray CT imaging observations of downward miscible displacements in a porous medium conducted under various velocities by Suekane and co-workers³⁶ are consistent with Dumoré's equation (5). In a configuration where the liquid at the top displaces the other denser, more viscous liquid at the bottom, these authors could observe a stable front at low velocity (U < U_m) and an unstable front with complex nonlinear interactions between fingers (tip-splitting, coalescence, shielding, and so on) at high velocity (U > U_m). Hereafter, we limit our analysis to the initiation of viscous fingering.

Considering that $d\rho/d\mu$ is equal to $d\rho/dy$ divided by $d\mu/dy$, where y is the H₂ molar content of the H₂/cushion gas mixture of interest, appropriate models are needed for the density and viscosity of H₂/CG mixtures as a function of y, temperature T, and pressure P. We consider a range of T and P from 298 to 373 K and 5 to 25 MPa representative of various storage conditions, from shallow cold (offshore) reservoirs to deep and hot (onshore) reservoirs.

We start by deriving those two quantities and the resulting limit velocity $U_m$ within the low-density approximation, which provides fairly accurate predictions for all H₂/CG mixtures of interest but H₂/CO₂ mixtures, as discussed in the Appendix and later in this section.

The low-density limit for the density of H₂/CG mixtures is nothing but the ideal gas (IG) law:

$$ \rho = \frac{P(yM_{H_2} + (1 - y)M_{CG})}{RT}, \qquad (6) $$

where $M_{H_2} = 2 \times 10^{-3}$ kg/mol is the molar weight of H₂ and $M_{CG}$ is the molar weight of the cushion gas, $R = 8.3145$ J/mol·K is the gas constant, P is the pressure (in Pascal), and T is the temperature (in kelvin). This approximation is equivalent to assuming a compressibility factor equal to 1, which is true to within less than 10% for H₂, CH₄, and N₂ and their mixtures within the T and P ranges of interest.

The low-density limit for the viscosity $\mu$ of gas mixtures as a function of the viscosities of their constituents, here $\mu_{H_2}$ and $\mu_{CG}$, can be found in standard textbooks²⁹

$$ \mu = \frac{y\mu_{H_2}}{y + (1 - y)\phi_{H_2,CG}} + \frac{(1 - y)\mu_{CG}}{1 - y + y\phi_{CG,H_2}}. \qquad (7) $$

where the exchange term

$$ \Phi_{CG,H_2} = \frac{\left[ 1 + M_r^{1/2} \left( \frac{M_{H_2}}{M_{CG}} \right)^{\frac{1}{4}} \right]^2}{\left[ 8 \left( 1 + \frac{M_{CG}}{M_{H_2}} \right)^{\frac{1}{2}} \right]^2}. \qquad (8) $$

These equations are those of the original Wilke's model.³⁷ An alternative exchange term has been proposed by Herning and Zipperer (HZ),³⁸

$$ \Phi_{CG,H_2} = \left( \frac{M_{H_2}}{M_{CG}} \right)^{1/2} \qquad (9) $$

where $\Phi_{H_2,CG}$ is obtained by interchanging CG and H₂ in the above two expressions and replacing the mobility ratio $M_r$ by $1/M_r$ in Eqn 8. $M_r = \mu_{CG}/\mu_{H_2}$ spans the following intervals: 1.35–2.1 (for CG = CH₄, 1.7–1.9 (for CG = gaseous CO₂), and 2.1–2.5 (for CG = N₂) in the T and P ranges 298–373 K and 5–25 MPa. $M_r > 1$ means that there is a high potential for viscous fingering in configurations where gravity can be neglected ($\alpha \approx 90^{\circ}$).

An analytical expression is obtained for $d\rho/d\mu$ from Eqns 6–8, which is a monotonically decreasing function of y. The minimum value of $d\rho/d\mu$ is obtained for $y \to 1$, meaning that the upper, H₂-rich region of the mixing zone is where fingering is triggered when the velocity reaches or exceeds the following critical value:

$$ U_m = \frac{kgP\cos(\alpha)}{RT} \frac{M_{CG} - M_{H_2}}{\mu_{CG}/\phi_{CG,H_2} - \mu_{H_2}\phi_{H_2,CG}}, \qquad (10) $$

where $U_m$ is a fraction of the value in the absence of mixing (Eqn 4):

$$ U_m/U_{im} = \frac{M_r - 1}{M_r/\phi_{CG,H_2} - \phi_{H_2,CG}}. \qquad (11) $$

This ratio depends only on the densities, molar weights, and viscosities of the two gases, here H₂ and the CG. Equations 10 and 11 are applicable only with CGs whose densities can be approximated by the IG law (Eqn 6) and whose viscosities are obtained by Wilke's equation (7), otherwise Eqn 5 must be used in combination with appropriate models for the mixture's viscosity and density – see the Appendix. The above approach and expressions can easily be extended to any couple of WG and CG that fit the IG law and Wilke's equation for their densities and viscosities, by
exhibit large viscosity variations while density changes only slightly. Later stages and the nonlinear evolution of the instability have been investigated by means of various theoretical and numerical approaches,³⁵ in which the density and viscosity profiles are shown to play an important role. The recent X-ray CT imaging observations of downward miscible displacements in a porous medium conducted under various velocities by Suekane and co-workers³⁶ are consistent with Dumoré's equation (5). In a configuration where the liquid at the top displaces the other denser, more viscous liquid at the bottom, these authors could observe a stable front at low velocity (U < U_m) and an unstable front with complex nonlinear interactions between fingers (tip-splitting, coalescence, shielding, and so on) at high velocity (U > U_m). Hereafter, we limit our analysis to the initiation of viscous fingering.

Considering that $d\rho/d\mu$ is equal to $d\rho/dy$ divided by $d\mu/dy$, where y is the H₂ molar content of the H₂/cushion gas mixture of interest, appropriate models are needed for the density and viscosity of H₂/CG mixtures as a function of y, temperature T, and pressure P. We consider a range of T and P from 298 to 373 K and 5 to 25 MPa representative of various storage conditions, from shallow cold (offshore) reservoirs to deep and hot (onshore) reservoirs.

We start by deriving those two quantities and the resulting limit velocity $U_m$ within the low-density approximation, which provides fairly accurate predictions for all H₂/CG mixtures of interest but H₂/CO₂ mixtures, as discussed in the Appendix and later in this section.

The low-density limit for the density of H₂/CG mixtures is nothing but the ideal gas (IG) law:

$$ \rho = \frac{P(yM_{H_2} + (1 - y)M_{CG})}{RT}, \quad (6) $$

where $M_{H_2} = 2 \times 10^{-3}$ kg/mol is the molar weight of H₂ and $M_{CG}$ is the molar weight of the cushion gas, $R = 8.3145$ J/mol·K is the gas constant, P is the pressure (in Pascal), and T is the temperature (in kelvin). This approximation is equivalent to assuming a compressibility factor equal to 1, which is true to within less than 10% for H₂, CH₄, and N₂ and their mixtures within the T and P ranges of interest.

The low-density limit for the viscosity $\mu$ of gas mixtures as a function of the viscosities of their constituents, here $\mu_{H_2}$ and $\mu_{CG}$, can be found in standard textbooks²⁹

$$ \mu = \frac{y\mu_{H_2}}{y + (1 - y)\phi_{H_2,CG}} + \frac{(1 - y)\mu_{CG}}{1 - y + y\phi_{CG,H_2}}. \quad (7) $$

where the exchange term

$$ \Phi_{CG,H_2} = \frac{\left[ 1 + M_r^{1/2} \left( \frac{M_{H_2}}{M_{CG}} \right)^{\frac{1}{4}} \right]^2}{\left[ 8 \left( 1 + \frac{M_{CG}}{M_{H_2}} \right)^{\frac{1}{2}} \right]^2}. \quad (8) $$

These equations are those of the original Wilke's model.³⁷ An alternative exchange term has been proposed by Herning and Zipperer (HZ),³⁸

$$ \Phi_{CG,H_2} = \left( \frac{M_{H_2}}{M_{CG}} \right)^{1/2} \quad (9) $$

where $\Phi_{H_2,CG}$ is obtained by interchanging CG and H₂ in the above two expressions and replacing the mobility ratio $M_r$ by $1/M_r$ in Eqn 8. $M_r = \mu_{CG}/\mu_{H_2}$ spans the following intervals: 1.35–2.1 (for CG = CH₄, 1.7–1.9 (for CG = gaseous CO₂), and 2.1–2.5 (for CG = N₂) in the T and P ranges 298–373 K and 5–25 MPa. $M_r > 1$ means that there is a high potential for viscous fingering in configurations where gravity can be neglected ($\alpha \approx 90^{\circ}$).

An analytical expression is obtained for $d\rho/d\mu$ from Eqns 6–8, which is a monotonically decreasing function of y. The minimum value of $d\rho/d\mu$ is obtained for $y \to 1$, meaning that the upper, H₂-rich region of the mixing zone is where fingering is triggered when the velocity reaches or exceeds the following critical value:

$$ U_m = \frac{kgP\cos(\alpha)}{RT} \frac{M_{CG} - M_{H_2}}{\mu_{CG}/\phi_{CG,H_2} - \mu_{H_2}\phi_{H_2,CG}}, \quad (10) $$

where $U_m$ is a fraction of the value in the absence of mixing (Eqn 4):

$$ U_m/U_{im} = \frac{M_r - 1}{M_r/\phi_{CG,H_2} - \phi_{H_2,CG}}. \quad (11) $$

This ratio depends only on the densities, molar weights, and viscosities of the two gases, here H₂ and the CG. Equations 10 and 11 are applicable only with CGs whose densities can be approximated by the IG law (Eqn 6) and whose viscosities are obtained by Wilke's equation (7), otherwise Eqn 5 must be used in combination with appropriate models for the mixture's viscosity and density – see the Appendix. The above approach and expressions can easily be extended to any couple of WG and CG that fit the IG law and Wilke's equation for their densities and viscosities, by
<table><caption>**Table 1.** Values of *U*<sub>m</sub>/*U*<sub>im</sub> calculated within the low-density (IG+Wilke) approximation for H<sub>2</sub> downwardly displacing CH<sub>4</sub> or N<sub>2</sub>.</caption><tr><td>Temperature (K)</td><td>298</td><td>298</td><td>373</td><td>373</td></tr><tr><td>Pressure (MPa)</td><td>5</td><td>25</td><td>5</td><td>25</td></tr><tr><td>H<sub>2</sub> + CH<sub>4</sub></td><td>0.173</td><td>0.340</td><td>0.175</td><td>0.268</td></tr><tr><td>H<sub>2</sub> + N<sub>2</sub></td><td>0.196</td><td>0.242</td><td>0.196</td><td>0.224</td></tr></table>

<table><caption>**Table 2.** Maximum velocity *U*<sub>m</sub> (m/day) for a stable downward displacement of CH<sub>4</sub> or N<sub>2</sub> by H<sub>2</sub> calculated within the low-density (IG+Wilke) approximation. Reservoir permeability: *k* = 1 Darcy. These values should be multiplied by the number of Darcy of the reservoir and, for tilted fronts, by cos(α).</caption><tr><td>Temperature (K)</td><td>298</td><td>298</td><td>373</td><td>373</td></tr><tr><td>Pressure (MPa)</td><td>5</td><td>25</td><td>5</td><td>25</td></tr><tr><td>H<sub>2</sub> + CH<sub>4</sub></td><td>1.45</td><td>3.75</td><td>1.94</td><td>2.97</td></tr><tr><td>H<sub>2</sub> + N<sub>2</sub></td><td>0.88</td><td>3.07</td><td>1.98</td><td>2.26</td></tr></table>

substituting the density, molar weight, and viscosity of the WG at storage conditions for those of H<sub>2</sub> in the above two equations. If the gases are not pure components but a mixture of components, then the densities, molar weights, and viscosities of the corresponding mixtures must be used in these equations, again taking care to use Eqns 10 and 11 only when appropriate.

Tables 1 and 2 provide values of *U*<sub>m</sub>/*U*<sub>im</sub> and *U*<sub>m</sub> for H<sub>2</sub> downwardly displacing CH<sub>4</sub> or N<sub>2</sub> as obtained within the low-density (IG+Wilke) approximation.

For downward displacements of CH<sub>4</sub> or N<sub>2</sub> by H<sub>2</sub> (corresponding to H<sub>2</sub> injection), the velocity limits *U*<sub>m</sub> lie in the range of a few meters per day. These values are considerably higher than practical velocities of a few centimeters/day, including when the angle α is high and approaches 90° or the reservoir has permeability in the tens of mDarcy. As a consequence, no fingering is to be expected, except where the mixing front tends to be locally vertical, that is, α ≈ 90° and, indeed, where α > 90°: then the front is unconditionally unstable.

As a practical example, injecting about 100,000 Nm<sup>3</sup> (or 0.82 tons) of H<sub>2</sub> per day (which is equivalent to about 2000 m<sup>3</sup> per day at reservoir pressure *P* ≈ 5 MPa) in an aquifer spanning an area of 1 km<sup>2</sup> with 10 wells corresponds to a Darcy velocity *U* ≈ 2 cm/day, which is about two orders of magnitude below the velocity limit *U*<sub>m</sub> of incipient fingering (see Table 2).

The above reasoning holds for mixing fronts positioned at some distance from the well. When injecting H<sub>2</sub> for the first time in a gas-bearing reservoir, the mixing front is initially located in the near wellbore region, where fluid velocities are high and, therefore, an unstable mixing front is expected. Viscous fingering could thus be an issue when injecting H<sub>2</sub> for the first time in a gas-bearing reservoir. This difficulty could be mitigated in part by controlling the injection rate in the early stages of the first injection sequence to obey the criterion for fluid velocities determined above and/or starting the injection process with CG-rich (or H<sub>2</sub>-poor) mixtures and gradually increasing the H<sub>2</sub> content in the injected fluid as the front moves away from the well and, when it is far enough, switch to 100% H<sub>2</sub> (y = 1). Grading the viscosity of a liquid displacing a more viscous liquid has been shown to mitigate viscous fingering,<sup>39</sup> while transverse dispersion attenuates finger growth provided the corresponding Peclet number, that is, fluid velocity, is small enough.<sup>40,41</sup> In subsequent injection sequences, the mixing front should be located far enough from the well – at a distance that, again, can be estimated by using the results of this section – which is equivalent to considering that the bottom part of the H<sub>2</sub> column serves as a CG for the working H<sub>2</sub>.<sup>4</sup>

More accurate models than the above IG and Wilke models are available to calculate *dρ/dμ* and its minimum value under the H<sub>2</sub>/CG mixtures and the *T* and *P* conditions of interest. Densities are better predicted by the GERG-2008 equation<sup>42</sup>, and more appropriate viscosity models are presented and discussed in the Appendix. Using those models, we find velocity limits for H<sub>2</sub>/CH<sub>4</sub> and H<sub>2</sub>/N<sub>2</sub> mixtures that are very close to those obtained with the IG+Wilke models and listed in Tables 1 and 2 (see the Appendix).

However, the IG+Wilke models fail to predict the densities and viscosities of H<sub>2</sub>/CO<sub>2</sub> mixtures under most of the *T* and *P* conditions of interest, especially for CO<sub>2</sub>-rich mixtures (y small). Experimentally, adding only 2 mol% H<sub>2</sub> to pure CO<sub>2</sub> is shown to lower the density by as much as 25%,<sup>43</sup> a feature that is much better accounted for by the GERG-2008 equation.<sup>42</sup> Likewise, Eqns 7 and 8 do not account for the viscosities of CO<sub>2</sub>-rich H<sub>2</sub>+CO<sub>2</sub> mixtures, as it is easily checked in comparison with recent experimental data<sup>44</sup> (see the Appendix). The models presented in the Appendix, such as the SuperTRAPP model (or the extended corresponding state [ECS] model), the Dean-Stiel-Brokaw (DSB) model (Simulis
<table><caption>**Table 1.** Values of *U*<sub>m</sub>/*U*<sub>im</sub> calculated within the low-density (IG+Wilke) approximation for H<sub>2</sub> downwardly displacing CH<sub>4</sub> or N<sub>2</sub>.</caption><tr><td>Temperature (K)</td><td>298</td><td>298</td><td>373</td><td>373</td></tr><tr><td>Pressure (MPa)</td><td>5</td><td>25</td><td>5</td><td>25</td></tr><tr><td>H<sub>2</sub> + CH<sub>4</sub></td><td>0.173</td><td>0.340</td><td>0.175</td><td>0.268</td></tr><tr><td>H<sub>2</sub> + N<sub>2</sub></td><td>0.196</td><td>0.242</td><td>0.196</td><td>0.224</td></tr></table>

<table><caption>**Table 2.** Maximum velocity *U*<sub>m</sub> (m/day) for a stable downward displacement of CH<sub>4</sub> or N<sub>2</sub> by H<sub>2</sub> calculated within the low-density (IG+Wilke) approximation. Reservoir permeability: *k* = 1 Darcy. These values should be multiplied by the number of Darcy of the reservoir and, for tilted fronts, by cos(α).</caption><tr><td>Temperature (K)</td><td>298</td><td>298</td><td>373</td><td>373</td></tr><tr><td>Pressure (MPa)</td><td>5</td><td>25</td><td>5</td><td>25</td></tr><tr><td>H<sub>2</sub> + CH<sub>4</sub></td><td>1.45</td><td>3.75</td><td>1.94</td><td>2.97</td></tr><tr><td>H<sub>2</sub> + N<sub>2</sub></td><td>0.88</td><td>3.07</td><td>1.98</td><td>2.26</td></tr></table>

substituting the density, molar weight, and viscosity of the WG at storage conditions for those of H<sub>2</sub> in the above two equations. If the gases are not pure components but a mixture of components, then the densities, molar weights, and viscosities of the corresponding mixtures must be used in these equations, again taking care to use Eqns 10 and 11 only when appropriate.

Tables 1 and 2 provide values of *U*<sub>m</sub>/*U*<sub>im</sub> and *U*<sub>m</sub> for H<sub>2</sub> downwardly displacing CH<sub>4</sub> or N<sub>2</sub> as obtained within the low-density (IG+Wilke) approximation.

For downward displacements of CH<sub>4</sub> or N<sub>2</sub> by H<sub>2</sub> (corresponding to H<sub>2</sub> injection), the velocity limits *U*<sub>m</sub> lie in the range of a few meters per day. These values are considerably higher than practical velocities of a few centimeters/day, including when the angle α is high and approaches 90° or the reservoir has permeability in the tens of mDarcy. As a consequence, no fingering is to be expected, except where the mixing front tends to be locally vertical, that is, α ≈ 90° and, indeed, where α > 90°: then the front is unconditionally unstable.

As a practical example, injecting about 100,000 Nm<sup>3</sup> (or 0.82 tons) of H<sub>2</sub> per day (which is equivalent to about 2000 m<sup>3</sup> per day at reservoir pressure *P* ≈ 5 MPa) in an aquifer spanning an area of 1 km<sup>2</sup> with 10 wells corresponds to a Darcy velocity *U* ≈ 2 cm/day, which is about two orders of magnitude below the velocity limit *U*<sub>m</sub> of incipient fingering (see Table 2).

The above reasoning holds for mixing fronts positioned at some distance from the well. When injecting H<sub>2</sub> for the first time in a gas-bearing reservoir, the mixing front is initially located in the near wellbore region, where fluid velocities are high and, therefore, an unstable mixing front is expected. Viscous fingering could thus be an issue when injecting H<sub>2</sub> for the first time in a gas-bearing reservoir. This difficulty could be mitigated in part by controlling the injection rate in the early stages of the first injection sequence to obey the criterion for fluid velocities determined above and/or starting the injection process with CG-rich (or H<sub>2</sub>-poor) mixtures and gradually increasing the H<sub>2</sub> content in the injected fluid as the front moves away from the well and, when it is far enough, switch to 100% H<sub>2</sub> (y = 1). Grading the viscosity of a liquid displacing a more viscous liquid has been shown to mitigate viscous fingering,<sup>39</sup> while transverse dispersion attenuates finger growth provided the corresponding Peclet number, that is, fluid velocity, is small enough.<sup>40,41</sup> In subsequent injection sequences, the mixing front should be located far enough from the well – at a distance that, again, can be estimated by using the results of this section – which is equivalent to considering that the bottom part of the H<sub>2</sub> column serves as a CG for the working H<sub>2</sub>.<sup>4</sup>

More accurate models than the above IG and Wilke models are available to calculate *dρ/dμ* and its minimum value under the H<sub>2</sub>/CG mixtures and the *T* and *P* conditions of interest. Densities are better predicted by the GERG-2008 equation<sup>42</sup>, and more appropriate viscosity models are presented and discussed in the Appendix. Using those models, we find velocity limits for H<sub>2</sub>/CH<sub>4</sub> and H<sub>2</sub>/N<sub>2</sub> mixtures that are very close to those obtained with the IG+Wilke models and listed in Tables 1 and 2 (see the Appendix).

However, the IG+Wilke models fail to predict the densities and viscosities of H<sub>2</sub>/CO<sub>2</sub> mixtures under most of the *T* and *P* conditions of interest, especially for CO<sub>2</sub>-rich mixtures (y small). Experimentally, adding only 2 mol% H<sub>2</sub> to pure CO<sub>2</sub> is shown to lower the density by as much as 25%,<sup>43</sup> a feature that is much better accounted for by the GERG-2008 equation.<sup>42</sup> Likewise, Eqns 7 and 8 do not account for the viscosities of CO<sub>2</sub>-rich H<sub>2</sub>+CO<sub>2</sub> mixtures, as it is easily checked in comparison with recent experimental data<sup>44</sup> (see the Appendix). The models presented in the Appendix, such as the SuperTRAPP model (or the extended corresponding state [ECS] model), the Dean-Stiel-Brokaw (DSB) model (Simulis
<figure><img src="image_4.png" /><figcaption>Figure 2. Values of $kgd\rho/d\mu$ as a function of the H₂ molar content y of the CO₂ + H₂ mixture for various temperature and pressure conditions. Densities are those of the GERG-2008 model,⁴² and viscosities are determined using the SuperTRAPP model (see the Appendix). The inset shows the region of high H₂ molar content (y → 1), where $kgd\rho/d\mu$ reaches its global minimum, for all T and P conditions but T = 298 K and P = 25 MPa.</figcaption></figure>

Thermodynamics) or Wilke's model with a HZ exchange coefficient (Eqn 9), provide viscosity values in much better agreement with experimental values.

Figure 2 displays the values of $kgd\rho/d\mu$ as a function of the H₂ molar fraction y in H₂/CO₂ mixtures at four different sets of temperatures and pressures, as obtained using the GERG-2008 model for densities and the SuperTRAPP model for viscosities. $kgd\rho/d\mu$ decreases monotonically with increasing y over the interval [0,1] only for the lowest pressure (5 MPa) and highest temperature (373 K). For the other T and P conditions, $d\rho/d\mu$ exhibits a maximum in the interval [0,1] and two local minima, for y = 0 and y = 1. For the lowest temperature (T = 298 K) and highest pressure (P = 25 MPa) investigated, the global minimum of $kgd\rho/d\mu$ is reached for y = 0, meaning that, when fluid velocity reaches or slightly exceeds the limit $U_m$, the fingering process is initiated at the bottom (i.e., CO₂-rich) tail of the mixing front, similarly to what happens at miscible liquid/liquid fronts.³⁴ In all other cases, the instability is initiated in the upper

<table><caption>Table 3. Maximum velocities $U_m$ (m/day) for a stable downward displacement of CO₂ by H₂ within the low-density approximation (IG+Wilke) and the GERG-2008 model for densities and SuperTRAPP model for viscosities (see the Appendix).</caption><tr><td>Temperature (K)</td><td>298</td><td>298</td><td>373</td><td>373</td></tr><tr><td>Pressure (MPa)</td><td>5</td><td>25</td><td>5</td><td>25</td></tr><tr><td>IG+Wilke</td><td>1.12</td><td>3.3</td><td>0.77</td><td>2.7</td></tr><tr><td>GERG-2008+SuperTRAPP</td><td>1.04</td><td>3.71</td><td>0.76</td><td>3.61</td></tr></table>

region of the mixing zone, where H₂ is the fluid saturating the porous medium.

Table 3 lists the velocity limits $U_m$ for a stable downward displacement of CO₂ by H₂ as calculated within the IG+Wilke models and the GERG-2008+SuperTRAPP models for a reservoir with permeability k = 1 Darcy. As discussed above, the incipient instability occurs on the H₂-rich side of the front (y = 1) in all cases, except with the
<figure><img src="image_5.png" /><figcaption>Figure 2. Values of $kgd\rho/d\mu$ as a function of the H₂ molar content y of the CO₂ + H₂ mixture for various temperature and pressure conditions. Densities are those of the GERG-2008 model,⁴² and viscosities are determined using the SuperTRAPP model (see the Appendix). The inset shows the region of high H₂ molar content (y → 1), where $kgd\rho/d\mu$ reaches its global minimum, for all T and P conditions but T = 298 K and P = 25 MPa.</figcaption></figure>

Thermodynamics) or Wilke's model with a HZ exchange coefficient (Eqn 9), provide viscosity values in much better agreement with experimental values.

Figure 2 displays the values of $kgd\rho/d\mu$ as a function of the H₂ molar fraction y in H₂/CO₂ mixtures at four different sets of temperatures and pressures, as obtained using the GERG-2008 model for densities and the SuperTRAPP model for viscosities. $kgd\rho/d\mu$ decreases monotonically with increasing y over the interval [0,1] only for the lowest pressure (5 MPa) and highest temperature (373 K). For the other T and P conditions, $d\rho/d\mu$ exhibits a maximum in the interval [0,1] and two local minima, for y = 0 and y = 1. For the lowest temperature (T = 298 K) and highest pressure (P = 25 MPa) investigated, the global minimum of $kgd\rho/d\mu$ is reached for y = 0, meaning that, when fluid velocity reaches or slightly exceeds the limit $U_m$, the fingering process is initiated at the bottom (i.e., CO₂-rich) tail of the mixing front, similarly to what happens at miscible liquid/liquid fronts.³⁴ In all other cases, the instability is initiated in the upper

<table><caption>Table 3. Maximum velocities $U_m$ (m/day) for a stable downward displacement of CO₂ by H₂ within the low-density approximation (IG+Wilke) and the GERG-2008 model for densities and SuperTRAPP model for viscosities (see the Appendix).</caption><tr><td>Temperature (K)</td><td>298</td><td>298</td><td>373</td><td>373</td></tr><tr><td>Pressure (MPa)</td><td>5</td><td>25</td><td>5</td><td>25</td></tr><tr><td>IG+Wilke</td><td>1.12</td><td>3.3</td><td>0.77</td><td>2.7</td></tr><tr><td>GERG-2008+SuperTRAPP</td><td>1.04</td><td>3.71</td><td>0.76</td><td>3.61</td></tr></table>

region of the mixing zone, where H₂ is the fluid saturating the porous medium.

Table 3 lists the velocity limits $U_m$ for a stable downward displacement of CO₂ by H₂ as calculated within the IG+Wilke models and the GERG-2008+SuperTRAPP models for a reservoir with permeability k = 1 Darcy. As discussed above, the incipient instability occurs on the H₂-rich side of the front (y = 1) in all cases, except with the
GERG-2008+SuperTRAPP models at 298 K and 25 MPa where it occurs at the CO₂-rich, that is, liquid-like, side of the front (y = 0): the agreement between the two values of Uc obtained within the two models is fortuitous.

A description of the size of those fingers and how they develop once they have been initiated is beyond the scope of this work. This development strongly depends on the density and viscosity profiles along the mixing front³⁵ as well as on the dispersive properties of the porous medium.²² To the best of our knowledge, the few experimental data available on how fingers develop when the instability sets in have been acquired with miscible liquid/liquid couples,³⁶ not with miscible gas/gas or gas/liquid systems.

A general comment is in order here. Density varies with viscosity in very different manners in contrasted gas/gas and liquid/liquid mixtures, and this has consequences for the velocity range that ensures front stability and mixing limited to dispersion. As a rule, two liquids, one viscous and dense (e.g., an oil) and the other lighter and less viscous (e.g., a solvent) exhibit a limited density difference (typically a few tens of percent difference), yet their viscosities may differ by orders of magnitude. The densities of H₂ and the CGs considered here differ by one order of magnitude or more, but their viscosity ratio rarely exceeds 2.5 (see the above: the maximum mobility ratio M<sub>r</sub> is 2.5, except when the CG is dense CO₂). When composition varies in H₂/CG mixtures, the values of the derivatives *dρ/dμ* – and their minima when composition varies – are therefore expected to be much larger than those of mixtures of two miscible liquids (e.g., oil/solvent mixtures). Compared to miscible liquid/liquid displacement fronts, H₂/CG fronts are, therefore, stable and piston like over a much larger range of H₂ velocities or injection rates and/or over a much larger range of reservoir permeabilities and inclination angles α. This is one of the most important results of this study, and a strong argument in favor of storing H₂ on top of a denser, more viscous gas acting as a CG. Dispersive mixing (see the previous subsection) is expected to be the only cause of the spreading of the front, except where the H₂-CG front is close to being vertical (i.e., α ≈ 90°).

In practice, the permeability field is not homogeneous and the question arises as to the effects of permeability heterogeneities on the miscible displacements of interest here. This subject has been investigated for miscible liquid/liquid systems by

various numerical simulation methods, most often in a linear drive configuration: one fluid is injected into a reservoir layer through a line of (injection) wells that displaces the other fluid to another line of (production) wells on the other side of the layer.⁴⁵⁻⁴⁷ The conclusions of those studies are that, for small enough variances and correlation lengths of the permeability field, the onset of fingering and finger growth is not affected by permeability heterogeneities, whereas finger growth interacts strongly with the permeability field when this field is strongly heterogeneous.⁴⁵⁻⁴⁷ The focus in this section is on how viscous and gravity forces compete to generate either a stable mixing front or fingering, for, respectively, low and high ratios of the gravity to viscous forces. The ratio of the viscous to gravity forces is one of the dimensionless numbers that control the flow patterns of H₂ displacing a CG.¹⁴ We suggest to define this ratio, or at least its fluid-dependent part, by using the minimum value of *dρ/dμ* over all possible mixtures of the two fluids.

Lastly, the hypothesis has been made of first-contact miscibility, that is, complete miscibility of H₂ and CG, irrespective of their proportion. It is worth pointing out here that first-contact miscibility is not achieved in CO₂-rich H₂ + CO₂ mixtures in the ranges 293–303 K and 4–7 MPa.⁴⁸ These mixtures and conditions are not addressed in this work.

## Caprock integrity and storage capacities

In this section, we first recall how buoyancy and the interfacial properties of the seal rock/brine/gas system are involved in the gas-trapping process beneath the seal rock – usually referred to as structural or stratigraphic trapping. A brief overview is then given of the experimental efforts undertaken to characterize those interfacial properties. Then, we examine the storage scheme of Fig. 1, with H₂ placed on top of a CG with limited mixing between the two gases (see the previous section), and draw consequences as to its safety and capacity. We also discuss why and how the placement of H₂ on top of an existing gas storage or a gas reservoir (e.g., a depleted hydrocarbon reservoir) reinforces the capillary-sealing efficiency of the seal rock and renders the storage of *both* gases more secure while allowing some amount of H₂ to be stored together with the gas already present.

Even though it is not further discussed and explored in this work, it is worth pointing out that interfacial or
GERG-2008+SuperTRAPP models at 298 K and 25 MPa where it occurs at the CO₂-rich, that is, liquid-like, side of the front (y = 0): the agreement between the two values of U_c obtained within the two models is fortuitous.

A description of the size of those fingers and how they develop once they have been initiated is beyond the scope of this work. This development strongly depends on the density and viscosity profiles along the mixing front³⁵ as well as on the dispersive properties of the porous medium.²² To the best of our knowledge, the few experimental data available on how fingers develop when the instability sets in have been acquired with miscible liquid/liquid couples,³⁶ not with miscible gas/gas or gas/liquid systems.

A general comment is in order here. Density varies with viscosity in very different manners in contrasted gas/gas and liquid/liquid mixtures, and this has consequences for the velocity range that ensures front stability and mixing limited to dispersion. As a rule, two liquids, one viscous and dense (e.g., an oil) and the other lighter and less viscous (e.g., a solvent) exhibit a limited density difference (typically a few tens of percent difference), yet their viscosities may differ by orders of magnitude. The densities of H₂ and the CGs considered here differ by one order of magnitude or more, but their viscosity ratio rarely exceeds 2.5 (see the above: the maximum mobility ratio *M*<sub>r</sub> is 2.5, except when the CG is dense CO₂). When composition varies in H₂/CG mixtures, the values of the derivatives *dρ/dμ* – and their minima when composition varies – are therefore expected to be much larger than those of mixtures of two miscible liquids (e.g., oil/solvent mixtures). Compared to miscible liquid/liquid displacement fronts, H₂/CG fronts are, therefore, stable and piston like over a much larger range of H₂ velocities or injection rates and/or over a much larger range of reservoir permeabilities and inclination angles α. This is one of the most important results of this study, and a strong argument in favor of storing H₂ on top of a denser, more viscous gas acting as a CG. Dispersive mixing (see the previous subsection) is expected to be the only cause of the spreading of the front, except where the H₂-CG front is close to being vertical (i.e., α ≈ 90°).

In practice, the permeability field is not homogeneous and the question arises as to the effects of permeability heterogeneities on the miscible displacements of interest here. This subject has been investigated for miscible liquid/liquid systems by

various numerical simulation methods, most often in a linear drive configuration: one fluid is injected into a reservoir layer through a line of (injection) wells that displaces the other fluid to another line of (production) wells on the other side of the layer.⁴⁵⁻⁴⁷ The conclusions of those studies are that, for small enough variances and correlation lengths of the permeability field, the onset of fingering and finger growth is not affected by permeability heterogeneities, whereas finger growth interacts strongly with the permeability field when this field is strongly heterogeneous.⁴⁵⁻⁴⁷ The focus in this section is on how viscous and gravity forces compete to generate either a stable mixing front or fingering, for, respectively, low and high ratios of the gravity to viscous forces. The ratio of the viscous to gravity forces is one of the dimensionless numbers that control the flow patterns of H₂ displacing a CG.¹⁴ We suggest to define this ratio, or at least its fluid-dependent part, by using the minimum value of *dρ/dμ* over all possible mixtures of the two fluids.

Lastly, the hypothesis has been made of first-contact miscibility, that is, complete miscibility of H₂ and CG, irrespective of their proportion. It is worth pointing out here that first-contact miscibility is not achieved in CO₂-rich H₂ + CO₂ mixtures in the ranges 293–303 K and 4–7 MPa.⁴⁸ These mixtures and conditions are not addressed in this work.

## Caprock integrity and storage capacities

In this section, we first recall how buoyancy and the interfacial properties of the seal rock/brine/gas system are involved in the gas-trapping process beneath the seal rock – usually referred to as structural or stratigraphic trapping. A brief overview is then given of the experimental efforts undertaken to characterize those interfacial properties. Then, we examine the storage scheme of Fig. 1, with H₂ placed on top of a CG with limited mixing between the two gases (see the previous section), and draw consequences as to its safety and capacity. We also discuss why and how the placement of H₂ on top of an existing gas storage or a gas reservoir (e.g., a depleted hydrocarbon reservoir) reinforces the capillary-sealing efficiency of the seal rock and renders the storage of *both* gases more secure while allowing some amount of H₂ to be stored together with the gas already present.

Even though it is not further discussed and explored in this work, it is worth pointing out that interfacial or
capillary forces or, more precisely, the ratios of the capillary to viscous and to buoyancy forces, that is, the so-called capillary and Bond numbers, have a strong impact on flow patterns in porous media: as an example, in the presence of water capillarity prevents the (nonwetting) gas to enter the small permeabilities zones at low enough flow rates (see for details in Wang et al.¹⁴).

## Capillary trapping of a single gas phase beneath a seal rock

The gas phase beneath the seal rock is overpressured with respect to the brine, which forms an irreducible phase in the reservoir and fully saturates both the aquifer and the seal rock. The buoyancy pressure is related to the gas column height h above the contact with the aquifer (water-gas contact [WGC]) and to the difference between phase densities, ρ<sub>w</sub> and ρ<sub>g</sub> (assumed below to be constant), as follows:

$$P_b = (\rho_w - \rho_g)gh. \qquad (12)$$

The assumption of constant gas density is valid as long as the column height does not exceed a few tens of meters – the storage scenario considered in practical applications and in this work – except if the T-P path along this column crosses the liquid/vapor phase envelope or approaches the critical point where compressibility is high: this excludes CO₂ in some cold reservoirs and very near the critical point of CO₂, P ≈ 7.38 MPa and T ≈ 304 K.

The trapping of the gas phase beneath the seal rock is ensured as long as the buoyancy pressure does not exceed the capillary entry (or displacement) pressure of the gas phase into the seal rock, equal from the Laplace equation to

$$P_{ce,g} = P_g - P_w = \frac{2\gamma_{w,g}\cos(\theta)}{r}, \qquad (13)$$

where r is a representative of the largest pore thresholds in the seal rock, γ<sub>w,g</sub> is the IFT between the brine and gas phases and θ is the contact angle at the triple (mineral/brine/gas) line in a drainage process, that is, gas displacing brine from the substrate. Stated otherwise, the storage capacity or maximum gas column height is determined by the competition between two adverse effects: buoyancy that drives the gas upwards, and capillarity that prevents the gas from breaking through the brine-saturated seal rock. The

maximum column height, noted H, is obtained by equating the two pressures Eqns 12 and 13<sup>49,50</sup>

$$H = \frac{P_{ce,g}}{g(\rho_w - \rho_g)} = \frac{2\gamma_{w,g}\cos(\theta)}{r(\rho_w - \rho_g)g}. \qquad (14)$$

The capillary entry (or displacement) pressure is determined experimentally or estimated using measurements conducted with a convenient gas, such as N₂, or using mercury intrusion measurements, and then using the equation P<sub>ce,H₂</sub> = P<sub>ce</sub>γ<sub>w,H₂</sub>cos(θ)/γ'cos(θ'), where P<sub>ce</sub> is the entry pressure of N₂ or liquid mercury, γ' and θ' are the IFT between brine and N₂ or between liquid mercury and its vapor, and the contact angle is that in the wetting phase (water or mercury vapor). The contact angle has a typical value of ≈ 40° for mercury.⁵¹

Clearly, the storage capacity or maximum storage height increases with brine/gas IFT γ<sub>w,g</sub> and cos(θ), that is, with the water-wet character of the porous rock, and it decreases with increasing buoyancy, that is, with increasing density difference between the brine and gas phases.

## Relevant interfacial properties

A considerable amount of experimental data, including an increasing number of experiments in silico (mostly molecular dynamics simulations), is available on the interfacial properties that control the capillary-sealing efficiency of seal rocks with respect to the gases of interest: the water (or brine)/gas IFT γ<sub>w,g</sub> and the contact angle θ<sub>w,g</sub> at the triple (mineral/brine/gas) line. Only a brief review with the relevant information is presented here, with a focus on how these two properties are modified when the hydrogen is contaminated with the CG, whether CO₂, CH₄, or N₂.

The IFT between two immiscible phases is related to the affinities between the molecules making up these phases. IFTs γ<sub>w,g</sub> between water or brine and a gas phase are correlated to the water/gas mutual solubility:⁵² the more water soluble the gas, the larger the IFT decrease when pressure increases or, equivalently (from Gibbs adsorption equation), the higher the adsorption of the gas onto water (or brine). For given T, P, and brine salinity, the following ranking

$$\gamma_{w,H_2} > \gamma_{w,N_2} > \gamma_{w,CH_4} > \gamma_{w,CO_2} \qquad (15)$$

is indeed observed. Massoudi and King⁵³ were the first and, to the best of our knowledge, the only researchers to have measured in the same experimental campaign
capillary forces or, more precisely, the ratios of the capillary to viscous and to buoyancy forces, that is, the so-called capillary and Bond numbers, have a strong impact on flow patterns in porous media: as an example, in the presence of water capillarity prevents the (nonwetting) gas to enter the small permeabilities zones at low enough flow rates (see for details in Wang et al.¹⁴).

## Capillary trapping of a single gas phase beneath a seal rock

The gas phase beneath the seal rock is overpressured with respect to the brine, which forms an irreducible phase in the reservoir and fully saturates both the aquifer and the seal rock. The buoyancy pressure is related to the gas column height h above the contact with the aquifer (water-gas contact [WGC]) and to the difference between phase densities, $ρ_w$ and $ρ_g$ (assumed below to be constant), as follows:

$$P_b = (\rho_w - \rho_g)gh. \qquad (12)$$

The assumption of constant gas density is valid as long as the column height does not exceed a few tens of meters – the storage scenario considered in practical applications and in this work – except if the T-P path along this column crosses the liquid/vapor phase envelope or approaches the critical point where compressibility is high: this excludes CO₂ in some cold reservoirs and very near the critical point of CO₂, $P \approx 7.38$ MPa and $T \approx 304$ K.

The trapping of the gas phase beneath the seal rock is ensured as long as the buoyancy pressure does not exceed the capillary entry (or displacement) pressure of the gas phase into the seal rock, equal from the Laplace equation to

$$P_{ce,g} = P_g - P_w = \frac{2\gamma_{w,g}\cos(\theta)}{r}, \qquad (13)$$

where r is a representative of the largest pore thresholds in the seal rock, $γ_{w,g}$ is the IFT between the brine and gas phases and $θ$ is the contact angle at the triple (mineral/brine/gas) line in a drainage process, that is, gas displacing brine from the substrate. Stated otherwise, the storage capacity or maximum gas column height is determined by the competition between two adverse effects: buoyancy that drives the gas upwards, and capillarity that prevents the gas from breaking through the brine-saturated seal rock. The

maximum column height, noted H, is obtained by equating the two pressures Eqns 12 and 13⁴⁹,⁵⁰

$$H = \frac{P_{ce,g}}{g(\rho_w - \rho_g)} = \frac{2\gamma_{w,g}\cos(\theta)}{r(\rho_w - \rho_g)g}. \qquad (14)$$

The capillary entry (or displacement) pressure is determined experimentally or estimated using measurements conducted with a convenient gas, such as N₂, or using mercury intrusion measurements, and then using the equation $P_{ce,H_2} = P_{ce}\gamma_{w,H_2}\cos(\theta)/\gamma'\cos(\theta')$, where $P_{ce}$ is the entry pressure of N₂ or liquid mercury, $γ'$ and $θ'$ are the IFT between brine and N₂ or between liquid mercury and its vapor, and the contact angle is that in the wetting phase (water or mercury vapor). The contact angle has a typical value of ≈ 40° for mercury.⁵¹

Clearly, the storage capacity or maximum storage height increases with brine/gas IFT $γ_{w,g}$ and $\cos(\theta)$, that is, with the water-wet character of the porous rock, and it decreases with increasing buoyancy, that is, with increasing density difference between the brine and gas phases.

## Relevant interfacial properties

A considerable amount of experimental data, including an increasing number of experiments in silico (mostly molecular dynamics simulations), is available on the interfacial properties that control the capillary-sealing efficiency of seal rocks with respect to the gases of interest: the water (or brine)/gas IFT $γ_{w,g}$ and the contact angle $θ_{w,g}$ at the triple (mineral/brine/gas) line. Only a brief review with the relevant information is presented here, with a focus on how these two properties are modified when the hydrogen is contaminated with the CG, whether CO₂, CH₄, or N₂.

The IFT between two immiscible phases is related to the affinities between the molecules making up these phases. IFTs $γ_{w,g}$ between water or brine and a gas phase are correlated to the water/gas mutual solubility:⁵² the more water soluble the gas, the larger the IFT decrease when pressure increases or, equivalently (from Gibbs adsorption equation), the higher the adsorption of the gas onto water (or brine). For given T, P, and brine salinity, the following ranking

$$\gamma_{w,H_2} > \gamma_{w,N_2} > \gamma_{w,CH_4} > \gamma_{w,CO_2} \qquad (15)$$

is indeed observed. Massoudi and King⁵³ were the first and, to the best of our knowledge, the only researchers to have measured in the same experimental campaign
water/gas IFTs for the above suite of gases under similar T and P. Subsequent studies have been concerned with only a subset of those gases and/or have gathered and examined IFTs from various sources.⁵⁴ Brine/H₂ IFTs lie in the range of ≈50–75 mN/m,⁵³⁻⁵⁸ whereas under similar T, P, and brine salinity, brine/CO₂ IFTs are about twice smaller, in the range of ≈20–35 mN/m.⁴⁹,⁵⁹,⁶⁰ Brine/CH₄ and brine/N₂ IFTs have values intermediate between brine/CO₂ IFTs and brine/H₂ IFTs, the latter being slightly smaller than the former.⁶¹,⁶² Correlations and simple models are available for brine–H₂ IFTs⁵⁴,⁵⁵ and brine–CO₂ IFTs⁶³,⁶⁴ as a function of T, P, and brine salinity (see also the recent critical review by Mouallem et al.⁶⁵). When H₂ is contaminated by the CG, there is indeed a decrease of the IFT, which is more pronounced with CO₂ than with the two other gases.⁵⁵,⁶⁶⁻⁶⁸ As a rule, brine/gas IFTs vary monotonously (albeit nonlinearly) with gas composition.⁶⁹

The information on contact angles and wetting characteristics of mineral or rock/brine/gas systems is much more abundant than that on brine/gas IFTs. In addition to T, P, and brine salinity, an important role is played by the mineral or rock, whose tension with brine or gas is not measured but can be evaluated indirectly.⁷⁰⁻⁷³ As a rule, H₂ does not significantly alter the water-wet character (i.e., contact angles remain low) of most minerals, even at high P,⁵⁶,⁷³,⁷⁴ whereas pressurized CO₂ may induce a significant increase in contact angles – up to values close to 90° (neutral-wet) at high P,⁷⁵,⁷⁶ especially when some (minute amount of) organic acid is present.⁷⁷ CO₂ added to H₂ alters the water-wet character of minerals or rock substrates more than N₂ and CH₄. Even though a few data are still lacking and/or exhibit contradictory trends, the cosines of the contact angles appear to rank similarly to brine/gas IFTs (see Eqn 15).

Brine–gas IFTs and contact angles can also be inferred from the capillary pressures obtained from two-phase (water and H₂) core floods or mercury injection experiments and are found to be consistent with the directly measured IFTs and angles.³³

An immediate consequence of those experimental data (see Eqn 13) is the increase of the capillary entry pressure and, therefore, the enhancement of storage safety when H₂ is injected on top of an existing CO₂, CH₄, or N₂ storage or, more generally, when a gas having larger brine/gas IFT and not altering the water-wet character of the seal rock is added in a reservoir on top of another gas – e.g., when adding

methane on top of CO₂. If the hydrogen in contact with the seal rock is not pure but contaminated by the CG, then the IFT is still higher and the contact angle is lower than the IFT and contact angles that prevailed before hydrogen injection. This point is further discussed below (see the subsection ‘Enhancement of storage safety by adding H₂ at the top of an existing gas column’).

## H₂ storage capacities with various CGs and CG column heights

The base case, against which configurations with CGs other than H₂ are compared, is that of an entire gas column consisting of H₂, part of it serving as a CG that is indeed not counted in the storage capacity. The reservoir is permeable enough (k larger than a few tens of mdarcy) to disregard the difference between the WGC (the maximum depth where free gas is present) and the free water level, defined as the depth where the pressure in the water phase is equal to the (extrapolated) pressure in the other, nonwetting phase (i.e., the gas phase). Given these assumptions, the maximum height of H₂ that can be stored safely in the aquifer is given by Eqn 14, or

$$H = \frac{P_{ce,H_2}}{g(\rho_w - \rho_{H_2})} \quad (16)$$

When as depicted in Fig. 1 a CG other than H₂ is used, the height of H₂ can be calculated using the following equation:

$$H = \frac{P_{ce,H_2} - (\rho_w - \rho_{CG})gh_{CG}}{(\rho_w - \rho_{H_2})g} \quad (17)$$

where h<sub>CG</sub> refers to the CG column height. Equation 17 is similar to Eqn 16 in which the capillary entry pressure P<sub>ce,H₂</sub> is replaced by an effective capillary entry pressure equal to P<sub>ce,H₂</sub> minus the buoyancy pressure exerted by the CG: (ρ<sub>w</sub> − ρ<sub>CG</sub>)gh<sub>CG</sub>. If the CG is made up of H₂, then Eqn 17 is equivalent to Eqn 16 with H replaced by H + h<sub>CG</sub>, the total height of the H₂ column.

It is clear from the above equation that, for larger CG densities, the buoyancy pressure exerted by the CG, (ρ<sub>w</sub> − ρ<sub>CG</sub>)gh<sub>CG</sub>, is smaller, and therefore the effective capillary entry pressure and H₂ storage height are larger: in other words, a given height h<sub>CG</sub> of a denser CG allows a larger quantity (or height H) of H₂ to be stored safely. This approach can also be easily extended to any WG/CG couple, including when the WG and CG are mixtures.
water/gas IFTs for the above suite of gases under similar T and P. Subsequent studies have been concerned with only a subset of those gases and/or have gathered and examined IFTs from various sources.⁵⁴ Brine/H₂ IFTs lie in the range of ≈50–75 mN/m,⁵³⁻⁵⁸ whereas under similar T, P, and brine salinity, brine/CO₂ IFTs are about twice smaller, in the range of ≈20–35 mN/m.⁴⁹,⁵⁹,⁶⁰ Brine/CH₄ and brine/N₂ IFTs have values intermediate between brine/CO₂ IFTs and brine/H₂ IFTs, the latter being slightly smaller than the former.⁶¹,⁶² Correlations and simple models are available for brine–H₂ IFTs⁵⁴,⁵⁵ and brine–CO₂ IFTs⁶³,⁶⁴ as a function of T, P, and brine salinity (see also the recent critical review by Mouallem et al.⁶⁵). When H₂ is contaminated by the CG, there is indeed a decrease of the IFT, which is more pronounced with CO₂ than with the two other gases.⁵⁵,⁶⁶⁻⁶⁸ As a rule, brine/gas IFTs vary monotonously (albeit nonlinearly) with gas composition.⁶⁹

The information on contact angles and wetting characteristics of mineral or rock/brine/gas systems is much more abundant than that on brine/gas IFTs. In addition to T, P, and brine salinity, an important role is played by the mineral or rock, whose tension with brine or gas is not measured but can be evaluated indirectly.⁷⁰⁻⁷³ As a rule, H₂ does not significantly alter the water-wet character (i.e., contact angles remain low) of most minerals, even at high P,⁵⁶,⁷³,⁷⁴ whereas pressurized CO₂ may induce a significant increase in contact angles – up to values close to 90° (neutral-wet) at high P,⁷⁵,⁷⁶ especially when some (minute amount of) organic acid is present.⁷⁷ CO₂ added to H₂ alters the water-wet character of minerals or rock substrates more than N₂ and CH₄. Even though a few data are still lacking and/or exhibit contradictory trends, the cosines of the contact angles appear to rank similarly to brine/gas IFTs (see Eqn 15).

Brine–gas IFTs and contact angles can also be inferred from the capillary pressures obtained from two-phase (water and H₂) core floods or mercury injection experiments and are found to be consistent with the directly measured IFTs and angles.³³

An immediate consequence of those experimental data (see Eqn 13) is the increase of the capillary entry pressure and, therefore, the enhancement of storage safety when H₂ is injected on top of an existing CO₂, CH₄, or N₂ storage or, more generally, when a gas having larger brine/gas IFT and not altering the water-wet character of the seal rock is added in a reservoir on top of another gas – e.g., when adding

methane on top of CO₂. If the hydrogen in contact with the seal rock is not pure but contaminated by the CG, then the IFT is still higher and the contact angle is lower than the IFT and contact angles that prevailed before hydrogen injection. This point is further discussed below (see the subsection ‘Enhancement of storage safety by adding H₂ at the top of an existing gas column’).

## H₂ storage capacities with various CGs and CG column heights

The base case, against which configurations with CGs other than H₂ are compared, is that of an entire gas column consisting of H₂, part of it serving as a CG that is indeed not counted in the storage capacity. The reservoir is permeable enough (k larger than a few tens of mdarcy) to disregard the difference between the WGC (the maximum depth where free gas is present) and the free water level, defined as the depth where the pressure in the water phase is equal to the (extrapolated) pressure in the other, nonwetting phase (i.e., the gas phase). Given these assumptions, the maximum height of H₂ that can be stored safely in the aquifer is given by Eqn 14, or

$$H = \frac{P_{ce,H_2}}{g(\rho_w - \rho_{H_2})} \qquad (16)$$

When as depicted in Fig. 1 a CG other than H₂ is used, the height of H₂ can be calculated using the following equation:

$$H = \frac{P_{ce,H_2} - (\rho_w - \rho_{CG})gh_{CG}}{(\rho_w - \rho_{H_2})g} \qquad (17)$$

where h<sub>CG</sub> refers to the CG column height. Equation 17 is similar to Eqn 16 in which the capillary entry pressure P<sub>ce,H₂</sub> is replaced by an effective capillary entry pressure equal to P<sub>ce,H₂</sub> minus the buoyancy pressure exerted by the CG: (ρ<sub>w</sub> − ρ<sub>CG</sub>)gh<sub>CG</sub>. If the CG is made up of H₂, then Eqn 17 is equivalent to Eqn 16 with H replaced by H + h<sub>CG</sub>, the total height of the H₂ column.

It is clear from the above equation that, for larger CG densities, the buoyancy pressure exerted by the CG, (ρ<sub>w</sub> − ρ<sub>CG</sub>)gh<sub>CG</sub>, is smaller, and therefore the effective capillary entry pressure and H₂ storage height are larger: in other words, a given height h<sub>CG</sub> of a denser CG allows a larger quantity (or height H) of H₂ to be stored safely. This approach can also be easily extended to any WG/CG couple, including when the WG and CG are mixtures.
<figure><img src="image_5.png" /><figcaption>Figure 3. Maximum hydrogen storage heights H at various reservoir pressures from 5 to 30 MPa (or reservoir depths from 500 to 3000 m) and T = 298 K (blue: cold reservoir), 323 K (green: warm reservoir), and 372 K (red: hot reservoir), obtained from Eqns 17 and 18 assuming a capillary entry pressure of the seal rock P<sub>ce,H<sub>2</sub></sub> ≈ 1.5 MPa at 298 K and 5 MPa. Brine/H<sub>2</sub> IFTs are those measured by Chow et al.<sup>55</sup> Densities are obtained from the NIST (webbook.nist.gov). The CG height is h<sub>CG</sub> = 40 m. Densities and IFTs are considered constant over the gas column for each reservoir pressure (or depth) and temperature.</figcaption></figure>

As an example, Fig. 3 displays the maximum column heights H of H₂ that can be safely stored on top of a CG with height h<sub>CG</sub> = 40 m, beneath a seal rock with entry pressure P<sub>ce,H<sub>2</sub></sub> ≈ 1.5 MPa at T = 298 K and P = 5 MPa. This is a conservative value, 3–5 times smaller than the observed breakthrough pressures into carbonate-rich caprocks of CO₂ or nitrogen,⁷⁸ which have less favorable interfacial properties than H₂. The maximum storage heights (or storage capacities) are calculated for various reservoir pressures P from 5 to 30 MPa (or storage depths from 500 to 3000 m) and temperatures T = 298 K (cold reservoirs), 323 K (warm reservoirs), and 372 K (hot reservoirs) by means of Eqn 17, where the capillary entry pressure at T and P is inferred from that at T<sub>ref</sub> = 298 K and P<sub>ref</sub> = 5 MPa as follows:

$$P_{ce,H_2}(T, P) = P_{ce,H_2}(T_{ref}, P_{ref}) \frac{\gamma_{w,H_2}(T, P)}{\gamma_{w,H_2}(T_{ref}, P_{ref})}. \quad (18)$$

In the context of this work, the contact angle θ is assumed to be independent of T and P. If the dependence of θ with T and P is known, then the capillary entry pressure at T and P should be corrected as follows:

$$P_{ce,H_2}(T, P) = P_{ce,H_2}(T_{ref}, P_{ref}) \times \frac{\gamma_{w,H_2}(T, P) \cos[\theta(T, P)]}{\gamma_{w,H_2}(T_{ref}, P_{ref}) \cos[\theta(T_{ref}, P_{ref})]} \quad (19)$$

Figure 3 shows the impact of the density difference between the aqueous phase and CG on H₂ storage capacity. When H₂ is used as a CG, the density difference is the largest and therefore the H₂ storage capacity is the smallest. When a dense CG such as CO₂ is used, storage capacity is significantly higher, especially under low or moderate temperatures when CO₂ is a liquid or a supercritical fluid.

The same trends (data not shown) are observed when expressing the storage capacity in terms of energy per area unit (e.g., A = 1 km²) instead of using the H₂ storage capacity (H) directly:

$$\text{Energy per area } E = \rho_{H_2} eH \quad (20)$$

where e ≈ 120 kJ/g is the energy density of H₂.⁷⁹

## Enhancement of storage safety by adding H₂ at the top of an existing gas column

When a gas (such as H₂) is injected at the crest of a reservoir already filled with a gas (CO₂, CH₄, N₂, or their mixtures), the injected gas becomes the gas in contact with the seal rock. As a rule, H₂ has a larger capillary entry pressure into the seal rock than most (if not all) of the gases cited above, that is, P<sub>ce,H<sub>2</sub></sub> = λP<sub>ce,g</sub> with λ > 1. Due to the enhanced capillary entry pressure, an amount (or height H) of H₂ can be stored in addition to the gas present without any risk of leakage by capillary failure of the seal rock. The height of the H₂ column that can be stored safely on top of the gas column with height h<sub>g</sub> is obtained by combining Eqns 17 and 14:

$$H = \frac{(\lambda - 1)(\rho_w - \rho_g)}{(\rho_w - \rho_{H_2})} h_g \quad (21)$$

As an example, consider H₂ being injected on top of a column of supercritical CO₂ with density ρ<sub>CO₂</sub> ≈ 500 kg/m³. From the known IFTs and contact angles (see the subsection Relevant interfacial properties), the coefficient λ relating the two entry pressures, P<sub>ce,H₂</sub> and P<sub>ce,CO₂</sub>, is in the range of 2. The height of H₂ that can be
<figure><img src="image_6.png" /><figcaption>Figure 3. Maximum hydrogen storage heights H at various reservoir pressures from 5 to 30 MPa (or reservoir depths from 500 to 3000 m) and T = 298 K (blue: cold reservoir), 323 K (green: warm reservoir), and 372 K (red: hot reservoir), obtained from Eqns 17 and 18 assuming a capillary entry pressure of the seal rock $P_{ce,H_2} \approx 1.5$ MPa at 298 K and 5 MPa. Brine/H$_2$ IFTs are those measured by Chow et al.⁵⁵ Densities are obtained from the NIST (webbook.nist.gov). The CG height is $h_{CG} = 40$ m. Densities and IFTs are considered constant over the gas column for each reservoir pressure (or depth) and temperature.</figcaption></figure>

As an example, Fig. 3 displays the maximum column heights H of H$_2$ that can be safely stored on top of a CG with height $h_{CG} = 40$ m, beneath a seal rock with entry pressure $P_{ce,H_2} \approx 1.5$ MPa at T = 298 K and P = 5 MPa. This is a conservative value, 3–5 times smaller than the observed breakthrough pressures into carbonate-rich caprocks of CO$_2$ or nitrogen,⁷⁸ which have less favorable interfacial properties than H$_2$. The maximum storage heights (or storage capacities) are calculated for various reservoir pressures P from 5 to 30 MPa (or storage depths from 500 to 3000 m) and temperatures T = 298 K (cold reservoirs), 323 K (warm reservoirs), and 372 K (hot reservoirs) by means of Eqn 17, where the capillary entry pressure at T and P is inferred from that at T$_{ref}$ = 298 K and P$_{ref}$ = 5 MPa as follows:

$$P_{ce,H_2}(T, P) = P_{ce,H_2}(T_{ref}, P_{ref}) \frac{\gamma_{w,H_2}(T, P)}{\gamma_{w,H_2}(T_{ref}, P_{ref})} \quad (18)$$

In the context of this work, the contact angle $\theta$ is assumed to be independent of T and P. If the dependence of $\theta$ with T and P is known, then the capillary entry pressure at T and P should be corrected as follows:

$$P_{ce,H_2}(T, P) = P_{ce,H_2}(T_{ref}, P_{ref}) \times \frac{\gamma_{w,H_2}(T, P) \cos[\theta(T, P)]}{\gamma_{w,H_2}(T_{ref}, P_{ref}) \cos[\theta(T_{ref}, P_{ref})]} \quad (19)$$

Figure 3 shows the impact of the density difference between the aqueous phase and CG on H$_2$ storage capacity. When H$_2$ is used as a CG, the density difference is the largest and therefore the H$_2$ storage capacity is the smallest. When a dense CG such as CO$_2$ is used, storage capacity is significantly higher, especially under low or moderate temperatures when CO$_2$ is a liquid or a supercritical fluid.

The same trends (data not shown) are observed when expressing the storage capacity in terms of energy per area unit (e.g., A = 1 km$^2$) instead of using the H$_2$ storage capacity (H) directly:

$$\text{Energy per area } E = \rho_{H_2} eH \quad (20)$$

where $e \approx 120$ kJ/g is the energy density of H$_2$.⁷⁹

## Enhancement of storage safety by adding H$_2$ at the top of an existing gas column

When a gas (such as H$_2$) is injected at the crest of a reservoir already filled with a gas (CO$_2$, CH$_4$, N$_2$, or their mixtures), the injected gas becomes the gas in contact with the seal rock. As a rule, H$_2$ has a larger capillary entry pressure into the seal rock than most (if not all) of the gases cited above, that is, $P_{ce,H_2} = \lambda P_{ce,g}$ with $\lambda > 1$. Due to the enhanced capillary entry pressure, an amount (or height H) of H$_2$ can be stored in addition to the gas present without any risk of leakage by capillary failure of the seal rock. The height of the H$_2$ column that can be stored safely on top of the gas column with height $h_g$ is obtained by combining Eqns 17 and 14:

$$H = \frac{(\lambda - 1)(\rho_w - \rho_g)}{(\rho_w - \rho_{H_2})} h_g \quad (21)$$

As an example, consider H$_2$ being injected on top of a column of supercritical CO$_2$ with density $\rho_{CO_2} \approx 500$ kg/m$^3$. From the known IFTs and contact angles (see the subsection Relevant interfacial properties), the coefficient $\lambda$ relating the two entry pressures, $P_{ce,H_2}$ and $P_{ce,CO_2}$, is in the range of 2. The height of H$_2$ that can be
safely stored is at least half the column height $h_{CO_2}$ of the CO₂ already present.

Similar conclusions can be drawn when injecting H₂ at the top of a CH₄-bearing reservoir or CH₄ at the top of a CO₂-bearing reservoir (the coefficient λ is not as large as that for the H₂/CO₂ couple, though). If the added gas is not H₂, the density $ρ_{H_2}$ has to be replaced in the above equation by the density of that gas.

All calculations in this section assume hydrostatic pressures, that is, fluids at rest. Fluid pressure deviates from the law of hydrostatics where fluid velocities are significant – near the bottom hole in injection/ withdrawal sequences. This induces the phenomenon of coning, that is, the deviation from horizontality of the nearby H₂/CG or H₂/water front in the absence of a CG.⁸⁰ When injecting H₂ at the crest of a gas-bearing reservoir (see Fig. 1), attention should be paid to the bottom-hole pressure buildup, as it is a proxy of the pressure buildup beneath the seal rock. If the well was previously used to inject the CG and a record of bottom-hole pressures is available, then this record should be used to determine the maximum bottom-hole pressure not to be exceeded in the ensuing H₂ injection process to avoid failure of the seal rock.

While the risk of capillary failure of the caprock is diminished by adding at the crest of an existing gas reservoir a lighter gas that renders the caprock more sealing (i.e., with a higher capillary entry pressure), the risk of mechanical failure of that rock is indeed not eliminated.²³ The latter risk increases because the pressure at the top of the gas column increases: an analysis of this risk is beyond the scope of this work.

## Conclusions and prospects

There are many benefits in storing hydrogen on top of another gas already present in a geological reservoir, particularly (but not only) when this gas is the CO₂ permanently stored as part of a CCS operation. These benefits are primarily economical: CAPEX costs are strongly reduced as there is no need to invest in the selection and characterization of a new reservoir, in the building of new infrastructures (wells, conduits, etc) and in a large volume of H₂ cushion. This study has examined two key issues, which condition the success of this storage scheme: (i) the extent of mixing between the WG ($H_2$) and the CG and (ii) storage safety and capacity.

(i) The mixing of H₂ with CGs such as CO₂, methane, nitrogen, and their mixtures is less problematic than anticipated: because of the large density contrast but limited viscosity contrast between H₂ and the CG, the mixing front remains stable over a large range of injection rates or velocities, reservoir permeabilities and inclinations of the front with respect to the horizontal. If the “bubble” of injected gas is expanded around the well close enough to the top seal, then a stable mixing front and therefore mixing limited to dispersion are expected.

(ii) As to storage safety and capacity:

- On the one hand, installing H₂ at the crest of an existing gas storage reinforces the sealing ability of the caprock. This benefit is the largest for CO₂-bearing reservoirs as CO₂ has the lowest IFTs with brines and alters the most water-wet character of caprock minerals. Some significant amounts of H₂ can be stored without compromising the safety of both storages. More generally, any lighter gas added on top of CO₂ (or another gas), provided it has a higher entry pressure into the seal rock, enhances the storage safety of the two gases.

- On the other hand, the trade-off between buoyancy and the capillary properties of the seal rock determines H₂ storage capacities. Denser CGs reduce buoyancy at the crest of the reservoir and thus increase storage capacities in the following order: CH₄, N₂, and CO₂. The use of CO₂ as CG for H₂ storage in aquifers or depleted gas reservoirs offers the most significant benefits: not only does it save a large amount of H₂, but it also enhances the safety of storage by reducing the buoyant pressure transmitted to the crest of the reservoir. The use of CO₂ as CG is particularly advantageous due to its high density, relatively close to that of water or brine, especially under the low or moderate temperatures and high pressures encountered in deep offshore reservoirs.

Further research is needed in various directions. One direction consists of better understanding and characterizing the role of heterogeneities in the mixing process – the results presently in this paper are obtained for homogeneous or weakly heterogeneous reservoirs. Another direction should address means to
safely stored is at least half the column height $h_{CO_2}$ of the CO₂ already present.

Similar conclusions can be drawn when injecting H₂ at the top of a CH₄-bearing reservoir or CH₄ at the top of a CO₂-bearing reservoir (the coefficient λ is not as large as that for the H₂/CO₂ couple, though). If the added gas is not H₂, the density ρH₂ has to be replaced in the above equation by the density of that gas.

All calculations in this section assume hydrostatic pressures, that is, fluids at rest. Fluid pressure deviates from the law of hydrostatics where fluid velocities are significant – near the bottom hole in injection/ withdrawal sequences. This induces the phenomenon of coning, that is, the deviation from horizontality of the nearby H₂/CG or H₂/water front in the absence of a CG.⁸⁰ When injecting H₂ at the crest of a gas-bearing reservoir (see Fig. 1), attention should be paid to the bottom-hole pressure buildup, as it is a proxy of the pressure buildup beneath the seal rock. If the well was previously used to inject the CG and a record of bottom-hole pressures is available, then this record should be used to determine the maximum bottom-hole pressure not to be exceeded in the ensuing H₂ injection process to avoid failure of the seal rock.

While the risk of capillary failure of the caprock is diminished by adding at the crest of an existing gas reservoir a lighter gas that renders the caprock more sealing (i.e., with a higher capillary entry pressure), the risk of mechanical failure of that rock is indeed not eliminated.²³ The latter risk increases because the pressure at the top of the gas column increases: an analysis of this risk is beyond the scope of this work.

## Conclusions and prospects

There are many benefits in storing hydrogen on top of another gas already present in a geological reservoir, particularly (but not only) when this gas is the CO₂ permanently stored as part of a CCS operation. These benefits are primarily economical: CAPEX costs are strongly reduced as there is no need to invest in the selection and characterization of a new reservoir, in the building of new infrastructures (wells, conduits, etc) and in a large volume of H₂ cushion. This study has examined two key issues, which condition the success of this storage scheme: (i) the extent of mixing between the WG (H₂) and the CG and (ii) storage safety and capacity.

(i) The mixing of H₂ with CGs such as CO₂, methane, nitrogen, and their mixtures is less problematic than anticipated: because of the large density contrast but limited viscosity contrast between H₂ and the CG, the mixing front remains stable over a large range of injection rates or velocities, reservoir permeabilities and inclinations of the front with respect to the horizontal. If the “bubble” of injected gas is expanded around the well close enough to the top seal, then a stable mixing front and therefore mixing limited to dispersion are expected.

(ii) As to storage safety and capacity:

- On the one hand, installing H₂ at the crest of an existing gas storage reinforces the sealing ability of the caprock. This benefit is the largest for CO₂-bearing reservoirs as CO₂ has the lowest IFTs with brines and alters the most water-wet character of caprock minerals. Some significant amounts of H₂ can be stored without compromising the safety of both storages. More generally, any lighter gas added on top of CO₂ (or another gas), provided it has a higher entry pressure into the seal rock, enhances the storage safety of the two gases.

- On the other hand, the trade-off between buoyancy and the capillary properties of the seal rock determines H₂ storage capacities. Denser CGs reduce buoyancy at the crest of the reservoir and thus increase storage capacities in the following order: CH₄, N₂, and CO₂. The use of CO₂ as CG for H₂ storage in aquifers or depleted gas reservoirs offers the most significant benefits: not only does it save a large amount of H₂, but it also enhances the safety of storage by reducing the buoyant pressure transmitted to the crest of the reservoir. The use of CO₂ as CG is particularly advantageous due to its high density, relatively close to that of water or brine, especially under the low or moderate temperatures and high pressures encountered in deep offshore reservoirs.

Further research is needed in various directions. One direction consists of better understanding and characterizing the role of heterogeneities in the mixing process – the results presently in this paper are obtained for homogeneous or weakly heterogeneous reservoirs. Another direction should address means to
mitigate near-well enhanced mixing in the early stages of an H₂ injection process in a gas-bearing reservoir. Reservoir/well configurations different from the one considered in this work – one vertical well crossing a laterally – and top-bounded reservoir at its crest – should also be examined, as well as the mechanical risks due to increased buoyancy pressure.

## Acknowledgements

This work has benefited from numerous discussions with Fernanda Veloso (BRGM), Farid Smai (BRGM), Pierre Chiquet (TEREGA), Alexandre Dos Santos (TEREGA), Roland Masson (University of Nice) and Sylvain Thibeau (TotalEnergies). The authors express their gratitude to TEREGA for their generous support to help funding the PhD of SB.

## Conflict of interest statement

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.

## References

1. Hiesl A, Ajanovic A, Haas R. On current and future economics of electricity storage. Greenhouse Gases Sci Technol. 2020; 10 (6): 1176-92. https://doi.org/10.1002/ghg.2030

2. Heinemann N, Wilkinson M, Adie K, Edlmann K, Thaysen EM, Hassanpouryouzband A, Haszeldine RS. Cushion Gas in hydrogen storage—a costly CAPEX or a valuable resource for energy crises? Hydrogen. 2022;3(4):550–63. https://doi.org/10.3390/hydrogen3040035

3. Scafidi J, Schirrer L, Vervoort I, Heinemann N. An open-source tool for the calculation of field deliverability and cushion-gas requirements in volumetric gas-reservoir storage sites. Geol Soc, Lond, Spl Publ. 2023; 528 (1): SP528-2022-71. https://doi.org/10.1144/SP528-2022-71

4. Carrière JF, Fasanino G, Tek MR. Mixing in underground storage reservoirs. 1985:SPE-14202-MS. https://doi.org/10.2118/14202-MS

5. Walker CJ, Huff RV. Feasibility of inert-gas cushions in gas storages. U.S. Dept. of the Interior, Bureau of Mines; 1964.

6. Misra B, Foh S, Shikari Y, Berry R, Labaune F. The use of inert base gas in underground natural gas storage. Paper presented at the SPE Unconventional Resources Conference/Gas Technology Symposium, Dallas, TX, 1988.

7. Oldenburg CM. Carbon dioxide as cushion gas for natural gas storage. Energy Fuels. 2003; 17 (1): 240–46. https://doi.org/10.1021/ef020162b

8. Oldenburg CM, Pan L. Utilization of CO₂ as cushion gas for porous media compressed air energy storage. Greenhouse Gases Sci Technol. 2013; 3 (2): 124–35. https://doi.org/10.1002/ghg.1332

9. Tremosa J, Jakobsen R, Le Gallo Y. Assessing and modeling hydrogen reactivity in underground hydrogen storage: a review and models simulating the Lobodice town gas storage. Front Energy Res. 2023; 11: 1145978. https://doi.org/10.3389/fenrg.2023.1145978

10. Gholami R. Hydrogen storage in geological porous media: Solubility, mineral trapping, H₂S generation and salt precipitation. J Storage Mater. 2023; 59: 106576. https://doi.org/10.1016/j.est.2022.106576

11. Haddad P, Ranchou-Peyruse M, Guignard M, et al. Geological storage of hydrogen in deep aquifers—an experimental multidisciplinary study. Energy Environ Sci. 2022; 15 (8): 3400–3415.

12. Thaysen EM, Armitage T, Slabon L, Hassanpouryouzband A, Edlmann K. Microbial risk assessment for underground hydrogen storage in porous rocks. Fuel. 2023; 352: 128852. https://doi.org/10.1016/j.fuel.2023.128852

13. Heinemann N, Scafidi J, Pickup G, Thaysen EM, Hassanpouryouzband A, Wilkinson M, et al. Hydrogen storage in saline aquifers: The role of cushion gas for injection and production. Int J Hydrogen Energy. 2021; 46 (79): 39284–96. https://doi.org/10.1016/j.ijhydene.2021.09.174

14. Wang G, Pickup G, Sorbie K, Mackay E. Scaling analysis of hydrogen flow with carbon dioxide cushion gas in subsurface heterogeneous porous media. Int J Hydrogen Energy. 2022; 47 (3): 1752–64. https://doi.org/10.1016/j.ijhydene.2021.10.224

15. Kanaani M, Sedaee B, Asadian-Pakfar M. Role of cushion gas on underground hydrogen storage in depleted oil reservoirs. J Storage Mater. 2022; 45: 103783. https://doi.org/10.1016/j.est.2021.103783

16. Saeed M, Jadhawar P. Optimizing underground hydrogen storage in aquifers: the impact of cushion gas type. Int J Hydrogen Energy. 2023; 52: 1537–49. https://doi.org/10.1016/j.ijhydene.2023.08.352

17. Zhao Q, Wang Y, Chen C. Numerical simulation of the impact of different cushion gases on underground hydrogen storage in aquifers based on an experimentally-benchmarked equation-of-state. Int J Hydrogen Energy. 2024;50(pt D):495–511. https://doi.org/10.1016/j.ijhydene.2023.07.262

18. Bo Z, Hörning S, Underschultz JR, Garnett A, Hurter S. Effects of geological heterogeneity on gas mixing during underground hydrogen storage (UHS) in braided-fluvial reservoirs. Fuel. 2024; 357: 129949. https://doi.org/10.1016/j.fuel.2023.129949

19. Hassanpouryouzband A, Joonaki E, Edlmann K, Haszeldine RS. Offshore geological storage of hydrogen: Is this our best option to achieve net-zero? ACS Energy Lett. 2021;6(6):2181-86. https://doi.org/10.1021/acsenergylett.1c00845

20. Aftab A, Hassanpouryouzband A, Xie Q, Machuca LL, Sarmadivaleh M. Toward a fundamental understanding of geological hydrogen storage. Ind Eng Chem Res. 2022; 61 (9): 3233-53. https://doi.org/10.1021/acs.iecr.1c04380

21. Chuoke R, van Meurs P, Poel van der C. The Instability of Slow, Immiscible, Viscous Liquid-Liquid Displacements in Permeable Media. Trans AIME. 1959; 216 (01): 188-94. https://doi.org/10.2118/1141-G

22. Homsy GM. Viscous fingering in porous media. Annu Rev Fluid Mech. 1987; 19 (1): 271-311. https://doi.org/10.1146/annurev.fl.19.010187.001415

23. Broseta D. Assessing seal rock integrity for CO₂ geological storage purposes. In: Pijaudier-Cabot G, Pereira JM, editors.
mitigate near-well enhanced mixing in the early stages of an H₂ injection process in a gas-bearing reservoir. Reservoir/well configurations different from the one considered in this work – one vertical well crossing a laterally – and top-bounded reservoir at its crest – should also be examined, as well as the mechanical risks due to increased buoyancy pressure.

## Acknowledgements

This work has benefited from numerous discussions with Fernanda Veloso (BRGM), Farid Smai (BRGM), Pierre Chiquet (TEREGA), Alexandre Dos Santos (TEREGA), Roland Masson (University of Nice) and Sylvain Thibeau (TotalEnergies). The authors express their gratitude to TEREGA for their generous support to help funding the PhD of SB.

## Conflict of interest statement

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.

## References

1. Hiesl A, Ajanovic A, Haas R. On current and future economics of electricity storage. Greenhouse Gases Sci Technol. 2020; 10 (6): 1176-92. https://doi.org/10.1002/ghg.2030

2. Heinemann N, Wilkinson M, Adie K, Edlmann K, Thaysen EM, Hassanpouryouzband A, Haszeldine RS. Cushion Gas in hydrogen storage—a costly CAPEX or a valuable resource for energy crises? Hydrogen. 2022;3(4):550–63. https://doi.org/10.3390/hydrogen3040035

3. Scafidi J, Schirrer L, Vervoort I, Heinemann N. An open-source tool for the calculation of field deliverability and cushion-gas requirements in volumetric gas-reservoir storage sites. Geol Soc, Lond, Spl Publ. 2023; 528 (1): SP528-2022-71. https://doi.org/10.1144/SP528-2022-71

4. Carrière JF, Fasanino G, Tek MR. Mixing in underground storage reservoirs. 1985:SPE-14202-MS. https://doi.org/10.2118/14202-MS

5. Walker CJ, Huff RV. Feasibility of inert-gas cushions in gas storages. U.S. Dept. of the Interior, Bureau of Mines; 1964.

6. Misra B, Foh S, Shikari Y, Berry R, Labaune F. The use of inert base gas in underground natural gas storage. Paper presented at the SPE Unconventional Resources Conference/Gas Technology Symposium, Dallas, TX, 1988.

7. Oldenburg CM. Carbon dioxide as cushion gas for natural gas storage. Energy Fuels. 2003; 17 (1): 240–46. https://doi.org/10.1021/ef020162b

8. Oldenburg CM, Pan L. Utilization of CO₂ as cushion gas for porous media compressed air energy storage. Greenhouse Gases Sci Technol. 2013; 3 (2): 124–35. https://doi.org/10.1002/ghg.1332

9. Tremosa J, Jakobsen R, Le Gallo Y. Assessing and modeling hydrogen reactivity in underground hydrogen storage: a review and models simulating the Lobodice town gas storage. Front Energy Res. 2023; 11: 1145978. https://doi.org/10.3389/fenrg.2023.1145978

10. Gholami R. Hydrogen storage in geological porous media: Solubility, mineral trapping, H₂S generation and salt precipitation. J Storage Mater. 2023; 59: 106576. https://doi.org/10.1016/j.est.2022.106576

11. Haddad P, Ranchou-Peyruse M, Guignard M, et al. Geological storage of hydrogen in deep aquifers—an experimental multidisciplinary study. Energy Environ Sci. 2022; 15 (8): 3400–3415.

12. Thaysen EM, Armitage T, Slabon L, Hassanpouryouzband A, Edlmann K. Microbial risk assessment for underground hydrogen storage in porous rocks. Fuel. 2023; 352: 128852. https://doi.org/10.1016/j.fuel.2023.128852

13. Heinemann N, Scafidi J, Pickup G, Thaysen EM, Hassanpouryouzband A, Wilkinson M, et al. Hydrogen storage in saline aquifers: The role of cushion gas for injection and production. Int J Hydrogen Energy. 2021; 46 (79): 39284–96. https://doi.org/10.1016/j.ijhydene.2021.09.174

14. Wang G, Pickup G, Sorbie K, Mackay E. Scaling analysis of hydrogen flow with carbon dioxide cushion gas in subsurface heterogeneous porous media. Int J Hydrogen Energy. 2022; 47 (3): 1752–64. https://doi.org/10.1016/j.ijhydene.2021.10.224

15. Kanaani M, Sedaee B, Asadian-Pakfar M. Role of cushion gas on underground hydrogen storage in depleted oil reservoirs. J Storage Mater. 2022; 45: 103783. https://doi.org/10.1016/j.est.2021.103783

16. Saeed M, Jadhawar P. Optimizing underground hydrogen storage in aquifers: the impact of cushion gas type. Int J Hydrogen Energy. 2023; 52: 1537–49. https://doi.org/10.1016/j.ijhydene.2023.08.352

17. Zhao Q, Wang Y, Chen C. Numerical simulation of the impact of different cushion gases on underground hydrogen storage in aquifers based on an experimentally-benchmarked equation-of-state. Int J Hydrogen Energy. 2024;50(pt D):495–511. https://doi.org/10.1016/j.ijhydene.2023.07.262

18. Bo Z, Hörning S, Underschultz JR, Garnett A, Hurter S. Effects of geological heterogeneity on gas mixing during underground hydrogen storage (UHS) in braided-fluvial reservoirs. Fuel. 2024; 357: 129949. https://doi.org/10.1016/j.fuel.2023.129949

19. Hassanpouryouzband A, Joonaki E, Edlmann K, Haszeldine RS. Offshore geological storage of hydrogen: Is this our best option to achieve net-zero? ACS Energy Lett. 2021;6(6):2181-86. https://doi.org/10.1021/acsenergylett.1c00845

20. Aftab A, Hassanpouryouzband A, Xie Q, Machuca LL, Sarmadivaleh M. Toward a fundamental understanding of geological hydrogen storage. Ind Eng Chem Res. 2022; 61 (9): 3233-53. https://doi.org/10.1021/acs.iecr.1c04380

21. Chuoke R, van Meurs P, Poel van der C. The Instability of Slow, Immiscible, Viscous Liquid-Liquid Displacements in Permeable Media. Trans AIME. 1959; 216 (01): 188-94. https://doi.org/10.2118/1141-G

22. Homsy GM. Viscous fingering in porous media. Annu Rev Fluid Mech. 1987; 19 (1): 271-311. https://doi.org/10.1146/annurev.fl.19.010187.001415

23. Broseta D. Assessing seal rock integrity for CO₂ geological storage purposes. In: Pijaudier-Cabot G, Pereira JM, editors.
Geomechanics in CO₂ storage facilities. Hoboken, New Jersey: John Wiley & Sons, Ltd.; 2013. p. 3–20.
https://doi.org/10.1002/9781118577424.ch1

24. Legatski MW, Katz DL. Dispersion coefficients for gases flowing in consolidated porous media. SPE. 1967; 7 (01): 43–53. https://doi.org/10.2118/1594-PA

25. Flowers TC, Hunt JR. Viscous and gravitational contributions to mixing during vertical brine transport in water-saturated porous media. Water Resour Res. 2007; 43 (1): W01407.

26. Alkindi A, Al-Wahaibi Y, Bijeljic B, Muggeridge A. Investigation of longitudinal and transverse dispersion in stable displacements with a high viscosity and density contrast between the fluids. J Contam Hydrol. 2011; 120-121: 170–83.
https://doi.org/10.1016/j.jconhyd.2010.06.006

27. Yang K, Kobeissi S, Ling N, Li M, Esteban L, May EF, Johns ML. Measurement of hydrogen dispersion in rock cores using benchtop NMR. Int J Hydrogen Energy. 2023; 48 (45): 17251–60. https://doi.org/10.1016/j.ijhydene.2023.01.197

28. Giddings JC, Seager SL. Method for the rapid determination of diffusion coefficients. theory and application. Ind Eng Chem Fundam. 1962; 1 (4): 277–83.

29. Poling BE, Prausnitz JM, O'Connell JP. The properties of gases and liquids. 5th ed. New York: McGraw-Hill; 2001.

30. Gelhar LW, Welty C, Rehfeldt KR. A critical review of data on field-scale dispersion in aquifers. Water Resour Res. 1992; 28 (7): 1955–74. https://doi.org/10.1029/92WR00607

31. Paterson L. The implications of fingering in underground hydrogen storage. Int J Hydrogen Energy. 1983; 8 (1): 53–59.

32. Christie MA. High-resolution simulation of unstable flows in porous media. SPE Reservoir Eng. 1989;4(3): 297–303.
https://doi.org/10.2118/16005-PA

33. Yekta AE, Manceau JC, Gaboreau S, Pichavant M, Audigane P. Determination of hydrogen-water relative permeability and capillary pressure in sandstone: application to underground hydrogen injection in sedimentary formations. Transp Porous Media. 2018; 122 (2): 333–56.
https://doi.org/10.1007/s11242-018-1004-7

34. Dumore J. Stability Considerations in downward miscible displacements. SPE J. 1964; 4 (4): 356–62.
https://doi.org/10.2118/961-PA

35. Manickam O, Homsy GM. Fingering instabilities in vertical miscible displacement flows in porous media. J Fluid Mech. 1995;288:75–102.
https://doi.org/10.1017/S0022112095001078

36. Mahardika MA, She Y, Koe T, Patmonoaji A, Nasir M, Matsushita S, Suekane T. Competition of gravity and viscous forces in miscible vertical displacement in a three-dimensional porous medium. Phys Fluids. 2022; 34 (7): 073102.
https://doi.org/10.1063/5.0090387

37. Wilke C. A viscosity equation for gas mixtures. J Chem Phys. 1950; 18 (4): 517–19.

38. Herning F, Zipperer L. Calculation of the viscosity of technical gas mixtures from the viscosity of the individual gases. Gas Wasserfach. 1936; 79: 69.

39. Sarma HK, Maini BB. An experimental evaluation of viscosity grading for controlling fingering in miscible displacements. J Can Pet Technol. 1993;32(01).

40. Perkins TK, Johnston OC, Hoffman RN. Mechanics of Viscous Fingering in Miscible Systems. SPE J. 1965; 5 (04): 301–17.
https://doi.org/10.2118/1229-PA

41. Gardner J, Ypma J. An investigation of phase-behavior/macroscopic-bypassing interaction in CO2 flooding. SPE J. 1984; 24 (05): 508–20.

42. Kunz O, Wagner W. The GERG-2008 Wide-range equation of state for natural gases and other mixtures: an expansion of GERG-2004. J Chem Eng Data. 2012; 57 (11): 3032–91.
https://doi.org/10.1021/je300655b

43. Sanchez-Vicente Y, Drage TC, Poliakoff M, Ke J, George MW. Densities of the carbon dioxide+hydrogen, a system of relevance to carbon capture and storage. Int J Greenhouse Gas Control. 2013; 13: 78–86.
https://doi.org/10.1016/j.ijggc.2012.12.002

44. Khosravi B, Betken B, Jakobsen JP, Løvseth SW, Span R. Viscosity measurements of CO₂-rich; CO₂ + N₂ and CO₂ + H₂ mixtures in gas or supercritical phase at temperatures between 273 and 473 K and pressures up to 8.7 MPa. Fluid Phase Equilib. 2022; 560: 113519.
https://doi.org/10.1016/j.fluid.2022.113519

45. Araktingi UG, Orr FM. Viscous fingering in heterogeneous porous media. SPE Adv Technol Sr. 1993; 1 (1): 71–80.
https://doi.org/10.2118/18095-PA

46. Tchelepi HA, Orr J. Interaction of viscous fingering, permeability heterogeneity, and gravity segregation in three dimensions. SPE Reservoir Eng. 1994; 9 (04): 266–71.
https://doi.org/10.2118/25235-PA

47. Zhou D, Fayers F, Orr Jr F. Scaling of multiphase flow in simple heterogeneous porous media. Paper SPE-27833-MS presented at the SPE/DOE Improved Oil Recovery Symposium, April 17–20, 1994.
https://doi.org/10.2118/27833-MS

48. Fandiño O, Trusler JM, Vega-Maza D. Phase behavior of (CO₂+H₂) and (CO₂+N₂) at temperatures between (218.15 and 303.15)K at pressures up to 15MPa. Int J Greenhouse Gas Control. 2015; 36: 78–92.
https://doi.org/10.1016/j.ijggc.2015.02.018

49. Chiquet P, Daridon JL, Broseta D, Thibeau S. CO₂/water interfacial tensions under pressure and temperature conditions of CO₂ geological storage. Energy Convers Manage. 2007; 48 (3): 736–44. https://doi.org/10.1016/j.enconman.2006.09.011

50. Naylor M, Wilkinson M, Haszeldine R. Calculation of CO₂ column heights in depleted gas fields from known pre-production gas column heights. Mar Pet Geol. 2011; 28 (5): 1083–93. https://doi.org/10.1016/j.marpetgeo.2010.10.005

51. Aveyard R, Clint J, Neumann A. The measurement of surface tensions and contact angles at elevated pressures and temperatures. J Colloid Interface Sci. 1973; 42 (2): 363–71.

52. Chen Z, Yang D. Correlations/estimation of equilibrium interfacial tension for methane/CO₂-water/brine systems based on mutual solubility. Fluid Phase Equilib. 2019; 483: 197–208. https://doi.org/10.1016/j.fluid.2018.11.037

53. Massoudi R, King Jr A. Effect of pressure on the surface tension of water. Adsorption of low molecular weight gases on water at 25 deg. J Phys Chem. 1974; 78 (22): 2262–66.

54. Hosseini M, Fahimpour J, Ali M, Keshavarz A, Iglauer S. H₂-brine interfacial tension as a function of salinity, temperature, and pressure; implications for hydrogen geo-storage. J Petr Sci Eng. 2022; 213: 110441.
https://doi.org/10.1016/j petrol.2022.110441

55. Chow YF, Maitland GC, Trusler JM. Interfacial tensions of (H₂O + H₂) and (H₂O + CO₂ + H₂) systems at temperatures
Geomechanics in CO₂ storage facilities. Hoboken, New Jersey: John Wiley & Sons, Ltd.; 2013. p. 3-20.
https://doi.org/10.1002/9781118577424.ch1

24. Legatski MW, Katz DL. Dispersion coefficients for gases flowing in consolidated porous media. SPE. 1967; 7 (01): 43-53. https://doi.org/10.2118/1594-PA

25. Flowers TC, Hunt JR. Viscous and gravitational contributions to mixing during vertical brine transport in water-saturated porous media. Water Resour Res. 2007; 43 (1): W01407.

26. Alkindi A, Al-Wahaibi Y, Bijeljic B, Muggeridge A. Investigation of longitudinal and transverse dispersion in stable displacements with a high viscosity and density contrast between the fluids. J Contam Hydrol. 2011; 120-121: 170-83.
https://doi.org/10.1016/j.jconhyd.2010.06.006

27. Yang K, Kobeissi S, Ling N, Li M, Esteban L, May EF, Johns ML. Measurement of hydrogen dispersion in rock cores using benchtop NMR. Int J Hydrogen Energy. 2023; 48 (45): 17251-60. https://doi.org/10.1016/j.ijhydene.2023.01.197

28. Giddings JC, Seager SL. Method for the rapid determination of diffusion coefficients. theory and application. Ind Eng Chem Fundam. 1962; 1 (4): 277-83.

29. Poling BE, Prausnitz JM, O'Connell JP. The properties of gases and liquids. 5th ed. New York: McGraw-Hill; 2001.

30. Gelhar LW, Welty C, Rehfeldt KR. A critical review of data on field-scale dispersion in aquifers. Water Resour Res. 1992; 28 (7): 1955-74. https://doi.org/10.1029/92WR00607

31. Paterson L. The implications of fingering in underground hydrogen storage. Int J Hydrogen Energy. 1983; 8 (1): 53-59.

32. Christie MA. High-resolution simulation of unstable flows in porous media. SPE Reservoir Eng. 1989;4(3): 297-303.
https://doi.org/10.2118/16005-PA

33. Yekta AE, Manceau JC, Gaboreau S, Pichavant M, Audigane P. Determination of hydrogen-water relative permeability and capillary pressure in sandstone: application to underground hydrogen injection in sedimentary formations. Transp Porous Media. 2018; 122 (2): 333-56.
https://doi.org/10.1007/s11242-018-1004-7

34. Dumore J. Stability Considerations in downward miscible displacements. SPE J. 1964; 4 (4): 356-62.
https://doi.org/10.2118/961-PA

35. Manickam O, Homsy GM. Fingering instabilities in vertical miscible displacement flows in porous media. J Fluid Mech. 1995;288:75-102.
https://doi.org/10.1017/S0022112095001078

36. Mahardika MA, She Y, Koe T, Patmonoaji A, Nasir M, Matsushita S, Suekane T. Competition of gravity and viscous forces in miscible vertical displacement in a three-dimensional porous medium. Phys Fluids. 2022; 34 (7): 073102.
https://doi.org/10.1063/5.0090387

37. Wilke C. A viscosity equation for gas mixtures. J Chem Phys. 1950; 18 (4): 517-19.

38. Herning F, Zipperer L. Calculation of the viscosity of technical gas mixtures from the viscosity of the individual gases. Gas Wasserfach. 1936; 79: 69.

39. Sarma HK, Maini BB. An experimental evaluation of viscosity grading for controlling fingering in miscible displacements. J Can Pet Technol. 1993;32(01).

40. Perkins TK, Johnston OC, Hoffman RN. Mechanics of Viscous Fingering in Miscible Systems. SPE J. 1965; 5 (04): 301-17.
https://doi.org/10.2118/1229-PA

41. Gardner J, Ypma J. An investigation of phase-behavior/macroscopic-bypassing interaction in CO2 flooding. SPE J. 1984; 24 (05): 508-20.

42. Kunz O, Wagner W. The GERG-2008 Wide-range equation of state for natural gases and other mixtures: an expansion of GERG-2004. J Chem Eng Data. 2012; 57 (11): 3032-91.
https://doi.org/10.1021/je300655b

43. Sanchez-Vicente Y, Drage TC, Poliakoff M, Ke J, George MW. Densities of the carbon dioxide+hydrogen, a system of relevance to carbon capture and storage. Int J Greenhouse Gas Control. 2013; 13: 78-86.
https://doi.org/10.1016/j.ijggc.2012.12.002

44. Khosravi B, Betken B, Jakobsen JP, Løvseth SW, Span R. Viscosity measurements of CO₂-rich; CO₂ + N₂ and CO₂ + H₂ mixtures in gas or supercritical phase at temperatures between 273 and 473 K and pressures up to 8.7 MPa. Fluid Phase Equilib. 2022; 560: 113519.
https://doi.org/10.1016/j.fluid.2022.113519

45. Araktingi UG, Orr FM. Viscous fingering in heterogeneous porous media. SPE Adv Technol Sr. 1993; 1 (1): 71-80.
https://doi.org/10.2118/18095-PA

46. Tchelepi HA, Orr J. Interaction of viscous fingering, permeability heterogeneity, and gravity segregation in three dimensions. SPE Reservoir Eng. 1994; 9 (04): 266-71.
https://doi.org/10.2118/25235-PA

47. Zhou D, Fayers F, Orr Jr F. Scaling of multiphase flow in simple heterogeneous porous media. Paper SPE-27833-MS presented at the SPE/DOE Improved Oil Recovery Symposium, April 17-20, 1994.
https://doi.org/10.2118/27833-MS

48. Fandiño O, Trusler JM, Vega-Maza D. Phase behavior of (CO₂+H₂) and (CO₂+N₂) at temperatures between (218.15 and 303.15)K at pressures up to 15MPa. Int J Greenhouse Gas Control. 2015; 36: 78-92.
https://doi.org/10.1016/j.ijggc.2015.02.018

49. Chiquet P, Daridon JL, Broseta D, Thibeau S. CO₂/water interfacial tensions under pressure and temperature conditions of CO₂ geological storage. Energy Convers Manage. 2007; 48 (3): 736-44. https://doi.org/10.1016/j.enconman.2006.09.011

50. Naylor M, Wilkinson M, Haszeldine R. Calculation of CO₂ column heights in depleted gas fields from known pre-production gas column heights. Mar Pet Geol. 2011; 28 (5): 1083-93. https://doi.org/10.1016/j.marpetgeo.2010.10.005

51. Aveyard R, Clint J, Neumann A. The measurement of surface tensions and contact angles at elevated pressures and temperatures. J Colloid Interface Sci. 1973; 42 (2): 363-71.

52. Chen Z, Yang D. Correlations/estimation of equilibrium interfacial tension for methane/CO₂-water/brine systems based on mutual solubility. Fluid Phase Equilib. 2019; 483: 197-208. https://doi.org/10.1016/j.fluid.2018.11.037

53. Massoudi R, King Jr A. Effect of pressure on the surface tension of water. Adsorption of low molecular weight gases on water at 25 deg. J Phys Chem. 1974; 78 (22): 2262-66.

54. Hosseini M, Fahimpour J, Ali M, Keshavarz A, Iglauer S. H₂-brine interfacial tension as a function of salinity, temperature, and pressure; implications for hydrogen geo-storage. J Petr Sci Eng. 2022; 213: 110441.
https://doi.org/10.1016/j petrol.2022.110441

55. Chow YF, Maitland GC, Trusler JM. Interfacial tensions of (H₂O + H₂) and (H₂O + CO₂ + H₂) systems at temperatures
of (298-448) K and pressures up to 45 MPa. Fluid Phase Equilib. 2018; 475: 37-44.
https://doi.org/10.1016/j.fluid.2018.07.022

56. Higgs S, Da Wang Y, Sun C, Ennis-King J, Jackson SJ, Armstrong RT, Mostaghimi P. In-situ hydrogen wettability characterisation for underground hydrogen storage. Int J Hydrogen Energy. 2022; 47 (26): 13062-75.
https://doi.org/10.1016/j.ijhydene.2022.02.022

57. Al-Mukainah H, Al-Yaseri A, Yekeen N, Hamad JA, Mahmoud M. Wettability of shale-brine-H2 system and H2-brine interfacial tension for assessment of the sealing capacities of shale formations during underground hydrogen storage. Energy Rep. 2022; 8: 8830-43.
https://doi.org/10.1016/j.egyr.2022.07.004

58. Omrani S, Ghasemi M, Singh M, Mahmoodpour S, Zhou T, Babaei M, Niasar V. Interfacial tension-temperature-pressure-salinity relationship for the hydrogen-brine system under reservoir conditions: integration of molecular dynamics and machine learning. Langmuir. 2023; 39 (36): 12680-91.
https://doi.org/10.1021/acs.langmuir.3c01424

59. Pereira LM, Chapoy A, Burgass R, Tohidi B. Interfacial tension of CO2+brine systems: Experiments and predictive modeling. Adv Water Resour. 2017; 103: 64-75.
https://doi.org/10.1016/j.adwatres.2017.02.015

60. Mutailipu M, Jiang L, Fu J, Wang Z, Yu T, Lu Z, Liu Y. Effects of Na+, K+, Ca2+, and Mg2+ cations on CO2-brine interfacial tension under offshore storage conditions. Greenhouse Gases Sci Technol. 2018; 8 (4): 762-80.
https://doi.org/10.1002/ghg.1787

61. Wiegand G, Franck E. Interfacial tension between water and non-polar fluids up to 473 K and 2800 bar. Ber Bunsen Ges Phys Chem. 1994; 98 (6): 809-17.

62. Yan W, Zhao GY, Chen GJ, Guo TM. Interfacial tension of (methane + nitrogen) + water and (carbon dioxide + nitrogen) + water systems. J Chem Eng Data. 2001; 46 (6): 1544-48.
https://doi.org/10.1021/je0101505

63. Zhang C, Wang M. CO2/brine interfacial tension for geological CO2 storage: A systematic review. J Petro Sci Eng. 2023; 220: 111154.
https://doi.org/10.1016/j.petrol.2022.111154

64. Mutailipu M, Song Y, Yao Q, Liu Y, Martin Trusler J. Solubility and interfacial tension models for CO2-brine systems under CO2 geological storage conditions. Fuel. 2024; 357: 129712.
https://doi.org/10.1016/j.fuel.2023.129712

65. Mouallem J, Arif M, Raza A, et al. Critical review and meta-analysis of the interfacial tension of CO2-brine and H2-brine systems: implications for CO2 and H2 geo-storage. Fuel. 2024; 356: 129575.
https://doi.org/10.1016/j.fuel.2023.129575

66. Ali A, Cole DR, Striolo A. Cushion gas effects on clay-hydrogen-brine wettability at conditions relevant to underground gas storage. Int J Hydrogen Energy. 2024; 58: 668-77.
https://doi.org/10.1016/j.ijhydene.2024.01.151

67. Mirchi V, Dejam M, Alvarado V. Interfacial tension and contact angle measurements for hydrogen-methane mixtures/brine/oil-wet rocks at reservoir conditions. Int J Hydrogen Energy. 2022; 47 (82): 34963-75.
https://doi.org/10.1016/j.ijhydene.2022.08.056

68. Dalal Isfehani Z, Sheidae A, Hosseini M, Fahimpour J, Iglauer S, Keshavarz A. Interfacial tensions of (brine + H2 + CO2) systems at gas geo-storage conditions. J Mol Liq. 2023; 374: 121279.
https://doi.org/10.1016/j.molliq.2023.121279

69. Muhammed NS, Haq B, Al Shehri D. CO2 rich cushion gas for hydrogen storage in depleted gas reservoirs: insight on contact angle and surface tension. Int J Hydrogen Energy. 2024; 50: 1281-301.
https://doi.org/10.1016/j.ijhydene.2023.09.148

70. Arif M, Barifcani A, Iglauer S. Solid/CO2 and solid/water interfacial tensions as a function of pressure, temperature, salinity and mineral type: Implications for CO2-wettability and CO2 geo-storage. Int J Greenhouse Gas Control. 2016; 53: 263-73.
https://doi.org/10.1016/j.ijggc.2016.08.020

71. Pan B, Yin X, Iglauer S. Rock-fluid interfacial tension at subsurface conditions: Implications for H2, CO2 and natural gas geo-storage. Int J Hydrogen Energy. 2021; 46 (50): 25578-85.
https://doi.org/10.1016/j.ijhydene.2021.05.067

72. Ali M, Pan B, Yekeen N, et al. Assessment of wettability and rock-fluid interfacial tension of caprock: implications for hydrogen and carbon dioxide geo-storage. Int J Hydrogen Energy. 2022; 47 (30): 14104-20.
https://doi.org/10.1016/j.ijhydene.2022.02.149

73. Esfandyari H, Hosseini M, Ali M, Iglauer S, Haghighi M, Keshavarz A. Assessment of the interfacial properties of various mineral/hydrogen/water systems. J Storage Mater. 2023; 60: 106637.
https://doi.org/10.1016/j.est.2023.106637

74. Chen C, Xia J. A comparative study on transport and interfacial physics of H2/CO2/CH4 interacting with H2O and/or silica by molecular dynamics simulation. Phys Fluids. 2024; 36 (1): 016606.
https://doi.org/10.1063/5.0184754

75. Iglauer S, Pentland CH, Busch A. CO2 wettability of seal and reservoir rocks and the implications for carbon geo-sequestration. Water Resour Res. 2015; 51 (1): 729-74.
https://doi.org/10.1002/2014WR015553

76. Yekeen N, Padmanabhan E, Sevoo TA, Kanesen KA, Okunade OA. Wettability of rock/CO2/brine systems: a critical review of influencing parameters and recent advances. J Ind Eng Chem. 2020; 88: 1-28.
https://doi.org/10.1016/j.jiec.2020.03.021

77. Ali M, Jha NK, Al-Yaseri A, Zhang Y, Iglauer S, Sarmadivaleh M. Hydrogen wettability of quartz substrates exposed to organic acids; Implications for hydrogen geo-storage in sandstone reservoirs. J Pet Sci Eng. 2021; 207: 109081.
https://doi.org/10.1016/j.petrol.2021.109081

78. Tonnet N, Mouronval G, Chiquet P, Broseta D. Petrophysical assessment of a carbonate-rich caprock for CO2 geological storage purposes. Energy Procedia. 2011; 4: 5422-29. 10th International Conference on Greenhouse Gas Control Technologies https://doi.org/10.1016/j.egypro.2011.02.527

79. National Institute of Standards and Technology (NIST). Thermophysical properties of fluid systems. 2023. https://webbook.nist.gov/chemistry/fluid/ [last accessed March 21, 2023]

80. Oldenburg CM, Finsterle S, Trautz RC. Water upconing in underground hydrogen storage: sensitivity analysis to inform design of withdrawal. Transp Porous Media. 2023;151:1-30.
https://doi.org/10.1007/s11242-024-02066-z

81. Ely JF, Hanley H. Prediction of transport properties. 1. Viscosity of fluids and mixtures. Ind Eng Chem Fundam. 1981; 20 (4): 323-32.

82. Millat J, Dymond J, Castro dCN, Wakeham W. Transport properties of fluids. Cambridge, UK: Cambridge University Press; 1996.

83. Huber ML, Lemmon EW, Bell IH, McLinden MO. The NIST REFPROP database for highly accurate properties of
of (298-448) K and pressures up to 45 MPa. Fluid Phase Equilib. 2018; 475: 37-44.
https://doi.org/10.1016/j.fluid.2018.07.022

56. Higgs S, Da Wang Y, Sun C, Ennis-King J, Jackson SJ, Armstrong RT, Mostaghimi P. In-situ hydrogen wettability characterisation for underground hydrogen storage. Int J Hydrogen Energy. 2022; 47 (26): 13062-75.
https://doi.org/10.1016/j.ijhydene.2022.02.022

57. Al-Mukainah H, Al-Yaseri A, Yekeen N, Hamad JA, Mahmoud M. Wettability of shale-brine-H2 system and H2-brine interfacial tension for assessment of the sealing capacities of shale formations during underground hydrogen storage. Energy Rep. 2022; 8: 8830-43.
https://doi.org/10.1016/j.egyr.2022.07.004

58. Omrani S, Ghasemi M, Singh M, Mahmoodpour S, Zhou T, Babaei M, Niasar V. Interfacial tension-temperature-pressure-salinity relationship for the hydrogen-brine system under reservoir conditions: integration of molecular dynamics and machine learning. Langmuir. 2023; 39 (36): 12680-91.
https://doi.org/10.1021/acs.langmuir.3c01424

59. Pereira LM, Chapoy A, Burgass R, Tohidi B. Interfacial tension of CO2+brine systems: Experiments and predictive modeling. Adv Water Resour. 2017; 103: 64-75.
https://doi.org/10.1016/j.adwatres.2017.02.015

60. Mutailipu M, Jiang L, Fu J, Wang Z, Yu T, Lu Z, Liu Y. Effects of Na+, K+, Ca2+, and Mg2+ cations on CO2-brine interfacial tension under offshore storage conditions. Greenhouse Gases Sci Technol. 2018; 8 (4): 762-80.
https://doi.org/10.1002/ghg.1787

61. Wiegand G, Franck E. Interfacial tension between water and non-polar fluids up to 473 K and 2800 bar. Ber Bunsen Ges Phys Chem. 1994; 98 (6): 809-17.

62. Yan W, Zhao GY, Chen GJ, Guo TM. Interfacial tension of (methane + nitrogen) + water and (carbon dioxide + nitrogen) + water systems. J Chem Eng Data. 2001; 46 (6): 1544-48.
https://doi.org/10.1021/je0101505

63. Zhang C, Wang M. CO2/brine interfacial tension for geological CO2 storage: A systematic review. J Petro Sci Eng. 2023; 220: 111154.
https://doi.org/10.1016/j.petrol.2022.111154

64. Mutailipu M, Song Y, Yao Q, Liu Y, Martin Trusler J. Solubility and interfacial tension models for CO2-brine systems under CO2 geological storage conditions. Fuel. 2024; 357: 129712.
https://doi.org/10.1016/j.fuel.2023.129712

65. Mouallem J, Arif M, Raza A, et al. Critical review and meta-analysis of the interfacial tension of CO2-brine and H2-brine systems: implications for CO2 and H2 geo-storage. Fuel. 2024; 356: 129575.
https://doi.org/10.1016/j.fuel.2023.129575

66. Ali A, Cole DR, Striolo A. Cushion gas effects on clay-hydrogen-brine wettability at conditions relevant to underground gas storage. Int J Hydrogen Energy. 2024; 58: 668-77.
https://doi.org/10.1016/j.ijhydene.2024.01.151

67. Mirchi V, Dejam M, Alvarado V. Interfacial tension and contact angle measurements for hydrogen-methane mixtures/brine/oil-wet rocks at reservoir conditions. Int J Hydrogen Energy. 2022; 47 (82): 34963-75.
https://doi.org/10.1016/j.ijhydene.2022.08.056

68. Dalal Isfehani Z, Sheidae A, Hosseini M, Fahimpour J, Iglauer S, Keshavarz A. Interfacial tensions of (brine + H2 + CO2) systems at gas geo-storage conditions. J Mol Liq. 2023; 374: 121279.
https://doi.org/10.1016/j.molliq.2023.121279

69. Muhammed NS, Haq B, Al Shehri D. CO2 rich cushion gas for hydrogen storage in depleted gas reservoirs: insight on contact angle and surface tension. Int J Hydrogen Energy. 2024; 50: 1281-301.
https://doi.org/10.1016/j.ijhydene.2023.09.148

70. Arif M, Barifcani A, Iglauer S. Solid/CO2 and solid/water interfacial tensions as a function of pressure, temperature, salinity and mineral type: Implications for CO2-wettability and CO2 geo-storage. Int J Greenhouse Gas Control. 2016; 53: 263-73.
https://doi.org/10.1016/j.ijggc.2016.08.020

71. Pan B, Yin X, Iglauer S. Rock-fluid interfacial tension at subsurface conditions: Implications for H2, CO2 and natural gas geo-storage. Int J Hydrogen Energy. 2021; 46 (50): 25578-85.
https://doi.org/10.1016/j.ijhydene.2021.05.067

72. Ali M, Pan B, Yekeen N, et al. Assessment of wettability and rock-fluid interfacial tension of caprock: implications for hydrogen and carbon dioxide geo-storage. Int J Hydrogen Energy. 2022; 47 (30): 14104-20.
https://doi.org/10.1016/j.ijhydene.2022.02.149

73. Esfandyari H, Hosseini M, Ali M, Iglauer S, Haghighi M, Keshavarz A. Assessment of the interfacial properties of various mineral/hydrogen/water systems. J Storage Mater. 2023; 60: 106637.
https://doi.org/10.1016/j.est.2023.106637

74. Chen C, Xia J. A comparative study on transport and interfacial physics of H2/CO2/CH4 interacting with H2O and/or silica by molecular dynamics simulation. Phys Fluids. 2024; 36 (1): 016606.
https://doi.org/10.1063/5.0184754

75. Iglauer S, Pentland CH, Busch A. CO2 wettability of seal and reservoir rocks and the implications for carbon geo-sequestration. Water Resour Res. 2015; 51 (1): 729-74.
https://doi.org/10.1002/2014WR015553

76. Yekeen N, Padmanabhan E, Sevoo TA, Kanesen KA, Okunade OA. Wettability of rock/CO2/brine systems: a critical review of influencing parameters and recent advances. J Ind Eng Chem. 2020; 88: 1-28.
https://doi.org/10.1016/j.jiec.2020.03.021

77. Ali M, Jha NK, Al-Yaseri A, Zhang Y, Iglauer S, Sarmadivaleh M. Hydrogen wettability of quartz substrates exposed to organic acids; Implications for hydrogen geo-storage in sandstone reservoirs. J Pet Sci Eng. 2021; 207: 109081.
https://doi.org/10.1016/j.petrol.2021.109081

78. Tonnet N, Mouronval G, Chiquet P, Broseta D. Petrophysical assessment of a carbonate-rich caprock for CO2 geological storage purposes. Energy Procedia. 2011; 4: 5422-29. 10th International Conference on Greenhouse Gas Control Technologies https://doi.org/10.1016/j.egypro.2011.02.527

79. National Institute of Standards and Technology (NIST). Thermophysical properties of fluid systems. 2023. https://webbook.nist.gov/chemistry/fluid/ [last accessed March 21, 2023]

80. Oldenburg CM, Finsterle S, Trautz RC. Water upconing in underground hydrogen storage: sensitivity analysis to inform design of withdrawal. Transp Porous Media. 2023;151:1-30.
https://doi.org/10.1007/s11242-024-02066-z

81. Ely JF, Hanley H. Prediction of transport properties. 1. Viscosity of fluids and mixtures. Ind Eng Chem Fundam. 1981; 20 (4): 323-32.

82. Millat J, Dymond J, Castro dCN, Wakeham W. Transport properties of fluids. Cambridge, UK: Cambridge University Press; 1996.

83. Huber ML, Lemmon EW, Bell IH, McLinden MO. The NIST REFPROP database for highly accurate properties of
industrially important fluids. Ind Eng Chem Res. 2022; 61 (42): 15449-72.

84. Chichester JC, Huber ML. Documentation and assessment of transport property model for mixtures implemented in NIST REFPROP (Version 8.0). US Department of Commerce, Technology Administration, National Institute of Standards and Technology; 2008.

85. Dean DE, Stiel LI. The viscosity of nonpolar gas mixtures at moderate and high pressures. AIChE J. 1965;11(3):526-32.

86. Predvoditelev A. Some invariant quantities in the heat-conduction and viscosity theory. Zh Fiz Khim. 1948; 22 (3): 339-48.

87. Galliéro G, Boned C, Baylaucq A. Molecular dynamics study of the Lennard-Jones fluid viscosity: application to real fluids. Ind Eng Chem Res. 2005; 44 (17): 6963-72.

88. Brokaw RS. Predicting transport properties of dilute gases. Ind Eng Chem Process Des Dev. 1969; 8 (2): 240-53.

# Appendix A: Models of density versus viscosity behavior for gas mixtures

In this appendix, we start by introducing models for the viscosity of gas mixtures that have a larger range of validity than the low-pressure (or low-density) approximation (Wilke's model) presented and used in the subsection 'Rate conditions for viscous fingering avoidance (front stability)'. Subsequently, we compare the results of those models for the gas mixtures of interest in this work, focusing on how viscosity and density are related since this relation determines the critical front velocity (see Eqn 5). Densities are calculated by the GERG-2008 model,⁴² which also has a larger range of validity than the low-pressure, IG approximation.

Those models are then used to assess the validity of the low-pressure (IG+Wilke) approximation for the gas mixtures and conditions of interest, with a focus on how density and viscosity are related in the H₂-rich region. Finally, mixtures of H₂ and CO₂ are examined over the whole interval of compositions.

## A.1. Viscosities of gas mixtures

A variety of methods have been developed to accurately estimate the fluid viscosities of pure components and their mixtures. These methods present advantages and limitations generally linked to their type, either “empirical” requiring a lot of experimental data or “predictive” with a solid theoretical background and requiring little or no experimental data. Given the very limited available high-pressure viscosity data especially for the H₂ + CO₂ system, predictive models are the

most appropriate in this case. The best known predictive models used to estimate the viscosities of high-pressure gas mixtures are

* the TRAPP (Transport Property Prediction) -type model (Eqn A.1) developed by Ely and Hanley⁸¹ based on the ECS theory and applicable to systems containing hydrocarbons, permanent gases, CO₂, H₂S, and refrigerants. However, the original version (TRAPP) is not recommended for use at reduced temperatures (T/T_c) above 0.925, and in these conditions, its variant called SuperTRAPP (see Millat et al.⁸²), which has been further refined and implemented under the generic name “ECS model,” is preferable and is currently available in REFPROP (NIST's software). In REFPROP's SuperTRAPP model, a scaling factor has been added so that, at infinite dilution of each compound, the viscosity of the pure fluid is exactly reproduced.⁸³ The definition and calculation details of each term in Eqn A.1 can be found elsewhere in the literature.⁸⁴

$$ 
\begin{align}
\Delta\mu(T, \varrho, x) &= \mu(T, \varrho, x) - \mu^*(T, x) \tag{A.1} \\
&= \Delta\mu(T_{\text{ref}}, \varrho_{\text{ref}})F_{\mu}(T, \varrho, x). \nonumber
\end{align}
$$

* the model proposed by Dean, Stiel, and Brokaw, hereafter referred to as DSB (Eqn A.2),⁸⁵ which is simple to implement, applies to high-pressure pure fluids and their mixtures and is available in Simulis Thermodynamics (ProSim software). Using a semitheoretical approach, Predvoditelev⁸⁶ showed that the residual viscosity Δμ(T, ρ, x) is simply a function of density. Subsequent works have been built on this finding, including those by Stiel and co-workers, and more recently by Galliero et al.⁸⁷ The Dean-Stiel version implemented in Simulis Thermodynamics uses the Brokaw model⁸⁸ to estimate the viscosity of gas mixtures at low-to-moderate pressure (or dilute-gas viscosity) μ*(T, x).

$$ 
\begin{align}
\Delta\mu(T, \varrho, x) &= \mu - \mu^*(T, x) \nonumber \\
&= \frac{10.8 \times 10^{-5}(e^{1.439\rho_R} - e^{-1.11\rho_R^{1.858}})}{\xi}, \tag{A.2}
\end{align}
$$

where ρ_R is a reduced density and ξ a viscosity parameter

$$ \xi = \frac{T_c^{1/6}}{M^{1/2}P_c^{2/3}} \tag{A.3} $$
industrially important fluids. Ind Eng Chem Res. 2022; 61 (42): 15449-72.

84. Chichester JC, Huber ML. Documentation and assessment of transport property model for mixtures implemented in NIST REFPROP (Version 8.0). US Department of Commerce, Technology Administration, National Institute of Standards and Technology; 2008.

85. Dean DE, Stiel LI. The viscosity of nonpolar gas mixtures at moderate and high pressures. AIChE J. 1965;11(3):526-32.

86. Predvoditelev A. Some invariant quantities in the heat-conduction and viscosity theory. Zh Fiz Khim. 1948; 22 (3): 339-48.

87. Galliéro G, Boned C, Baylaucq A. Molecular dynamics study of the Lennard-Jones fluid viscosity: application to real fluids. Ind Eng Chem Res. 2005; 44 (17): 6963-72.

88. Brokaw RS. Predicting transport properties of dilute gases. Ind Eng Chem Process Des Dev. 1969; 8 (2): 240-53.

# Appendix A: Models of density versus viscosity behavior for gas mixtures

In this appendix, we start by introducing models for the viscosity of gas mixtures that have a larger range of validity than the low-pressure (or low-density) approximation (Wilke's model) presented and used in the subsection 'Rate conditions for viscous fingering avoidance (front stability)'. Subsequently, we compare the results of those models for the gas mixtures of interest in this work, focusing on how viscosity and density are related since this relation determines the critical front velocity (see Eqn 5). Densities are calculated by the GERG-2008 model,⁴² which also has a larger range of validity than the low-pressure, IG approximation.

Those models are then used to assess the validity of the low-pressure (IG+Wilke) approximation for the gas mixtures and conditions of interest, with a focus on how density and viscosity are related in the H₂-rich region. Finally, mixtures of H₂ and CO₂ are examined over the whole interval of compositions.

## A.1. Viscosities of gas mixtures

A variety of methods have been developed to accurately estimate the fluid viscosities of pure components and their mixtures. These methods present advantages and limitations generally linked to their type, either “empirical” requiring a lot of experimental data or “predictive” with a solid theoretical background and requiring little or no experimental data. Given the very limited available high-pressure viscosity data especially for the H₂ + CO₂ system, predictive models are the

most appropriate in this case. The best known predictive models used to estimate the viscosities of high-pressure gas mixtures are

* the TRAPP (Transport Property Prediction) -type model (Eqn A.1) developed by Ely and Hanley⁸¹ based on the ECS theory and applicable to systems containing hydrocarbons, permanent gases, CO₂, H₂S, and refrigerants. However, the original version (TRAPP) is not recommended for use at reduced temperatures (T/T_c) above 0.925, and in these conditions, its variant called SuperTRAPP (see Millat et al.⁸²), which has been further refined and implemented under the generic name “ECS model,” is preferable and is currently available in REFPROP (NIST's software). In REFPROP's SuperTRAPP model, a scaling factor has been added so that, at infinite dilution of each compound, the viscosity of the pure fluid is exactly reproduced.⁸³ The definition and calculation details of each term in Eqn A.1 can be found elsewhere in the literature.⁸⁴

$$ \begin{aligned} \Delta\mu(T, \varrho, x) &= \mu(T, \varrho, x) - \mu^*(T, x) \\ &= \Delta\mu(T_{\text{ref}}, \varrho_{\text{ref}})F_{\mu}(T, \varrho, x). \end{aligned} \quad (A.1) $$

* the model proposed by Dean, Stiel, and Brokaw, hereafter referred to as DSB (Eqn A.2),⁸⁵ which is simple to implement, applies to high-pressure pure fluids and their mixtures and is available in Simulis Thermodynamics (ProSim software). Using a semitheoretical approach, Predvoditelev⁸⁶ showed that the residual viscosity Δμ(T, ρ, x) is simply a function of density. Subsequent works have been built on this finding, including those by Stiel and co-workers, and more recently by Galliero et al.⁸⁷ The Dean-Stiel version implemented in Simulis Thermodynamics uses the Brokaw model⁸⁸ to estimate the viscosity of gas mixtures at low-to-moderate pressure (or dilute-gas viscosity) μ*(T, x).

$$ \begin{aligned} \Delta\mu(T, \varrho, x) &= \mu - \mu^*(T, x) \\ &= \frac{10.8 \times 10^{-5}(e^{1.439\rho_R} - e^{-1.11\rho_R^{1.858}})}{\xi}, \end{aligned} \quad (A.2) $$

where ρ_R is a reduced density and ξ a viscosity parameter

$$ \xi = \frac{T_c^{1/6}}{M^{1/2}P_c^{2/3}} \quad (A.3) $$
with *T*<sub>c</sub>, *P*<sub>c</sub>, and *M* the molar-averaged critical temperatures, pressures, and molecular weights of the mixture.

## A.2. Density versus viscosity of H₂-rich mixtures

The densities and viscosities of mixtures of H₂ with either CH₄, N₂, or CO₂ have been calculated for the temperatures and pressures of interest by using the low-pressure IG and Wilke models (presented in the main text) and by the SuperTRAPP or DSB models for the viscosities (see the above) combined with the GERG-2008 model for densities. Some of those calculations are shown in Fig. A.1, where the density is plotted as a function of viscosity for the three mixtures at *T* = 298 K and *P* = 5 MPa. Density increases along with viscosity as the CO₂ content is increased in the mixture: only H₂-rich mixtures (*0*.< *y**< *1*) are considered in Fig. A.1, as the minimum of *dρ/dμ* over the whole interval of mixture compositions is reached for *y* → 1 (100% H₂). The density versus viscosity curves calculated with the three models are very close to each other and, in particular, have nearly identical slopes *dρ/dμ* where both densities and viscosities reach their minimum value (for *y* = 1). These slopes multiplied by the product *kg* of reservoir permeability by the acceleration of gravity are the critical velocities *U*<sub>m</sub> (see Eqn 5), which coincide with the analytical expression, Eqn 10, and in particular with the velocities given in Table 2 for a medium with *k* = 1 darcy. Interestingly, the critical velocity of an H₂/CO₂ front at *T* = 298 K and *P* = 5 MPa can be predicted by using the IG+Wilke approach.

## A.3. Mixtures of H₂ with CO₂

How viscosities and densities are related in H₂/CO₂ mixtures must be examined over the whole interval of compositions, including the CO₂-rich region where the IG and Wilke approximations are not appropriate. We have evaluated the densities and viscosities of those mixtures using four modeling approaches: the three approaches mentioned above (IG+Wilke, GERG-2008+SuperTRAPP, GERG-2008+DSB) and another modeling approach, labeled IG+Wilke (HZ), in which densities are calculated using the IG law and viscosities Wilke's model with the HZ exchange term (Eqn 9). For those four approaches, the densities of those mixtures are plotted as a function of viscosity in Fig. A.2 over the whole compositional interval (from

<figure><img src="image_6.png" /><figcaption>(A) H₂ + CH₄ at 298 K and 5 MPa</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>(B) H₂ + N₂ at 298 K and 5 MPa</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>(C) H₂ + CO₂ at 298 K and 5 MPa</figcaption></figure>

<caption>**Figure A.1.** Density versus viscosity of H₂-rich (y between 0.9 and 1) gas mixtures at T = 298 K, P = 5 MPa.</caption>

*y* = 0 to *y* = 1) for six temperature and pressure conditions: *T* = [298, 323, 373] K, *P* = 5 and 25 MPa. The few experimental data on CO₂-rich mixtures available, acquired very recently by Khosravi et al.⁴⁴, are represented in that figure as well.

The IG law underestimates the density of CO₂-rich mixtures, and the IG+Wilke models give an unphysical viscosity increase when a small amount of H₂ is added to CO₂ at 5 MPa. Clearly, the approach based on the
with *T*<sub>c</sub>, *P*<sub>c</sub>, and *M* the molar-averaged critical temperatures, pressures, and molecular weights of the mixture.

## A.2. Density versus viscosity of H₂-rich mixtures

The densities and viscosities of mixtures of H₂ with either CH₄, N₂, or CO₂ have been calculated for the temperatures and pressures of interest by using the low-pressure IG and Wilke models (presented in the main text) and by the SuperTRAPP or DSB models for the viscosities (see the above) combined with the GERG-2008 model for densities. Some of those calculations are shown in Fig. A.1, where the density is plotted as a function of viscosity for the three mixtures at *T* = 298 K and *P* = 5 MPa. Density increases along with viscosity as the CO₂ content is increased in the mixture: only H₂-rich mixtures (*0*.< *y**< *1*) are considered in Fig. A.1, as the minimum of *dρ/dμ* over the whole interval of mixture compositions is reached for *y* → 1 (100% H₂). The density versus viscosity curves calculated with the three models are very close to each other and, in particular, have nearly identical slopes *dρ/dμ* where both densities and viscosities reach their minimum value (for *y* = 1). These slopes multiplied by the product *kg* of reservoir permeability by the acceleration of gravity are the critical velocities *U*<sub>m</sub> (see Eqn 5), which coincide with the analytical expression, Eqn 10, and in particular with the velocities given in Table 2 for a medium with *k* = 1 darcy. Interestingly, the critical velocity of an H₂/CO₂ front at *T* = 298 K and *P* = 5 MPa can be predicted by using the IG+Wilke approach.

## A.3. Mixtures of H₂ with CO₂

How viscosities and densities are related in H₂/CO₂ mixtures must be examined over the whole interval of compositions, including the CO₂-rich region where the IG and Wilke approximations are not appropriate. We have evaluated the densities and viscosities of those mixtures using four modeling approaches: the three approaches mentioned above (IG+Wilke, GERG-2008+SuperTRAPP, GERG-2008+DSB) and another modeling approach, labeled IG+Wilke (HZ), in which densities are calculated using the IG law and viscosities Wilke's model with the HZ exchange term (Eqn 9). For those four approaches, the densities of those mixtures are plotted as a function of viscosity in Fig. A.2 over the whole compositional interval (from

<figure><img src="image_7.png" /><figcaption>(A) H₂ + CH₄ at 298 K and 5 MPa</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>(B) H₂ + N₂ at 298 K and 5 MPa</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>(C) H₂ + CO₂ at 298 K and 5 MPa</figcaption></figure>

<caption>**Figure A.1.** Density versus viscosity of H₂-rich (y between 0.9 and 1) gas mixtures at T = 298 K, P = 5 MPa.</caption>

*y* = 0 to *y* = 1) for six temperature and pressure conditions: *T* = [298, 323, 373] K, *P* = 5 and 25 MPa. The few experimental data on CO₂-rich mixtures available, acquired very recently by Khosravi et al.⁴⁴, are represented in that figure as well.

The IG law underestimates the density of CO₂-rich mixtures, and the IG+Wilke models give an unphysical viscosity increase when a small amount of H₂ is added to CO₂ at 5 MPa. Clearly, the approach based on the
<figure><img src="image_9.png" /><figcaption>Figure A.2. Calculated density versus viscosity of H₂ + CO₂ mixtures at *T* = [298, 323, 373] K, *P* = 5 and 25 MPa. The black squares at 5 MPa and 298 and 323 K are experimental measurements by Khosravi et al.⁴⁴.</figcaption></figure>

GERG-2008 and SuperTRAPP models is the most appropriate. This approach, together with the GERG-2008 + DSB approach, predicts densities and viscosities that differ strongly from the predictions by the IG and Wilke models, especially at high pressure (25 MPa). There is an inflection point for a composition *y* between 0 and 1 (see Fig. A.2), where the derivative

$dρ/dμ$ is maximum and therefore $dρ/dμ$ has two minima, one for *y* = 0 and one for *y* = 1. The derivatives $dρ/dμ$ over the whole interval of compositions are displayed in Fig. 2 and discussed in the subsection 'Rate conditions for viscous fingering avoidance (front stability)' of the main text for the four extreme conditions, 298 and 373 K, 5 and 25 MPa.
<figure><img src="image_10.png" /><figcaption>Figure A.2. Calculated density versus viscosity of H₂ + CO₂ mixtures at *T* = [298, 323, 373] K, *P* = 5 and 25 MPa. The black squares at 5 MPa and 298 and 323 K are experimental measurements by Khosravi et al.⁴⁴.</figcaption></figure>

GERG-2008 and SuperTRAPP models is the most appropriate. This approach, together with the GERG-2008 + DSB approach, predicts densities and viscosities that differ strongly from the predictions by the IG and Wilke models, especially at high pressure (25 MPa). There is an inflexion point for a composition *y* between 0 and 1 (see Fig. A.2), where the derivative

$dρ/dμ$ is maximum and therefore $dρ/dμ$ has two minima, one for *y* = 0 and one for *y* = 1. The derivatives $dρ/dμ$ over the whole interval of compositions are displayed in Fig. 2 and discussed in the subsection 'Rate conditions for viscous fingering avoidance (front stability)' of the main text for the four extreme conditions, 298 and 373 K, 5 and 25 MPa.
<figure><img src="image_10.png" /></figure>

**Sabrine Ben Rhouma**

Sabrine Ben Rhouma defended her PhD in December 2023 from the University of Pau and Pays de l'Adour and the University of Nice. Her research focuses on hydrogen storage in aquifers using CO₂ as cushion gas with a focus on computational thermodynamics and reservoir

engineering.

<figure><img src="image_11.png" /></figure>

**Daniel Broseta**

Daniel Broseta has been a professor at the University of Pau and Pays de l'Adour since 2005, following 15 years as a research engineer with IFPEN. He has been working on various aspects of CO₂ geological storage and, more recently, on gas hydrates.

<figure><img src="image_12.png" /></figure>

**Salaheddine Chabab**

Salaheddine Chabab obtained his PhD from Ecole des Mines ParisTech in 2020 on the thermodynamics of underground gas storage. He has developed thermophysical calculation software, published several scientific articles and book chapters, and presented his work at several

international conferences. He has benefited from the expertise of several renowned laboratories at Heriot Wat University, DTU, and ENSTA Paris. He is currently leading a 5-year partnership project (junior chair HYDR) at UPPA on the impacts of hydrogen storage on deep aquifers.
<figure><img src="image_11.png" /></figure>

**Sabrine Ben Rhouma**

Sabrine Ben Rhouma defended her PhD in December 2023 from the University of Pau and Pays de l'Adour and the University of Nice. Her research focuses on hydrogen storage in aquifers using CO₂ as cushion gas with a focus on computational thermodynamics and reservoir

engineering.

<figure><img src="image_12.png" /></figure>

**Daniel Broseta**

Daniel Broseta has been a professor at the University of Pau and Pays de l'Adour since 2005, following 15 years as a research engineer with IFPEN. He has been working on various aspects of CO₂ geological storage and, more recently, on gas hydrates.

<figure><img src="image_13.png" /></figure>

**Salaheddine Chabab**

Salaheddine Chabab obtained his PhD from Ecole des Mines ParisTech in 2020 on the thermodynamics of underground gas storage. He has developed thermophysical calculation software, published several scientific articles and book chapters, and presented his work at several

international conferences. He has benefited from the expertise of several renowned laboratories at Heriot Wat University, DTU, and ENSTA Paris. He is currently leading a 5-year partnership project (junior chair HYDR) at UPPA on the impacts of hydrogen storage on deep aquifers.