Natural H2 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 H2 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 H2 concentration is variable with time. Here, the time and space evolution of H2 concentration is simulated using a 2-D advective-diffusive model of H2 transport in porous media. Several parameters have been tested as the depth and periodicity of the H2 point source (pulsed), bacterial H2 consumption and permeability heterogeneities of the soil. The radius of the structure is linked to the time spent by the H2 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, H2 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 H2 concentration values will be highly dependent on the time at which the measurement is performed. H2 monitoring through a sensor network is thus mandatory to characterize the H2 dynamics in the soil of fairy circles.

Des émissions de H2 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 H2 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 H2 du sol y est variable. Sur la base d’une modélisation du transport de H2 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 H2, la consommation de H2 dans le sol et aux hétérogénéités du sol. 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’H2. 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.

Introduction

H2 transfer in subsurface

Subsurface H2 gas concentrations far above atmospheric level and subsequent H2 flow towards the surface have been reported in an increasing number of publications. H2 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 H2-forming reactions and processes are likely to produce H2 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 H2 under abiotic conditions (e.g., Brunet, 2019 and references therein). H2 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. H2 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, Fe2+ rich rocks are present, as it has been suggested as a source for the continental H2 accumulations in Mali (Prinzhofer et al., 2018). Radiolytic H2 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 1011 moles per year. In addition, H2 may originate from mantle degassing (Larin et al., 2010) since H2 is the initial and major component of the solar system. Mechanoradical H2 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 H2/CH4 ratio in a gas field in the Jura front belt (Deronzier and Giouse, 2020) however its overall role for large scale H2 production remains unknown. Bacterial activity is also prone to produce H2 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 H2 transfer in continental settings characterized by a soil cover and which leads to the formation of fairy circles. In the absence of soil cover, H2 emanations are usually reported along faults. In the Oman ophiolite, for example, H2-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), H2 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 H2.

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

Shape of the H2 emitting structures in the presence of a soil layer

Due to their circular shape, the structures from which H2 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 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 H2 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 H2 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 H2 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.

Time-resolved H2-concentration measurements within the fairy circles

Up to now, even though the list of basins where H2 emanations have been reported is getting longer (Zgonnik, 2020), both time-resolved and spatial distribution measurements of H2 concentration in soils from fairy circles are still rare. However, the datasets cited in the previous section all show that the soil H2-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 H2 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 H2 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).

In Brazil, the monitored circular structure is a depression with a diameter close to 500 m. The published H2-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 H2 flow for about 6 h a day. H2 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 H2 concentration maximum was reached during the night. A sensor installed in the center of the structure for few months recorded only a few H2 concentration pulses but most of the 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 H2 concentration in the soil, near the surface, is changing with time along the day is shown by all data.

The effect, on H2 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 H2 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 H2 flows. Basically, the consumption of H2 by bacteria increases with the residence time of H2 in the soil whatever the transport mode, advective or diffusive.

There are only few published data of H2 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 H2 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 H2 emission do take place. A single zero H2 concentration (that we will note [H2] further in this paper) value at a given time cannot be interpreted as a lack of H2 flow in a dynamic system. Based on a numerical model of H2 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 H2 concentration; (3) the variation through time of the near surface signal. We will therefore discuss the causes of the apparent randomness of the H2 concentration recorded near the surface.

Transport of H2 across a soil

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 H2 consumption. The parameters, experimentally derived for the soils of Brazilian H2-emitting structure, have been used in simulations. Two types of calculation were performed, namely with and without H2 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 [H2] at which microbiological uptake kinetics was measured, so up to few thousandths of ppm. For the sake of simplicity, modelled H2 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 [H2] 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 H2 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.

Boundary conditions

Soil permeability: the average permeability used in most of the simulations reported here is 10 Da (or 10−11 m2), 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 H2 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 H2 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 H2, being both poorly permeable we will consider that the H2 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 H2 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 H2 leakage: Since all the available monitoring data show that the H2 flow is periodic we will mainly use that hypothesis. H2 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 H2 leakage will be modeled for comparison.

Amount of the H2 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 different ΔPfracture will be tested and the size of the area influenced by the H2 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 H2 consumption will not be discussed here in detail.

Results and discussion

Permanent flow of H2

Figure 2 shows the concentration in the soil above a permanent leakage of H2 located at 30 m. After a few hours, H2 reaches the surface and the horizontal concentration profile has a bell shape with a maximum located just above the leak point. Then, the H2 concentration profile spreads laterally. Therefore, the size of the area affected by the H2 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 H2 soil concentration has a “bell shape” distribution with its maximum located at the center of the structure

  • The requested time for H2 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 [H2] 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 H2 signal fits with the observed data described above.

A second set of runs has been performed adding the effect of H2 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. H2 consumption limits the radial expansion of the H2 emitting zone. With bacterial activity, after few days the permanent state is achieved, there is an equilibrium between the H2 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 H2 in the soil to reach the surface and the kinetics of the H2 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 [H2] values neither the variabilities of the [H2] versus time. In addition, one may note that without consumption of H2 in the soil the distribution of H2 is continuing and growth of the structure continues over years. In the literature, the soil H2 concentration is always variable but nobody noticed so far large changes in the structure size through time. Whatever its origin, bacterial or not, H2 consumption in the soils seems mandatory to explain this stability in case of continuous flow.

H2 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 H2 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 H2 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 [H2] 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 H2 concentration gradient, [H2]zforumla 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 [H2] in the center is zero.

The results of this modeling, even very simple, are a key to understand literature data where apparently random [H2] 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 H2 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.

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 H2 in the soil, H2 bearing area could be very large, a few hundreds of meters even for a near surface leakage. The scavenging of the H2 in the soil (e.g., by bacteria) will limit the size of the structure for a given point source flux and depth. The H2 travel distance which can continuously grow in the absence of H2 consumption becomes limited when H2 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 H2 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 H2 “bubble”, which is prone to vary. Since the time dependency of the pulse is the same (Fig. A1B), that increase of pressure of the H2 bubble means also than a higher quantity of H2 is released into the soil from the leak point.

The source being pulsed, the [H2] 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 [H2] will be displayed (Fig. 5), which corresponds to the maximum [H2] 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 H2 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 H2 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 H2 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 H2 pulse is more or less pressurized, that means also that the size increases if the quantity of the H2 released at each pulse is getting larger. It seems particularly interesting to explain the occurence of small structures within the large structures as observed in Carolina bay.

The fact that H2 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.

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 H2 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 H2 consumption by bacterial activity is considered.

Figure 6A shows a snapshot of the resulting H2 flow. In that case, H2 is first moving upward and laterally as previously but when H2 flow reaches the low permeability area, it gets around it taking the fastest way to the surface. The ring of high H2 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 [H2] near surface at 80 cm depth, which corresponds to the maximum [H2] 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 H2 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 H2 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 H2 is consumed before it reaches the surface but the “travel time” is only slightly affected by the bacteria. The consumption of H2 affects the H2 gradient, when it is strongly dependent on the permeability. In that case, after 2 days the signal is stable, the deep H2 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.

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 [H2] 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 H2 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 H2 in the soil is globally smaller since the amount of H2 released at the source is globally smaller. We do not know the exact quantity of H2 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 H2 concentration is stable as for a constant source. However, a constant [H2] at given location in time and with depth does not fit with the observed data.

Other potential processes that may influence the H2 concentration in the soil

As explained previously we modeled the transport of H2 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 H2 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 H2 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 H2 consumption by inorganic processes is not modeled, for instance adsorbtion of H2 in the shally soil in the center of the structure could also amplify the H2 signal attenuation phenomenon in the center of the structures but is not computed here.

External causes of H2 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 H2 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.

Conclusions

The present 2-D model of H2 transport (advection and diffusion) in a soil layer from a leak point of periodic [H2] 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 H2. 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 H2 in the soil.

  • A pulsed source allows to explain the spatial variability of the [H2] 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 H2 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 H2 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 [H2] does increase are both possible, as was observed in Oman. The modeling has stressed the shift between the [H2] 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 H2-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 H2 will be higher in the areas of high permeability. The H2 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 H2 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 [H2] 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 [H2] signal in order to deduce the source characteristics (depth, overpressure, [H2]) is not attainable unless many subsurface parameters are known; among them, soil thickness, permeability variation and H2 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.

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).  
(φci)t+(uci)=(D*ci)+R*,
graphic
(1)
 
(ρφ)t=(ρkμP),
graphic
(2)
 
u=kμ(Pρg),
graphic
(3)
 
ρ=PMRT,
graphic
(4)
 
(Pφ)t=(kμPP),
graphic
(5)
 
Dpor=φτDatm,
graphic
(6)
 
τ=φ1/3,
graphic
(7)
 
D*=Dpor+α|u|.
graphic
(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, ci is the concentration of H2 in a gas phase within porous media, R* is a rate of H2 consumption by bacteria, and u signifies the gas velocity, ϕ denotes the porosity, τ 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 Da is a diffusion coefficient of H2 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 H2 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 m2), 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.

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

This is an Open Access article distributed under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.