# Space and time distribution of subsurface H2 concentration in so-called “fairy circles”: Insight from a conceptual 2-D transport model

Andrey Myagkiy, Isabelle Moretti, Fabrice Brunet

► To cite this version:

Andrey Myagkiy, Isabelle Moretti, Fabrice Brunet. Space and time distribution of subsurface H2 concentration in so-called “fairy circles”: Insight from a conceptual 2-D transport model. Bulletin de la Société Géologique de France, 2020, 191 (1), 10.1051/bsgf/2020010. hal-02938377

HAL Id: hal-02938377

https://hal.science/hal-02938377

Submitted on 2 Mar 2021

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
# Space and time distribution of subsurface H₂ concentration in so-called “fairy circles”: Insight from a conceptual 2-D transport model

Andrey Myagkiy<fnref n="1" />,<fnref n="3" />, Isabelle Moretti<fnref n="2" />,<fnref n="4" />,* and Fabrice Brunet<fnref n="1" />

Received: 7 December 2019 / Accepted: 15 April 2020

**Abstract** – Natural H₂ emissions from the ground have now been measured in many places worldwide. These emissions can be localized on faults or be more diffuse in some sedimentary basins, usually of Proterozoic age. In such a case, emanation zones are often visible from aerial images or on high-resolution topographic maps since they correspond to slight depressions of circular to elliptic shape. Furthermore, the rounded depressions are covered with a scrubby vegetation which often contrasts with the surrounding vegetation. Although the emission structure displays a very regular shape, the distribution of H₂ concentration in the first meter of soil in such a structure does show a clear pattern. For example, the maximum concentration is almost never measured in the center of the structure and the few time-resolved data show that the soil H₂ concentration is variable with time. Here, the time and space evolution of H₂ concentration is simulated using a 2-D advective-diffusive model of H₂ transport in porous media. Several parameters have been tested as the depth and periodicity of the H₂ point source (pulsed), bacterial H₂ consumption and permeability heterogeneities of the soil. The radius of the structure is linked to the time spent by the H₂ in the soil that depends on the soil permeability, the depth of the gas leakage point and the pressure of the bubble. To account for field observations, the case of a shaly, less permeable, heterogeneity in the center of the structures has been modeled. It resulted in an increase of the concentration toward the rim of the structure and a close to zero signal in its center. If the deep signal is periodic with a frequency smaller than a few hours, H₂ concentration within the soil is almost constant; in other cases, the near surface concentration wave reflects the concentration periodicity of the source with a delay (in the range of 12 h for 30 m of soil) and so the near surface H₂ concentration values will be highly dependent on the time at which the measurement is performed. H₂ monitoring through a sensor network is thus mandatory to characterize the H₂ dynamics in the soil of fairy circles.

**Keywords:** natural H₂ / H₂ emission / soil micro-seeps / heterogeneous soil permeability

**Résumé** – Variabilité dans le temps et l'espace dans le sol des structures émettrices dits « ronds de sorcières »: apport d'un modèle conceptuel de transport 2D. Des émissions de H₂ du sol ont maintenant été notées à de nombreux endroits, elles peuvent être localisées sur des failles ou plus diffuses dans certains bassins sédimentaires, généralement Protérozoïques. Dans ce dernier cas, les émanations sont souvent visibles à partir d'images satellites ou de cartes topographiques à haute résolution puisqu'elles correspondent à de légères dépressions ; la végétation est affectée, souvent détruite et la forme de la structure est circulaire ou elliptique. Les mesures de la concentration de H₂ dans le sol, lorsque le dihydrogène s'échappe des bassins sédimentaires, semblent aléatoires. Le maximum n'est presque jamais localisé au centre de la structure, les quelques données continues en temps montrent que la concentration en H₂ du sol y est variable. Sur la base d'une modélisation du transport de H₂ dans le sol, nous proposons une explication à ces données, les caractéristiques de cette concentration qui pourrait paraître aléatoire sont principalement dues à l'émission non constante de H₂, la consommation de H₂ dans le sol et aux hétérogénéités du sol.

<footnote n="1">¹ Univ. Grenoble Alpes, Univ. Savoie Mont Blanc, CNRS, IRD, IFSTTAR, ISTerre, Grenoble, France</footnote>

<footnote n="2">² Engie SA, Courbevoie, France</footnote>

<footnote n="3">³ Storengy, Bois-Colombes, France</footnote>

<footnote n="4">⁴ E2S UPPA, Pau, France</footnote>

<footnote>*Corresponding author: isabelle.moretti@univ-pau.fr</footnote>
La taille de la structure émettrice est liée au temps passé par le l'hydrogène dans le sol, qui est influencé à la fois par la perméabilité du sol, par la profondeur de la fuite initiale de gaz et par la pression de la bulle. L'hétérogénéité des perméabilités des sols y entraîne une distribution non homogène de la concentration d'H₂. Les parties particulièrement moins perméables situées au centre des structures entraînent une augmentation de la concentration vers les limites de la structure et un signal proche de zéro au centre, comme observé sur le terrain. Si le signal source est régulier avec une fréquence horaire inférieure à quelques heures, un signal quasi constant peut être atteint. Dans les autres cas, le signal proche de la surface reflète le signal source en profondeur avec un retard (d'environ 12 h pour 30 m de sol) et les valeurs proches de la surface dépendent donc beaucoup de l'heure des mesures. Une surveillance permanente, de même que la connaissance de la perméabilité du sol, sont indispensables pour interpréter les données proches de la surface en termes de débit.

**Mots clés** : hydrogène natif / concentration d'H₂ dans les sols / émanations d'H₂ / hétérogénéité des sols

# 1 Introduction

## 1.1 H₂ transfer in subsurface

Subsurface H₂ gas concentrations far above atmospheric level and subsequent H₂ flow towards the surface have been reported in an increasing number of publications. H₂ gas concentrations at around 1 m below the surface have been measured in Russia and USA among others (Larin et al., 2014; Zgonnik et al., 2015) and even monitored in Brazil (Prinzhofer et al., 2019) on a time-resolved basis. The source of this hydrogen is still debated since different H₂-forming reactions and processes are likely to produce H₂ in the subsurface (Guélard et al., 2017). The hydrothermal oxidation by the seawater of ultramafic rocks exhumed along the mid oceanic ridges is one of the best-known redox process that produces H₂ under abiotic conditions (e.g., Brunet, 2019 and references therein). H₂ fluxes related to this phenomenon, called serpentinization, have been observed along the mid-Atlantic ridge Charlou et al. (2002, 2010) or the Pacific one (Welhan and Craig, 1979) where it takes place at high temperatures under hydrothermal conditions. H₂ emissions are also observed in convergence area where the oceanic lithosphere outcrops onshore such as in Oman, Philippines or New Caledonia (Neal and Stanger, 1983; Abrajano et al., 1988; Sano et al., 1993; Monnin et al., 2014; Deville and Prinzhofer, 2016; Vacquand et al., 2018). The same type of process based on the interaction between ultramafic rocks and water may occur at lower temperature and the fluid could be meteoric water (e.g., Mayhew et al., 2013; Miller et al., 2017).

Similar reactions of oxidation/hydration may also take place in continental settings if, for instance, Fe²⁺ rich rocks are present, as it has been suggested as a source for the continental H₂ accumulations in Mali (Prinzhofer et al., 2018). Radiolytic H₂ production in the Precambrian continental crust, which accounts for 72% of the continents, is another possibility (Sherwood et al., 2014) with a contribution which, for these last authors, could be comparable to the one of the oceanic domains, around 10¹¹ moles per year. In addition, H₂ may originate from mantle degassing (Larin et al., 2010) since H₂ is the initial and major component of the solar system. Mechanoradical H₂ generation along fault planes when they are creeping, or active, have also been observed (Klein et al., 2020); following these authors it can occur wherever there is comminution of silicate-bearing rocks. A frictional process has been proposed to explain the change in H₂/CH₄ ratio in a gas

field in the Jura front belt (Deronzier and Giouse, 2020) however its overall role for large scale H₂ production remains unknown. Bacterial activity is also prone to produce H₂ in specific environments (Conrad and Seiler, 1980), however this alternative has also not yet been quantified at the global scale.

Without attempting to discriminate between these possible sources, the present work focuses on the H₂ transfer in continental settings characterized by a soil cover and which leads to the formation of fairy circles. In the absence of soil cover, H₂ emanations are usually reported along faults. In the Oman ophiolite, for example, H₂-rich gases have the same exit points as alkaline spring waters, which are located along fault and shear zones (Neal and Stanger, 1983). In such a case, water and gases react immediately with the atmosphere resulting in carbonate deposits such as in Oman (Noel et al., 2018) and in New Caledonia (Deville and Prinzhofer, 2016). In case of water accumulation (river or lake), H₂ being poorly soluble in water at low pressure (Lopez-Lazaro et al., 2019 and reference therein), continuous degassing can be observed in the form of bubbles. They show a continuous but not permanent flow of H₂.

After a synthesis of the characteristics of H₂ emitting structures in terms of shape and H₂ concentrations in the soils, we explore, in the current study, through two-dimensional reaction-transport modeling, the situation where the H₂ migration toward the surface ends in a soil cover. Our model of H₂ subsurface distribution aims at accounting for the H₂ concentrations described as non-constant neither in time nor space reported worldwide in the circular-like shape structures.

## 1.2 Shape of the H₂ emitting structures in the presence of a soil layer

Due to their circular shape, the structures from which H₂ is leaking out to the atmosphere are often called “fairy circle”. We will use this wording even if the structures are not perfectly circular, and even sometimes perfectly elliptic as in South Carolina (Zgonnik et al., 2015). Figure 1 shows some examples of these structures.

Larin et al. (2014) have investigated a large number of them in the Voronezh Oblast, about 600 km south-east of Moscow. The size distribution deduced from the 562 measured circles shows variability from 3 km down to 100 m (Fig. 2 in Larin et al., 2014). Small structures are more frequent, and the authors noticed that above 1.5 km there is no more relationship
<figure><img src="image_2.png" /><figcaption>**Fig. 1.** Synthesis of the fairy circle emissions characteristics. A. A fairy circle in Brazil showing how the vegetation is affected. B. Variability on time over 5 weeks for a given sensor in the structure monitored in Brazil. See Prinzhofer et al. (2019) and Myagkiy et al. (2020) for details about that case. C. Daily variations for 6 sensors installed on the structure monitored in Brazil showing the 24 h period. D. Variability in size, of the emitting structures in USA (modified from Prinzhofer and Deville, 2015). E. Two profiles across a structure showing the spatial variability of H₂ concentration value between the center and the rim of the structure, Lake Podowoye, Russia. These values have not been taken at the same time but over 2 days. Modified from Larin et al. (2014). F. Profile showing an example of [H₂] concentration variation versus depth in USA (modified from Zgonnik et al., 2015).</figcaption></figure>

between size and frequency. The authors noted that if the center of the structures corresponds sometimes to a true depression, up to 5 m deep, often there is no real depression but just changes in the vegetation. When the center corresponds to a depression, wetlands are present. Small lakes and changes in the soil chemistry are noted, Larin et al. (2014) noticed that soils are sandier at the circle rim and more shaly in the center. Similarly, in South Carolina (USA), hundreds of depressions have been also reported, first by researchers working on wetlands and biodiversity (Semlitsch, 2000) and then by those studying H₂ emissions (Zgonnik et al., 2015). They noted, as in Russia, a large spectrum of sizes from 8 km down to 100 m for the major axis of the depressions which, in that case, are ellipses (Fig. 1D).

Larin et al. (2014) reported a structural alignment and suggested that the fairy circles are located above faults from where the H₂ leaks out. The crucial role of faults for fluid (liquid as gas) migration is well established and will not be discussed here. The general view is that pervasive circulation of the fluids in the geological medium is large when bulk porosity and permeability are large – that is the case in all carrier beds and aquifers – whereas as soon as the permeability is low (shale, metamorphic rock, granite ...) the fluids circulate mainly in the fractures and faults. The cross-section (Fig. 3 in Larin et al., 2014) shows that, in their case study, the Precambrian basement is located between 0 and 500 m depth. The fact that H₂ escapes from and along fractures from the basement is thus a possibility, but when reaching the overlying soil cover, upward transport will proceed within the soil medium which has usually a relatively “good” porosity and

permeability that are at least few orders of magnitude higher than that of the Precambrian basement.

When there is no soil cover, as in Oman, New Caledonia (Deville and Prinzhofer, 2016; Vacquand et al., 2018) but also Turkey in Chimeres area or in Philippines, the gas mainly seeps from fractures, these settings will not be discussed here. Roughly speaking the shape of these seeps is the shape of the fracture network and is also out of the scope of this study.

## 1.3 Time-resolved H₂-concentration measurements within the fairy circles

Up to now, even though the list of basins where H₂ emanations have been reported is getting longer (Zgonnik, 2020), both time-resolved and spatial distribution measurements of H₂ concentration in soils from fairy circles are still rare. However, the datasets cited in the previous section all show that the soil H₂-concentration is variable across the structure and that the maximum concentration is barely, if not never, in the center of the structure (Figs. 5 and 6, Larin et al., 2014; Prinzhofer et al., 2019).

The only published permanent monitoring of an H₂ concentration in the soil porosity of an emitting structure has been performed in the São Francisco basin in Brazil (Prinzhofer et al., 2019). A few other examples of published permanent H₂ monitoring data are either related to seismic risk and were carried out in the vicinity of active faults, as in San Andreas (Sato et al., 1986) and in Japan or for safety control during mining operations (Kola Peninsula, Russia, Firstov and Shirokov, 2005; Nivin et al., 2016).
<figure><img src="image_3.png" /><figcaption>**Fig. 2.** Continuous leakage of H₂ with the source at 30 m and without active bacteria in the soil. Soil gas permeability is 10 Da. A. Snapshots of an H₂ propagation. The color from blue to red indicates the [H₂] of the propagating gas normalized to the source [H₂]. B. [H₂] distribution as a function of time at 80-cm depth.</figcaption></figure>

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

<figure><img src="image_5.png" /><figcaption>**Fig. 3.** Continuous leakage of H₂ with the source at 30 m, considering the active H₂-consuming bacteria in the soil. Soil gas permeability is 10 Da. A. Snapshots of an H₂ propagation. B. [H₂] distribution as a function of time at 80-cm depth.</figcaption></figure>

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

In Brazil, the monitored circular structure is a depression with a diameter close to 500 m. The published H₂-concentration dataset covers more or less one and half month of record. It shows a daily wavelength consistent with a pulsed source of a still unclear origin. Almost all the sensors which are near the rim of the structure and installed at 80 cm depth, show a maximum at around noon and are recording H₂ flow for about

6 h a day. H₂ soil concentrations are measured each hour and were around 150 ppm during the week of August 2019 that is reported in Prinzhofer et al. (2019). Another sensor recorded higher values (up to 1200 ppm), it showed also a 24 h period but H₂ concentration maximum was reached during the night. A sensor installed in the center of the structure for few months recorded only a few H₂ concentration pulses but most of the
<figure><img src="image_7.png" /><figcaption>**Fig. 4.** Pulse emission of H₂ from the source at 30-m depth. A. Snapshots of an H₂ gas pulse propagation in time. B. [H₂%] distribution along a vertical line passing through the source point (see white line on the *t*=0 snapshot of Fig. 4A). Note that the [H₂] measured at a given depth will strongly vary over the pulse propagation cycle. C. [H₂] distribution at 80-cm depth as a function of time. Reference case: no active bacteria, K = 10 Da, surplus in gas pressure at source is 2 bar, Depth of source = 30 m.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**Fig. 5.** Maximum [H₂] envelop at 80-cm depth as a function of source pressure (0.2, 0.5 and 1 MPa) for a 30-meter source (A) without active bacteria and (B) with active bacteria.</figcaption></figure>

time had a signal close to 0 (Prinzhofer and Engie's Team personal communication). The 24 h period is not yet fully understood (Myagkiy et al., 2020) nor the gap, or delay, between the highly emitting area and the other ones. However, the fact that the H₂ concentration in the soil, near the surface, is changing with time along the day is shown by all data.

The effect, on H₂ diffusive transport, of bacterial activity has already been investigated on soils from the same structure of the São Francisco basin (Myagkiy et al., 2020). The kinetics of H₂ consumption by the bacteria present in the soil, which was calibrated in laboratory experiments, is proved to be very large and may reduce or even completely hide the deeper H₂
flows. Basically, the consumption of H₂ by bacteria increases with the residence time of H₂ in the soil whatever the transport mode, advective or diffusive.

There are only few published data of H₂ in the soil at various depths. In Oman, Zgonnik et al. (2015) recorded concentrations at 3 different depths of 30, 60 and 110 cm and found concentrations of 150, 375 and 3400 ppm respectively (S54, f1, 2 and 3, Tab. A1 in the cited paper). However, in the same article they report, also data showing an inverse gradient with H₂ concentration which decreases with depth (sample 51 f1 and f2). These data correspond to one-shot measurements without indication of sampling time so even if these variations of concentration with depth are real, they could merely be a snapshot of a situation that may evolve with time.

The monitoring carried out in Brazil has already shown that recording data versus time for a least few days is crucial not only to compute a flow but also even just to know whether H₂ emission do take place. A single zero H₂ concentration (that we will note [H₂] further in this paper) value at a given time cannot be interpreted as a lack of H₂ flow in a dynamic system. Based on a numerical model of H₂ transmission in the soil we will evaluate here which are the key parameters that affect the near surface signal and can account for (1) the size of these structures; (2) the presence of a ring with higher H₂ concentration; (3) the variation through time of the near surface signal. We will therefore discuss the causes of the apparent randomness of the H₂ concentration recorded near the surface.

# 2 Transport of H₂ across a soil

## 2.1 Physical governing equations and assumptions

A 2-D numerical advective-diffusive model of gas reactive transport in porous media was developed. For that purpose, the numerical simulator COMSOL Multiphysics v 5.2 was used. All the details concerning the equations and the parameters used in our computations are presented in Appendix A.

Bacterial activity, that was found to have a prominent role in the consumption, or production of gas in the upper part of the subsurface, was incorporated in calculations through the Michaelis–Menten kinetic law of H₂ consumption. The parameters, experimentally derived for the soils of Brazilian H₂-emitting structure, have been used in simulations. Two types of calculation were performed, namely with and without H₂ bacterial consumption, for a comparison. The details of its measurements can be found in Myagkiy et al. (2020). The model assumes a constant gradient in bacteria distribution; the maximum value of bacterial consumption activity is assumed to appear at the surface and reaches zero at the bottom of the model.

It should be noted that the numerical model was resolved in the range of [H₂] at which microbiological uptake kinetics was measured, so up to few thousandths of ppm. For the sake of simplicity, modelled H₂ concentrations in soil were expressed as percentage of the source concentration. The impact of concentration at source in the range of applicability of microbiological consumption law was found to be negligible compared to the other parameters. The [H₂] at the surface was set to zero in the model as the content of dihydrogen in atmosphere is negligible compared to the modeled ones.

The impact of water saturation degree in the soil is not regarded in this study but can have an important impact on hydrogen characteristic transport time in porous media (see e. g. Myagkiy et al., 2020). High water content in the soil reduces its effective porosity and restricts the diffusion and advection of H₂ into the soil. For the shallow subsurface environment up to 30 m down, as in our model, the pressure and temperature are too low to warrant a significant dissolution of hydrogen in water. The dihydrogen arising from deep source is most likely to exceed the limit of its water solubility at these conditions which may lead to formation of the bubbles. That can be one of the ways of transport of dihydrogen through an aquifer. This mode of transport may be described by formation of gas clusters (e.g. Li and Yortsos, 1995; Dominguez et al., 2000) and is out of the scope of this study, which highlights the transport in the vadose zone.

## 2.2 Boundary conditions

**Soil permeability:** the average permeability used in most of the simulations reported here is 10 Da (or 10⁻¹¹ m²), appropriate for a moderately permeable soil (Auer et al., 1996; Lazik et al., 2016), the influence of heterogeneity in this anisotropy and especially of lower permeabilities in clayey area will be discussed.

**Punctual deep leakage:** The modelling case involves a deep origin H₂ source; the gas is, therefore, supposed to vent upward into the overlying rocks and then undergoes gas-soil interactions in the soil layer. The source is represented by a source point at the bottom of the modelling domain, where the gas injects due to overpressure and the dihydrogen enters the porous soil matrix. This periodic concentration pulse of H₂ is obviously an ideal case. However, in the absence of constraints on the emitting zone itself, more complex emission type would be difficult to justify. The Precambrian or ophiolitic rocks, which are the more likely source rock for H₂, being both poorly permeable we will consider that the H₂ emission happens at a given point in the subsurface which corresponds to a fracture that allows, when the fluid pressure is sufficient the escape of the gas. Like a periodic geyser from which gases escape in geothermal area. Other geological contexts may lead to this kind of punctual leakage: the degassing at the top of an aquifer since the H₂ solubility is strongly decreasing with the pressure; leakage by natural hydraulic fracturing is also possible, such a microcrack may get open when local stress tensor changes or at the top of an accumulation if the seal in not good enough as observed above the HC gas accumulation (Sibson, 2003; Petrie, 2014), capillary seal break is another possibility as already proposed for methane leakage (Cartwright and Santamarina, 2015).

**Periodicity of the deep H₂ leakage:** Since all the available monitoring data show that the H₂ flow is periodic we will mainly use that hypothesis. H₂ monitoring data from Brazil and Kola Peninsula suggest a near 24 h period and a 6 h long flow. We will use these values in our model for the pulsed source. In addition, shorter pulse periods and even continuous H₂ leakage will be modeled for comparison.

**Amount of the H₂ leaking:** one can suspect that the size of the depression is related to the amount and/or pressure of dihydrogen that is degassing during one cycle. Therefore, the
<figure><img src="image_9.png" /><figcaption>**Fig. 6.** A. Snapshots of an H₂-bearing gas pulse propagation as for Figure 4 but with the introduction of a low permeability zone (0.01 Da) simulating a clay-rich soil heterogeneity at the center of the structure. Active bacteria are considered here. B. Maximum [H₂] envelop at 80-cm depth.</figcaption></figure>

different $\Delta P_{\text{fracture}}$ will be tested and the size of the area influenced by the H₂ flow compared to the fairly circles sizes.

**Permeability and porosity of the soils:** the transport of fluid in soil is highly dependent of its hydraulic characteristics. In the soil, as in the rock, values may be highly variable, of at least three orders of magnitude, from 10 Da (for gravel) to 0.01 Da (to shaly soil). The influence of this parameter and of its distribution within the soil of a same structure will be tested.

**Bacterial activity:** the impact of this parameter has already been discussed in a previous paper (Myagkiy et al., 2020), so we will incorporate this effect in the discussion but the impact of the variation of the quantity of bacteria and kinetic of the H₂ consumption will not be discussed here in detail.

# 3 Results and discussion

## 3.1 Permanent flow of H₂

Figure 2 shows the concentration in the soil above a permanent leakage of H₂ located at 30 m. After a few hours, H₂ reaches the surface and the horizontal concentration profile has a bell shape with a maximum located just above the leak point. Then, the H₂ concentration profile spreads laterally. Therefore, the size of the area affected by the H₂ flow increases with time. Whereas it is about 200 m large after around 2 months, it reaches 400 m after a bit more than a year.

One may note that:

- The H₂ soil concentration has a “bell shape” distribution with its maximum located at the center of the structure

- The requested time for H₂ to reach the surface is around 12 h (number which is obviously depending on the depth of the leakage and pressure of dihydrogen recharge from the deepest level). Above the leak point after 2 days, the [H₂]

remains constant in intensity over time. In the center, all the values are so the same, the maximum that could be reached, i.e. the same concentration as at the leakage point

- Progressively, the concentration of the soil is increasing laterally and the affected area, the fairy circle, is getting larger and larger through time. The process of the growing of the structure, is getting slower but continues, the size of the structure, its diameter, is about 50 m after 1 day, 200 m after 50 days and 400 m after 500 days.

Clearly, none of the characteristics of this H₂ signal fits with the observed data described above.

A second set of runs has been performed adding the effect of H₂ consumption in soil due to bacterial activity. The results are shown in Figure 3. In such a case, the kinetics of the consumption bound the growing of the fairy circle. H₂ consumption limits the radial expansion of the H₂ emitting zone. With bacterial activity, after few days the permanent state is achieved, there is an equilibrium between the H₂ emission and consumption, and the size of the structure has reached its maximum which is about 100 m of diameter for this 30-m deep source. The concentration along profile remains thus constant. The balance between the transport time of H₂ in the soil to reach the surface and the kinetics of the H₂ consumption is linked to soil characteristic and bacterial activity; it is absolutely not a one-to-one relationship.

In conclusion, none of these results fit with the literature data, they do not explain the ring of high [H₂] values neither the variabilities of the [H₂] versus time. In addition, one may note that without consumption of H₂ in the soil the distribution of H₂ is continuing and growth of the structure continues over years. In the literature, the soil H₂ concentration is always variable but nobody noticed so far large changes in the
<figure><img src="image_10.png" /><figcaption>Fig. 7. Concentration at the center and on the border of the structure in case of heterogeneity at the center with K = 0.01 Da and active bacteria. A. Emission every 24 h. B. Emission every 6 h.</figcaption></figure>

structure size through time. Whatever its origin, bacterial or
not, H₂ consumption in the soils seems mandatory to explain
this stability in case of continuous flow.

## 3.2 H₂ transport reference case pulsed source

A third set of runs has been performed using a pulsed
source at a given depth of 30 m, with bacterial activity, the
influence of the source depth will be discussed later on. The
H₂ pulse chosen as boundary condition at the leakage point is
12 h long and shows a maximum at *t*=6 h. Appendix A that
resumes the equations and hypotheses shows the concentration
of H₂ at the source point for the various runs. Figure 4A shows
the rise of the hydrogen from the source, as previously the
surface is reached above the leaking point after 12 h. Then, a
ring develops defining a high [H₂] area which gets slightly
larger until it reaches the maximum diameter of 100 m in the
case of a leak point situated at 30 m. Figure 4B shows the same
results but highlights the profile versus depth at different times.
Interestingly, one may observe that according to the time of the
day, the H₂ concentration gradient, $\frac{\partial[\text{H}_2]}{\partial z}$ at the near surface can
be negative (t0 + 18 h) or positive (6 h and 12 h) depending on
the measurement time.

Figure 4C presents the near surface emanations across the
structure at 4 different time-frames, the values are given at
80 cm depth in the center of the structure since it is roughly
where the monitoring or the measurement are done in the
literature. One may observe the increase of concentration at the
center of the structure until 12 h. Then, the peak splits up. The
ring forms between 12 and 18 h and then it slightly grows. It is
worth to notice that at 24 h (dashed red line), the ring (the two
picks in the profile) is perfectly shaped, the [H₂] in the center is
zero.

The results of this modeling, even very simple, are a key to
understand literature data where apparently random [H₂]
values are reported through the structures as summarized in the
introduction section. The fact that the source is pulsed induces

the ring formation and the inversion of the maximum between
the center of the structure and the rim. The delay, i.e. the travel
time, between the pulse maximum at the bottom of the studied
zone and the moment it reaches the surface explains these
versus time and versus depth H₂ concentration curves. The
mode of transport, advective and diffusive, and the soil
characteristics, essentially its permeability, have a strong
influence on the velocity of hydrogen transport as it will be
discussed later on.

## 3.3 Size of the structures

The first run with a continuous flow without active bacteria
or any other trapping mechanisms already shows that if there is
no consumption of H₂ in the soil, H₂ bearing area could be very
large, a few hundreds of meters even for a near surface leakage.
The scavenging of the H₂ in the soil (e.g., by bacteria) will limit
the size of the structure for a given point source flux and depth.
The H₂ travel distance which can continuously grow in the
absence of H₂ consumption becomes limited when H₂ is
consumed on its way to the surface. Accordingly, a deeper
source will result in smaller structures or even in no structure at
all if all the emitted H₂ is consumed on its way to the surface.
However, the literature data in Russia (Larin et al., 2014) as in
the USA (Zgonnik et al., 2015) show that in the same area, the
size of the structures can be variable. Considering the fact that
the thickness of the soil does not usually change too quickly,
alternative explanations have to be looked for to explain the
variability of the structure diameters. We, therefore, tested the
pressure of the H₂ “bubble”, which is prone to vary. Since the
time dependency of the pulse is the same (Fig. A1B), that
increase of pressure of the H₂ bubble means also than a higher
quantity of H₂ is released into the soil from the leak point.

The source being pulsed, the [H₂] signal measured at the
surface is time dependent. This time dependency has been
described previously so we will focus our analysis on the size
of the structure. Therefore, instead of a succession of
snapshots, the envelop of oscillating [H₂] will be displayed (Fig. 5), which corresponds to the maximum [H₂] measured on each point of the surface within the first 30 h (with and without bacterial activity).

With the source at 30 m and without bacteria, the circle is 100 m large after 1 day if the pressure is low, 2 bars, and can reach a diameter of 140 m if the H₂ is overpressured (10 bars) when escaping from the basement. These values are reduced to 80 m and 120 m (respectively) in the case of bacterial activity consuming H₂ in the subsurface.

The numerical relationship between size of the structure and depth of the leak point is depending also on the soil permeability and porosity that will be discussed later. However, the prominent observation is the link between the characteristics of the H₂ pulse and the size of the structures. In a given basin, for given soil characteristics and thickness, the size of the structures may vary if the H₂ pulse is more or less pressurized, that means also that the size increases if the quantity of the H₂ released at each pulse is getting larger. It seems particularly interesting to explain the occurrence of small structures within the large structures as observed in Carolina bay.

The fact that H₂ consumption in the soil has a strong impact on the size of the structure by limiting its growth has already been shown with the first runs with permanent leakage. The possibility of inverting the surface signal to infer the depth of the leakage point is therefore illusory, too many unknown, or poorly known, characteristics are impacting the size of the emitting zone. In addition to that, we only model the last step of the transport of the gas toward the surface, when the gas invaded the unconsolidated porous media, as evidenced on high resolution seismic data gas chimney (along fault or not) may allow efficient upward migration of gas from deeper level (Gay et al., 2006; Maia et al., 2016). The knowledge of the structural pattern and fault network and so seismic data are mandatory to see where the chimneys root.

## 3.4 Heterogeneity of the soil

Up to now we have considered a homogeneous soil of constant permeability and porosity in our model; practically it is almost never the case. The circle corresponding to slight depressions, fine particles are transported towards the center by rainwater where they accumulate to form a shaly deposit. It is what has been observed in Russia (Larin et al., 2010) as in Brazil (Prinzhofer et al., 2019). Shaly soils have a very different permeability than silty ones, up to 3 orders of magnitude (e.g. Terzaghi and Peck, 1968; Freeze and Cherry, 1979 in Tab. 2.2, p. 29). To test the influence of permeability heterogeneity on the soil H₂ concentration, we added to our model a low permeability shaly layer at the center, above the emitting point. The computation was performed considering that the low-permeability layer is 3 m thick and 50 m wide and has a permeability of 0.01 Da, i.e. 3 orders of magnitude lower than the surrounding soil, and two times smaller in porosity. In this run the H₂ consumption by bacterial activity is considered.

Figure 6A shows a snapshot of the resulting H₂ flow. In that case, H₂ is first moving upward and laterally as previously but when H₂ flow reaches the low permeability area, it gets around it taking the fastest way to the surface. The ring of high H₂ values already present due to temporal effect is now reinforced

at the border of the less permeable area. The maximum is no more at the center of the structure but shifted at the rim of the low permeability layer. This can be clearly seen in Figure 6B that represents the maximum envelop of [H₂] near surface at 80 cm depth, which corresponds to the maximum [H₂] measured at each point of this depth within the first 30 h.

The amplitudes of the signal versus time at the center of the structure and at the border of the shaly layer are shown in Figure 7A, at the same depth as in Figure 6B. It highlights the strong impact of the permeability on the time of H₂ transport. The maximum is reached near the surface on the border after almost 1 day, that means 18 h after the peak of the emission, when it reaches the center after 30 h so 24 h after the peak of the emission. These 6 h of “delay” results in a higher H₂ consumption in the soil and so the signal in the center is close to zero.

It is also worth to note that with bacterial activity 60% of the H₂ is consumed before it reaches the surface but the “travel time” is only slightly affected by the bacteria. The consumption of H₂ affects the H₂ gradient, when it is strongly dependent on the permeability. In that case, after 2 days the signal is stable, the deep H₂ flow reaches the surface after a few hours, up to almost one and a half day for the center or far from the shaly part and then flow and consumption neutralize each other. A steady state is reached.

## 3.5 Frequency of the pulse

As already mentioned, a near 24 h long periodicity has been observed within the few published monitoring data, in Russia, in mines, as well as in Brazil in a grass field. The explanation for that period is still debated but the modeling presented here shows that a continuous flow is not compatible with the field data. In addition, we tested the hypothesis of a higher pulse frequency with 4 pulses per days instead of one. In order to compare the results with the previous case all the other parameters are kept the same, the leak point is 30 m deep and there is a low-permeability layer, 50 m large and 3 m thick, in the center of modeled area. Figure 7B shows the [H₂] in the soil in the middle of the structure, so in the shaly part, and at the rim of this shaly part. As in the 24h hours frequency case, Figure 7A, in Figure 7B the H₂ front first reaches the rim, after about 36 h, and then the positions of maximum and minimum are noticed but they are poorly marked. Globally, after 2 days, there is an almost constant signal. It is even more pronounced in the center of the structure, the travel time is larger since the soil permeability is smaller but after 3 days the signal reaches a kind of steady state and the signal is near constant. The amount of H₂ in the soil is globally smaller since the amount of H₂ released at the source is globally smaller. We do not know the exact quantity of H₂ that is leaking in the various basins so we will not discuss so much the final concentrations. However it is worth to realize that even if the soil is just 30 m thick a high frequency pulse will result in almost a permanent signal and oscillation will not be detected by the near surface measurement. After a maximum of 3 days, the H₂ concentration is stable as for a constant source. However, a constant [H₂] at given location in time and with depth does not fit with the observed data.
## 3.6 Other potential processes that may influence the H₂ concentration in the soil

As explained previously we modeled the transport of H₂ near the surface in an undersaturated high permeability soil. The chosen approach has limitations and doesn't consider all the physical processes that may exist in the subsurface, but also in the atmosphere, and affect the H₂ flow.

HC gas seeps have been mainly studied offshore and are the marks of the presence of gas in the subsurface; the related pockmarks and/or mud volcano are very common and rather easily imaged with the high-resolution bathymetry and 3D seismic data now available (e.g. Maia et al., 2016 and reference inside). The images clearly show that the gas escape through narrow conduits from the reservoir, which very often, offshore, is the free gas zone below the hydrate layer. These conduits are vertical and may in some cases correspond to small faults even in young, and so uncompacted, sediments. Except onshore near the surface above the piezometric level, all sediments are water-saturated however seismic data as experiments show that the gas may invade unconsolidated porous media as soon as the overpressure overpass the effective stress (Fauria and Rempel, 2011). The hydrogen as the methane is poorly soluble in water at low temperature and low pressure (Lopez-Lazaro et al., 2019). Therefore, we simplify the model and consider only the free gas, not the dissolved one. We also considered that there is enough gas in this area of gas leakage to model it as one single phase and not as isolated bubbles in a water column. Simulation of individual bubbles, their migration and coalescence, could be done at pore scale (e.g. Mahabadi et al., 2018) but at basin scale it is obviously more challenging, and it is not the level of detail that we target in this study. In term of H₂ consumption in the soil, our model is also simple, we only take into account the near surface bacterial activity and the distribution of the bacteria versus depth is an exponential decrease. The H₂ consumption by inorganic processes is not modeled, for instance adsorption of H₂ in the shaly soil in the center of the structure could also amplify the H₂ signal attenuation phenomenon in the center of the structures but is not computed here.

External causes of H₂ flow variation are also neglected, in our computation the air pressure is set constant over the simulation. Effect of barometric pumping (Auer et al., 1996) is not considered on the advective H₂ transport in the soil. Water content in the soil, due to for instance the seasonality of the rainfall may also change the soil parameters and increase the complexity of the signal.

Despite of all these limitations, of which we are aware, some conclusions can be drawn.

## 4 Conclusions

The present 2-D model of H₂ transport (advection and diffusion) in a soil layer from a leak point of periodic [H₂] to the surface offers a simple tool to interpret features that are encountered in the field where active emissions of dihydrogen take place. Several parameters may influence the final concentration and the time and depth dependency of this H₂. The goal of this paper is not to adjust a set of parameters to perfectly fit to the data on a given case but rather to understand

what controls the 2D distribution and amplitude of hydrogen concentration.

- Observed shapes: fairy circles. Our modeling shows that a punctual hydrogen leakage at a given depth results in a circular structure at the surface. The fact that currently detected emitting structures are not growing and have almost stable shape over few years of observations is incompatible with a permanent flow over years and also strongly suggests the consumption of H₂ in the soil.

- A pulsed source allows to explain the spatial variability of the [H₂] in the soil as observed in Brazil (Prinzhofer et al., 2019) among other localities. The presented model is 2D, but one may note that the fact that the emitting structures are often aligned above accidents as in Russia (Larin et al., 2010) are fully compatible with a model where the H₂ is moving in the basement fractures, escape from the low permeability rock through pumped fracture opening along the main faults and then get transported through another mode discussed here within the soil.

- Temporal variability of the data: the measurements done in the same structure but with a gap of a few days or a few months have shown that the H₂ concentration is not constant and the monitoring done in Brazil have shown that the signal varies all over the day. Moreover, negative gradients such as decrease of signal with depth, and the positive ones where [H₂] does increase are both possible, as was observed in Oman. The modeling has stressed the shift between the [H₂] pulse and the near surface signal even if the soil is perfectly homogeneous. With “regular” hypothesis for the permeability of the soil and the distance to the source point up to few hours of time lag may be noted between the maxima. Few hours between the pulse and the higher concentration near the surface and few hours between the center of the fairy circle and the rim.

- The size of the circle is mainly a function of the transport time between the localized source and the surface. There is so no simple and univocal way to compute the depth of the source point from the size of the circle. For a given set of parameters, the circle will be larger when the source is deeper and has higher recharge pressure but if the permeability of the soil is lower for a given depth the circle will be smaller. Lower permeability leads to bigger contact time with H₂-consuming bacteria. In addition, if a ΔP exists the structure diameter will be also larger.

- Higher values on the rim and near zero in the center of the structure. If the soil is inhomogeneous, the flow of H₂ will be higher in the areas of high permeability. The H₂ flow is slow in the low permeability zone, as a result the gradient is getting smaller and the consumption by bacterial activity is higher. On the field, usually clay minerals fill the center of a depressed structure, a low permeability may form and the H₂ flow may vanish in the center and increase around the low-permeability layer. So pulse and low permeability zone in the center both explain the rim of high values but that the hypothesis of a low permeability on the center that allow to understand why the [H₂] signal remains low in the center which is above the source point.

One may also notice that the so far published data are not compatible with a high frequency pulsed source that lead
quickly, after few ten of meters of soil, to a constant signal near the surface. At this stage, it can be concluded that the inversion of the near surface [H₂] signal in order to deduce the source characteristics (depth, overpressure, [H₂]) is not attainable unless many subsurface parameters are known; among them, soil thickness, permeability variation and H₂ consumption potential are the most important. Permanent monitoring at different depths would also be interesting, for us even mandatory, to go further in the interpretation of the near-surface signal.

*Acknowledgement.* We thank ENGIE who supported this study and Andrey Myagkiy post-doctoral position. Complementary funding by CNRS/INSU through the TelluS Program is also acknowledged.

## References

Abrajano TA, Sturchio NC, Bohlke JK, Lyon GL, Poreda RJ, Stevens CM. 1988. Methane-hydrogen gas seeps, Zambales Ophiolite, Philippines: deep or shallow origin? *Chemical Geology* 71(1-3): 211-22. DOI: 10.1016/0009-2541(88)90116-7.

Auer LH, Rosenberg ND, Birdsell KH, Whitney EM. 1996. The effects of barometric pumping on contaminant transport. *Journal of Contaminant Hydrology* 24(2): 145-166.

Brunet F. 2019. Hydrothermal production of H₂ and magnetite from steel slags: a geo-inspired approach based on olivine serpentinization. *Front. Earth Sci.* 7. DOI: 10.3389/feart.2019.00017.

Cartwright J, Santamarina C. 2015. Seismic characteristics of fluid escape pipes in sedimentary basins: implications for pipe genesis. *Mar. Petroleum Geol.* 65: 126e140.

Charlou JL, Donval JP, Fouquet Y, Jean-Baptiste P, Holm N. 2002. Geochemistry of high H₂ and CH₄ vent fluids issuing from ultramafic rocks at the Rainbow hydrothermal field (36°14′N, MAR). *Chemical Geology* 191(4): 345-59.

Charlou J-L, Donval JP, Konn C, Ondreas H, Fouquet Y, Jean-Baptiste P, et al. 2010. High production and fluxes of H₂ and CH₄ and evidence of abiotic hydrocarbon synthesis by serpentinization in ultramafic-hosted hydrothermal systems on the mid-atlantic ridge. *Geophysical Monograph Series* 188. Washington DC: American Geophysical Union.

Conrad R, Seiler W. 1980. Contribution of hydrogen production by biological nitrogen fixation to the global hydrogen budget. *Journal of Geophysical Research Oceans* 85: 5493-5498.

Deronzier JF, Giouse H. 2020. Vaux-en-Bugey (Ain, France): the first gas field produced in France, providing learning lessons for natural hydrogen in the sub-surface? *BSGF-Earth Sci. Bull.* 5: 191. DOI: 10.1051/bsgf/2020005.

Deville E, Prinzhofer A. 2016. The origin of N₂-H₂-CH₄-rich natural gas seepages in ophiolitic context: A major and noble gases study of fluid seepages in New Caledonia. *Chem. Geol.* 440: 139-47.

Dominguez A, Bories S, Prat M. 2000. Gas cluster growth by solute diffusion in porous media. Experiments and automaton simulation on pore network. *International Journal of Multiphase Flow* 26(12): 1951-1979.

Fauria K, Rempel A. 2011. Gas invasion into water-saturated, unconsolidated porous media: implication for gas hydrate reservoirs. *EPSL* 312: 188-193. DOI: 10.1016/j.epsl.2011.09.042.

Firstov PP, Shirokov VA. 2005. Dynamics of molecular hydrogen and its relation to deformational processes at the Petropavlovsk-Kamchatskii geodynamic test site: evidence from observations in 1999-2003. *Geochemistry International.* 43(11): 1056-1064. [Translated from *Geokhimiya* 2005(11): 1151-1160].

Freeze RA, Cherry JA. 1979. Groundwater. NJ: Prentice Hall, 604 p.

Gay A, Lopez M, Cochonat P, Levaché D, Sermondadaz G, Seranne M. 2006. Evidences of early to late fluid migration from an Upper Miocene turbiditic channel revealed by 3D seismic coupled to geochemical sampling within seafloor pockmarks, Lower Congo Basin. V 23(3): 387-399. DOI: 10.1016/j.marpetgeo.2006.02.004.

Guélard J, Beaumont V, Rouchon V, Guyot F, Pillot D, et al. 2017. Natural H₂ in Kansas. Deep or shallow origin? *Geochem. Geophys. Geosyst.* 18(5): 1841-65.

Klein F, Tarnas J, Bach W. 2020. Abiotic sources of molecular hydrogen on Earth. *Elements* 16: 19-24.

Larin NV, Larin VN, Gorbatikov AV. 2010. Circular structures, caused by the deep seeping of hydrogen. Degassing of the Earth: Geotectonics, geodynamics, Deep fluids, oil and gas, hydrocarbons and life. Russian conference with international participation dedicated to the 100th anniversary of academician P.N. Kropotkin, 282 p. (in Russian)

Larin N, Zgonnik V, Rodina S, Deville E, Prinzhofer A, Larin VN. 2014. Natural molecular hydrogen seepage associated with surficial, rounded depressions on the European Craton in Russia. *Natural Resources Research* 24(3): 369-83.

Lazik D, Ebert S, Neumann PP, Bartholmai M. 2016. Characteristic length measurement of a subsurface gas anomaly - Amonitoring approach for heterogeneous flow path distributions. *International Journal of Greenhouse Gas Control* 47: 330-341.

Li X, Yortsos YC. 1995. Theory of multiple bubble growth in porous media by solute diffusion. *Chemical Engineering Science* 50(8): 1247-1271.

Lopez-Lazaro C, Bachau P, Moretti I, Ferrando N. 2019. Hydrogen solubility in aqueous NaCl solutions: from molecular simulation to equation of state. *BSGF*. DOI: 10.1051/bsgf/2019008.

Mahabadi N, Zheng X, Yun TS, van Paassen L, Jang J. 2018. Gas bubble migration and trapping in porous media: Pore-scale simulation. *Journal of Geophysical Research: Solid Earth* 123: 1060-1071. DOI: 10.1002/2017JB015331.

Maia A, Cartwright J, Andersen E. 2016. Shallow plumbing system inferred from spatial analysis of pockmark array. *Marine and Petroleum Geology* 77: 865-881. DOI: 10.1016/j.marpetgeo.2016.07.029.

Mayhew LE, Ellison ET, McCollom TM, Trainor TP, Templeton AS. 2013. Hydrogen generation from low-temperature water-rock reactions. *Nature Geoscience* 6: 478-484. DOI: 10.1038/ngeo1825.

Miller HM, Mayhew LE, Ellison ET, Kelemen P, Kubo M, Templeton AS. 2017. Low temperature hydrogen production during experimental hydration of partially-serpentinized dunite. *Geochemistry et Cosmochimica Acta* 209: 161-183. DOI: 10.1016/j.gca.2017.04.022.

Millington RJ, Quirk JP. 1961. Permeability of porous solids. *Transactions of the Faraday Society* 57: 1200-1207.

Monnin C, Chavagnac V, Boulart C, Ménez B, Gérard M, Gérard E, et al. 2014. Fluid chemistry of the low temperature hyperalkaline hydrothermal system of Prony Bay (New Caledonia). *Biogeo-sciences* 11(20): 5687-5706. DOI: 10.5194/bg-11-5687-2014.

Myagkiy A, Brunet F, Popov C, Krüger R, Guimarães H, Charlet L, et al. 2020. H₂ dynamics in the soil of an H₂-emitting zone (São Francisco Basin, Brazil): Microbial uptake quantification and reactive transport modelling. *Applied Geochemistry* 212. DOI: 10.1016/j.apgeochem.2019.104474.

Neal C, Stanger G. 1983. Hydrogen generation from mantle source rocks in Oman. *Earth and Planetary Science Letters* 66: 315-320.

Nivin VA, Pukha VV, Lovchikov AV, Rakhimov RG. 2016. Changes in the molecular hydrogen concentration in an underground mine
(Lovozero rare-metal deposit, Kola peninsula). *Doklady Earth Sciences* 471(2): 1261–1264.

Noel J, Godard M, Oliot E, I Martinez. 2018. Evidence of polygenetic carbon trapping in the Oman Ophiolite: Petro-structural, geochemical, and carbon and oxygen isotope study of the Wadi Dima harzburgite-hosted carbonates (Wadi Tayin massif, Sultanate of Oman). *Lithos* 323: 218–237.

Petrie E. 2014. Rock strength of caprock seal lithologies: Evidence for past seal failure, migration of fluids and the analysis of the reservoir seal interface in outcrop and the subsurface. Dissertation Thesis, Utah State University, Logan, Utah, 208 p.

Prinzhofer A, Deville E. 2015. Hydrogene naturel. La prochaine révolution énergétique ? Paris : Belin.

Prinzhofer A, Cisse CST, Diallo AB. 2018. Discovery of a large accumulation of natural hydrogen in Bourakebougou (Mali). *Int. J. Hydrogen Energy* 43: 19315–26.

Prinzhofer A, Moretti I, François J, Pacheco C, D'Agostino A, Werly J, et al. 2019. Natural hydrogen continuous emission from sedimentary basins: The example of a Brazilian H2-emitting structure. *International Journal of Hydrogen Energy* 44(12): 5676–5685. DOI: 10.1016/j.ijhydene.2019.01.119.

Sano Y, Urabe A, Wakita H, Wushiki H. 1993. Origin of hydrogen-nitrogen gas seeps, Oman. *Applied Geochemistry* 8: 1–8.

Sato M, Sutton AJ, McGee KA, Russel-Robinson S. 1986. Monitoring of hydrogen along the San Andreas and Calaveras Faults in Central California in 1980–1984. *Journal of Geophysical Research* 91(B12): 12315–12326.

Semlitsch RD. 2000. Size does matter: the value of small isolated wetlands. *Natl Wetl News* 22: 5–6.

Sherwood Lollar B, Onstott TC, Lacrampe-Couloume G, Ballentine CJ. 2014. The contribution of the Precambrian continental lithosphere to globa H₂ production. *Nature* 516: 379–382. DOI: 10.1038/nature14017.

Sibson RH. 2003. Brittle-fault controls on maximum sustainable overpressure in different tectonic regimes. *AAPG Bulletin* 87: 901–908. DOI: 10.1306/01290300181.

Terzaghi K, Peck RB. 1968. Soil mechanics in engineering practice. New York: John Wiley & Sons, 729 p.

Vacquand C, Deville E, Beaumont V, Guyot F, Sissmann O, Pillot D, et al. 2018. Reduced gas seepages in ophiolitic complexes: Evidences for multiple origins of the H₂-CH₄-N₂ gas mixtures. *Geochimica Cosmochimica Acta* 223: 437–461.

Welhan JA, Craig H. 1979. Methane and hydrogen in East Pacific rise hydrothermal fluids. *Geophys. Res. Lett.* 6: 829–831.

Zgonnik V. 2020. The occurrence and geoscience of natural hydrogen: A comprehensive review. *Earth Science Review* 203. DOI: 10.1016/j.earscirev.2020.103140.

Zgonnik V, Beaumont V, Deville E, Larin N, Pillot D, Farrell K. 2015. Evidences for natural hydrogen seepages associated with rounded subsident structures: the Carolina bays (Northern Carolina, USA). *Prog. Earth Planet. Sci.* 2: 31. DOI: 10.1186/s40645-015-0062-5.

## Appendix A Physical governing equations

Equations that govern the hydrodynamic system are represented by the equation of gas transport (1), continuity equation (2) and Darcy flow equation (3). Darcy's Law is combined with continuity equation and the ideal gas equation (4) becomes equation (5).

$$ \frac{\partial(\varphi c_i)}{\partial t} + \nabla \cdot (u c_i) = \nabla \cdot (D^* \nabla c_i) + R^*, \quad (1) $$

$$ \frac{\partial(\rho\varphi)}{\partial t} = -\nabla \cdot \left(\rho - \frac{k}{\mu}\nabla P\right), \quad (2) $$

$$ u = -\frac{k}{\mu}(\nabla P - \rho g), \quad (3) $$

$$ \rho = \frac{PM}{RT}, \quad (4) $$

$$ \frac{\partial(P\varphi)}{\partial t} = -\nabla \cdot \left(-\frac{k}{\mu}P\nabla P\right), \quad (5) $$

$$ D_{por} = \varphi \tau D_{atm}, \quad (6) $$

$$ \tau = \varphi^{1/3}, \quad (7) $$

$$ D^* = D_{por} + \alpha|u|. \quad (8) $$

For all-gas conditions as discussed in this paper, the tortuosity factor reduces to the porous medium value (Eq. (7)), according to the tortuosity correlation of Millington and Quirk (1961). Therefore, the system of equations is the following: where, $c_i$ is the concentration of H₂ in a gas phase within porous media, $R^*$ is a rate of H₂ consumption by bacteria, and $u$ signifies the gas velocity, $\phi$ denotes the porosity, $\tau$ is the tortuosity, $D^*$ is the sum of effective diffusion coefficient for the hydrogen gas in a porous media and a velocity-dependent dispersion term, so called mechanical dispersion, and $D_a$ is a diffusion coefficient of H₂ in air.

The diffusive transport in the gas was traditionally described by molecular diffusion. The diffusive flux was modeled with the use of Fick's law that has been optimized for porous media and takes into account porosity and tortuosity factors. The molecular mean free path of gas particles was assumed to be way larger than pore dimension and, thus, no influence of pore walls was considered. Preferential transport pathways in the soil, such as presence of cracks, worm channels or decayed and live root channels have not been considered in this model. The modelling case involves a deep origin H₂ source; the gas is, therefore, supposed to vent upward into the overlying sediments and then undergo the gas-soil interactions in the soil horizon. The source is represented by a source point at the bottom of the modelling domain, where the gas injects due to the pressure surplus and the dihydrogen enters the porous matrix of domain, the boundary conditions versus time is presented in Figure A1.

For numerical calculations reported here, the permeability used is 10 Da (or $10^{-11}$ m²), appropriate for a moderately permeable soil. Lower permeabilities should cause smaller effects.

The main physical parameters used in the model are listed in Table A1.
<figure><img src="image_11.png" /></figure>

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

<figure><img src="image_13.png" /><figcaption>**Fig. A1.** Signal at the source leaking point. A. Permanent source. B. Periodic source with one emission a day. C. Periodic source with four emissions a day.</figcaption></figure>

<table><caption>**Table A1.** Parameter values used in the model.</caption><tr><td>L</td><td>Layer thickness</td><td>30 m</td></tr><tr><td>K</td><td>Permeability</td><td>10 Da (1 Da = 10<sup>-12</sup> m<sup>2</sup>)</td></tr><tr><td>μ</td><td>Hydrogen viscosity</td><td>0.88e-5 Pa*s</td></tr><tr><td>φ</td><td>Porosity</td><td>0.3</td></tr><tr><td>P<sub>0</sub></td><td>Pressure</td><td>10<sup>5</sup> Pa</td></tr><tr><td>ρ</td><td>Hydrogen density</td><td>0.0832 kg/m<sup>3</sup></td></tr><tr><td>D<sub>atm</sub></td><td>Molecular diff. coeff.</td><td>7e-5 m<sup>2</sup>/s</td></tr><tr><td>ΔP<sub>fracture</sub></td><td>Surplus in gas pressure<br>at fracture</td><td>2 bar, 5 bar, 10 bar</td></tr></table>

Note: The numerical problem has been solved with parallel direct sparse solver (PARDISO), with relative and absolute tolerance of 0.01 and 0.001 respectively in all the simulation cases presented in this paper. The reference case has dimensions of 30*200 m with homogeneous permeability field. The size of simulation domains involved in the calculation is composed of 3794 nodes, which conducts correct convergence of the simulations.

**Cite this article as:** Myagkiy A, Moretti I, Brunet F. 2020. Space and time distribution of subsurface H₂ concentration in so-called “fairy circles”: Insight from a conceptual 2-D transport model, *BSGF - Earth Sciences Bulletin* 191: 13.