International Journal of Geophysics
Volume 2021 (2021), Article ID 6647834, 14 pages
https://doi.org/10.1155/2021/6647834

Review Article

Inferring Hydraulic Pressure Changes from Induced Seismicity Observations: Three Showcases from Geothermal Reservoirs

Christopher Koch,¹ Stefan Baisch,¹ Elmar Rothert,¹ and John Reinecker²

¹Q-con GmbH, 76887 Bad Bergzabern, Germany
²GeoThermal Engineering, 76135 Karlsruhe, Germany

Correspondence should be addressed to Christopher Koch; koch@q-con.de

Received 22 October 2020; Accepted 1 June 2021; Published 28 June 2021

Academic Editor: Veronica Pazzi

Copyright © 2021 Christopher Koch et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

We apply a recently developed approach for inferring in situ fluid pressure changes from induced seismicity observations to datasets from geothermal reservoirs at St. Gallen (Switzerland), Paralana (Australia), and Cooper Basin (Australia), respectively. The approach, referred to as seismohydraulic pressure mapping (SHPM), is based on mapping the seismic moment of induced earthquakes. Relative fluid pressure changes are inferred from the stress deficit of fracture patches slipping repeatedly. The SHPM approach was developed for the specific scenario, where induced earthquakes occur on a single, larger-scale plane with slip being driven by the regional stress field. We demonstrate that this scenario applies to the three datasets under investigation, indicating that geothermal systems in crystalline rock could typically be fault-dominated. For all datasets, individual earthquake source geometry could not be determined from source spectra due to the attenuation of the high signal frequencies. Instead, SHPM was applied assuming a constant stress drop in a circular crack model. Absolute values of inferred pressure change scale with the assumed stress drop while the spatiotemporal pattern of pressure changes remains similar even when varying stress drop by one order of magnitude. We demonstrate how the associated mismapping of seismic moment tends to average out when hypocentres are densely spaced. Our results indicate that SHPM could provide important information for calibrating numerical reservoir models.

1. Introduction

Fluid injection into the subsurface, as performed in different energy technologies, can be accompanied by induced seismicity. Already in the 1950s, fluid overpressure has been identified as a dominating mechanism causing seismicity.

While the conceptual understanding of the underlying physical processes has existed for decades [1], it turned out to be difficult to describe the earthquake process quantitatively, mostly due to the lack of in situ stress and pressure measurements.

For example, interpretations of the level of hydraulic overpressure at which earthquakes are triggered by waste water disposal in Oklahoma vary between <0.01 MPa [2] and 1 MPa [3], depending on model assumptions. Here, fundamental characteristics of fluid pressure propagation in the subsurface are not well understood, which may lead to nonoptimal seismic risk mitigation strategies like reducing the disposed fluid volume at those wells located within a certain distance to observed seismicity [4]. For geothermal reservoirs, understanding hydraulic processes in the subsurface quantitatively is paramount since these are directly linked to the systems’ economics.

In a recent study, Baisch [5] has proposed a new approach for inferring in situ fluid pressure changes from induced seismicity observations, which could be an important step towards calibrating numerical models. The approach, hereafter referred to as seismohydraulic pressure mapping (SHPM), is based on relative stress and pressure changes occurring between repeated slippage of the same patch of a fracture. Baisch [5] demonstrates the performance of the SHPM approach using induced seismicity data from an EGS reservoir (Enhanced Geothermal System) in the Cooper Basin, Australia. For this, he has made several assumptions which are closely linked to certain characteristics of the induced seismicity in the Cooper Basin, namely, that all earthquakes have occurred on the same, large-scale (planar) fault and share a common source mechanism, and the slip area of neighbouring events is overlapping. Although these conditions may seem to be very specific, Baisch [5] speculates that similar conditions may actually prevail at other locations, thus allowing application of the same method.

Investigation of this hypothesis is the focus of the current study, where we present three datasets meeting the assumptions formulated above. With these datasets, we seek to further study the performance and limitations of the SHPM approach. We infer relative pressure changes and discuss the plausibility of SHPM pressure changes based on physical principles. We do, however, not use SHPM observations for calibrating subsurface models. Thereby, we limit our scope to demonstrating the general applicability of the SHPM approach to these datasets. The aspect of model calibration by comparison with analytical solutions (e.g., [6]) or numerical models will be addressed in subsequent studies.

Our first dataset stems from a geothermal reservoir at St. Gallen (Switzerland) where fluid injection activities have caused a series of a few hundred earthquakes occurring on a single fault zone as indicated by the distribution of relative hypocentre locations [7]. Compared to the dataset used by Baisch [5], the number of earthquakes in this dataset is two orders of magnitude smaller and overlap of the slip area of neighbouring earthquakes is less pronounced. With this dataset, we aim to explore the performance of SHPM for small datasets with only few fracture patches slipping repeatedly.

Our second dataset stems from geothermal reservoir stimulation at Paralana (Australia). Different to the situation in the Cooper Basin, published hypocentre distributions indicate that seismicity has occurred on a complex (volumetric) fracture network. Through reprocessing of the data, it could be demonstrated that the apparent fracture complexity is a processing artefact making the dataset suitable for the application of the SHPM methodology.

With our third dataset, we aim to investigate the performance of SHPM in a scenario where initial stress conditions are strongly heterogeneous. The dataset stems from the same geothermal reservoir in the Cooper Basin that has already been studied by Baisch [5]. Our analysis, however, is based on data from a restimulation of the Habanero#1 well in 2005, where the initial stimulation has resulted in heterogeneous stress conditions prior to restimulation.

2. Materials and Methods

2.1. Methodology

The SHPM methodology is described in detail by Baisch [5] and can be summarized as follows [5].

For an arbitrary patch on a fault with coordinates , which is repeatedly activated in the course of continuous pore pressure increase, it can be demonstrated that relative changes of fluid pressure can be related to relative changes of shear- and normal-stresses and by

with denoting the time of seismic activation and denoting the static coefficient of friction.

Let denote the slip area of event . The relative stress changes in Equation (1) can be expressed by a contribution resulting from slip on initialized at and a contribution resulting from slip occurring outside between time and as well as from aseismic, poroelastic, and thermoelastic contributions:

Following the line of argumentation of Baisch [5], we assume that aseismic, poroelastic, and thermoelastic contributions are of secondary order in our data examples and can be ignored, i.e.,

Furthermore, we make the assumptions that fault cohesion can be neglected and that changes of normal stresses are of secondary order, in which case fluid pressure changes are solely related to changes of the shear stress (Equation (1)).

In the specific scenario where induced earthquakes occur on patches of the same, large-scale fault and share a common rake direction, stress interaction between different patches of the fault (the stress contribution at location resulting from an earthquake located outside ) can be modelled analytically using the slider block concept and Okada’s [8] semianalytical solutions for efficient calculation of relative stress changes (compare Figure 1 of Baisch [5]).

2.2. Geothermal Project St. Gallen (Switzerland)

2.2.1. Project Overview

The St. Gallen deep geothermal project is located in the Swiss Molasse Basin between Lake Constance and the Alps. A 3D seismic survey provides a detailed image of the geothermal target fault zone, a system of major NNE-SSW striking normal faults constituting the St. Gallen Fault Zone (SFZ) [9].

The first geothermal well (St. Gallen GT-1) was drilled in 2013 into fractured Mesozoic limestones to a total depth of 4,253 mBGL ([10], see also Figure 1(b)). Following a chemical stimulation in July 2013, a gas kick was detected that had to be controlled through injection of heavy mud leading to a series of induced earthquakes with magnitudes up to [7]. Well cleaning activities in mid-September 2013 before production tests were accompanied by a new increase in seismicity resulting in a project-wide, total amount of 347 locatable induced earthquakes [7].

(a)

(b)

(c)

Figure 1: (a) Time-encoded relative hypocentre locations in perspective view with best-fitting plane determined by linear regression. Trajectory of St. Gallen GT-1 is indicated by a black line. Red shaded area denotes presumably falsely located event cluster (see Induced Seismicity above). Coordinates are given with respect to the top of GT-1 (47.415200°N 9.328801°E). (b) The structural model of St. Gallen displaying the well path of the geothermal well (black line), top Palaeozoic (violet) and faults (as interpreted from the 3D seismic survey [9]) as well as hypocentres of the induced seismicity (from the Swiss Seismological Service (SED), red dots). View from northeast. (c) Earthquake source-centred, stereographic projection (lower hemisphere) of -wave polarity data for all 347 events. The beachball indicates the compound fault plane solution and is approximately similar to the largest event ( 3.5) determined by Diehl et al. [16] with strike/dip/rake-values of 124°/72°/174°, respectively.

Geothermal activities were monitored with a local seismic station network consisting of 9 surface seismometers (SGT01-05: broadband, SGT06-09: short-period), two accelerometers, and a single borehole geophone deployed at 205 m depth. After the 3.5 seismic event occurred in the morning of the 20th of July 2013, the monitoring network was extended by 3 additional short-period seismometers (SGT10-12; see Figure 2). In the following months until May 2014, several stations were removed and some stations were upgraded, relocated, and renamed (SGT13-17) [11].

Figure 2: St. Gallen station network in map view. Wellhead position of injection well GT-1 is indicated by a yellow star. Blue triangles denote surface stations (SGT01-SGT12); red triangle denotes the instrument deployed in shallow borehole (SGT00). The position of the two accelerometers and the new position of station SGT06, SGT07, and SGT09 after slight relocation in October 2013 (renamed SGT13-17) are not shown. The inset in the upper left shows a map of Switzerland with the location of the St. Gallen project site. Data: [13, 14]. Map data from OpenStreetMap (published under ODbL).

Given the comparatively small volume of only 649 m³ fluid injected during well control operations [10], the seismicity response was larger than could be expected by, e.g., the McGarr [12] model. It could be speculated that an additional fluid/gas volume from an (artesian) gas pocket could have contributed to the mass balance. This hypothesis, however, could not be tested so far due to the absence of reservoir pressure readings. SHPM provides the opportunity to further test this hypothesis.

2.2.2. Induced Seismicity

Relative hypocentre locations determined by Diehl et al. [7] indicate that the induced seismicity aligns along a subvertical fault, steeply dipping towards northwest [7].

It has to be noted that the small cluster of seismicity at the lower edge of the seismicity (red shaded area in Figure 1(a)) is interpreted to be an artefact presumably caused by a anomaly and is most likely part of the main cluster [7].

Although event-specific fault plane solutions are not well constrained, a compound fault plane solution determined following Baisch et al. [15] (Figure 1(c)) is reasonably consistent with the mechanism determined by Diehl et al. [7] for the 3.5 seismic event and also with the general strike of the faults mapped in this region. This single-fault interpretation meets the requirements for applying the SHPM methodology.

We attempted to determine source parameters by fitting -wave spectra with a model but found that corner frequencies are strongly depending on the assumed attenuation model.

2.3. Geothermal Project Paralana (Australia)

2.3.1. Project Overview

The Paralana geothermal project, located in the Poontana Basin in south-eastern Australia, was launched in 2007. Two geothermal wells (Paralana-1b and Paralana-2) were drilled to a depth of 1,760 mBSL and 3,960 mBSL, respectively [17]. Paralana-2 reached the Mesoproterozoic basement consisting of undifferentiated felsic porphyry [17] and intersected a naturally fractured system with an overpressure of ~27 MPa at a depth of approximately 3,640-3,820 mBSL [18]. After completion, the well was perforated at a depth of 3,639 to 3,645 mBSL [17].

In January 2011, a “Diagnostic Fracture Injection Test” was performed injecting a small amount of fluid (approx. 14 m³) into the Paralana-2 well [17]. During the 4-hour test, wellhead pressure increased to ~60 MPa and approximately 300 induced earthquakes with a maximum magnitude of 1.4 were detected [17, 18]. Seismicity was monitored with a 12-station network consisting of surface seismometers and geophones deployed in shallow boreholes and in Paralana-1b at a depth of approximately 1,760 mBSL [17, 18]. For further details of the seismic monitoring network, we refer to Hasting et al. [18].

In July 2011, a massive hydraulic stimulation was conducted through the well Paralana-2, injecting a total fluid volume of 3,100 m³ at a maximum flow rate of 27 l/s and wellhead pressure up to 62 MPa (Figure 3). For the massive stimulation, the seismic monitoring network was temporarily extended by an additional set of five surface seismometers and four accelerometers [18] (see Figure 4).

(a)

(b)

(c)

Figure 3: (a) Injection rate, (b) wellhead pressure, and (c) daily rate of induced events as a function of time during the main stimulation of Paralana-2. Time window is restricted to 10th-17th of July 2011 when the surface stations of the seismic monitoring network were in operation. Injection was terminated on 15th of July 2011 after injecting a total volume of approximately 3,100 m³ [17]. Seismic activity continued after the injection period. Hydraulic data taken from daily operation reports. Note that the time basis of hydraulic data reported by Albaric et al. [17] and Hasting et al. [18] seems to be partly inconsistent with the daily operation reports.

Figure 4: Paralana station network in map view. Wellhead position of injection well Paralana-2 is indicated by a yellow star. Blue triangles denote surface stations; red triangles denote instruments deployed in boreholes. The inset shows the location of the Paralana project site. Map data from OpenStreetMap (published under ODbL).

During hydraulic stimulation, approximately 11,000 seismic events were automatically detected [18]. Quality control of the processed data was limited to events with a magnitude resulting in 2,600 events occurring between 10th of July and 23rd of August [18].

Hypocentre locations obtained by Hasting et al. [18] exhibit pronounced spatial scattering. Noting that the hypocentre distribution of Hasting et al. [18] was biased by the false interpretation of a reflected phase, Albaric et al. [17] obtained a more focussed image of the hypocentre distribution (see Figure 9 in [17]), but the reprocessed hypocentre distribution still indicates a complex fracture network.

2.3.2. Data Processing

As part of the current study, we reprocessed the induced seismicity data of the massive hydraulic stimulation of the Paralana geothermal project. Reprocessing was restricted to those events considered by Hasting et al. [18] in the time window of 10th of July 2011 to 17th of July 2011 when the surface stations of the seismic monitoring network were in operation. For this time window, a total number of 1,757 induced events could be relocated.

Besides the issue of a reflected phase masking S-onsets at up to five borehole stations, quality control revealed that some instruments exhibit a time drift (Figure 5). To account for the issue of instrument time drift, hypocentres were relocated using differential travel times , which are to first order not sensitive to instrument time drift :

with denoting origin time, denoting the hypocentral distance, and , denoting seismic wave velocities. Rearranging Equation (4) yields:

Figure 5: Seismogram example of an induced earthquake (July 13th, 2011 03 : 23 : 27 UTC) recorded at stations B06 and B07 of the Paralana seismic network. Differential times between (red) and (green) phase onsets are not consistent with the absolute timing of the onset indicating instrument time drift in the order of 0.4 s.

The velocity term on the right-hand side of Equation (5) was calibrated by assuming that the first twenty events occurring during the main stimulation are located close to the flow exit at the wellbore as typically observed in EGS stimulations (e.g., [19]). The resulting station-dependent velocity terms are almost constant for these events indicating that they are indeed approximately colocated. Data from five borehole stations (B01, B02, B03, B05, and B08) had to be discarded as only -wave arrivals could be determined.

2.3.3. Induced Seismicity

Figure 6(a) shows the resulting event distribution, where hypocentres predominantly align along a subhorizontal plane-like structure at approximately 3,700 mBSL depth. Scattering of hypocentre locations around the best-fitting plane approximately follows a normal distribution with a standard deviation of 30 m (Figure 6(a)), indicating that the apparent vertical thickness of the hypocentre distribution might be dominated by hypocentre location errors and that earthquakes could actually have occurred on a single plane.

(a)

(b)

(c)

Figure 6: (a) Absolute hypocentre locations in perspective view with best-fitting plane determined by linear regression. The histogram on the right shows vertical distances of the absolute hypocentre locations to the best-fitting plane. Dotted curve indicates normal distribution with 0 mean and 30 m standard deviation. (b) Map view of fault plane solutions determined for the 40 strongest events. Fracture intersection of the Paralana-2 injection well is indicated by a yellow star. (c) Earthquake source-centred, stereographic projection (lower hemisphere) of -wave polarity data for the 40 strongest earthquakes. Beachball indicates compound fault plane solution. Coordinates are given with respect to the top of Paralana-2 (30.21287°S/139.72850°E).

Fault plane solutions were determined for the 40 strongest events based on -wave polarities and amplitude ratios following the approach of Baisch et al. [20]. Furthermore, a compound fault plane solution was determined for the same set of events following Baisch et al. [15]. The resulting fault plane solutions consistently indicate oblique thrust faulting along the plane outlined by the hypocentre distribution (Figures 6(a) and 6(b)), conforming with focal mechanisms determined by Albaric et al. [17].

This scenario of a planar structure where seismicity is driven by a common mechanism fulfills the general requirements for applying the SHPM approach.

Our attempt to determine source parameters by fitting -wave spectra with a model resulted in the same conclusion derived for the previous dataset, i.e., that corner frequency cannot be reliably resolved from source spectra due to the (unknown) impact of attenuation.

2.4. Cooper Basin Habanero#1 Restimulation (Australia)

2.4.1. Project Overview

Our third showcase stems from the same geothermal reservoir in the Cooper Basin where the SHPM method has already been applied [5].

Here, we study the 2005 restimulation of the geothermal well Habanero#1 which targets the same subhorizontal fault as the Habanero#4 well studied by Baisch [5]. Both Habanero#1 stimulations were performed prior to the Habanero#4 stimulation. The same data processing was applied to all Habanero datasets. For data processing of the initial Habanero#1 stimulation in 2003, we refer to Baisch et al. [19] and for the second stimulation in 2005 to Baisch et al. [15].

Figure 7 shows the hypocentre distribution associated with the Habanero#1 restimulation. As noted in Baisch et al. [15], seismicity did not start near the injection well but at the outer rim of the previously stimulated region. This is an expression of the so-called Kaiser effect [21], providing direct evidence of strongly heterogeneous stress conditions on the fault prior to restimulation.

(a)

(b)

(c)

Figure 7: (a) Hypocentre locations of the seismicity induced during the 2005 restimulation of well Habanero#1 (red star) in map view [15]. Each seismic event is displayed by a globe scaled to the event magnitude. Colour encoding denotes occurrence time according to the legend. Seismic activity associated with the initial stimulation of well Habanero#1 in 2003 [19] is indicated by grey dots. Contour line indicates 2003 seismicity. (b) Absolute hypocentre locations in perspective view with best-fitting plane determined by linear regression. (c) Earthquake source-centred, stereographic projection (lower hemisphere) of -wave polarity data of the dominating type 1 events, see [15]. Beachball indicates compound fault plane solution. Coordinates are given with respect to the top of Habanero#1 (27.81639°S/140.75972°E).

Consistent with previous findings, source parameters determined from -wave spectra strongly depend on the assumed attenuation model.

3. Results and Discussion

3.1. Inferred Pressure Changes, St. Gallen

Applying SHPM to this dataset closely follows the procedure described by Baisch [5]. In a first step, hypocentre locations were projected onto the best-fitting plane determined by least squares. Subsequently, projected hypocentres were rotated into the horizontal plane and aligned with the direction of the maximal horizontal stress () striking 160° [22].

Similar to observations for the Cooper Basin dataset reported by Baisch [5], corner frequencies and hence source radii could not be reliably resolved from source spectra for this dataset (see Materials and Methods). Therefore, mapping of the seismic moment as described in Baisch [5] was performed assuming a constant, event-independent stress drop.

Figure 8 shows the spatial distribution of the maximum fluid pressure changes inferred from repeated slip for assumed stress drop values of 0.1, 0.5, and 1.0 MPa and a coefficient of friction of , respectively.

(a)

(b)

(c)

Figure 8: Spatial distribution of the maximum cumulative fluid pressure changes inferred from repeated slip in top view of best-fitting plane. Pressure changes are stated in megapascals according to the colour map, which is saturated at 20 MPa. The 3D arrow indicates northern (green), western (red), and depth (blue) direction after projection and rotation of hypocentre locations. A constant stress drop of (a) 0.1 MPa, (b) 0.5 MPa, and (c) 1.0 MPa was assumed for determining fluid pressure changes. Contour line denotes region of main seismic activity.

Since our analysis is based on relative hypocentre locations, the wellbore location relative to the seismicity distribution is not exactly resolved and therefore not shown in Figure 8. We might speculate that the SHPM overpressure maximum at -750/800 (Figure 8(c)) is located close to the flow exit of the wellbore.

Interpretation of SHPM results shows that fluid pressure changes can be inferred in a limited reservoir area only with few sampling points in time. This directly results from the low earthquake density compared, e.g., to Baisch et al. [5]. Another consequence of the low earthquake density is the strong dependency of inferred fluid pressure changes on the assumed stress drop value. Pressure changes are generally smaller than 4 MPa when assuming a stress drop of 0.1 MPa and increase to up to 20 MPa if 1 MPa stress drop is assumed.

While the individual earthquake source geometry tends to averaging out in case of high earthquake density (as is the case in the Baisch et al. dataset [5]), it has a first-order impact on inferred pressure changes in the current dataset.

To discuss the plausibility of the SHPM results, inferred pressure changes were compared to calculated reservoir pressure at the wellbore (Figure 9). Downhole pressure in the open-hole section of the wellbore was calculated from pressure measurements at the wellhead while accounting for the weight of the fluid column using measured fluid density. This approximation is not exact, e.g., due to temperature-density and mixing effects in the wellbore and highly dynamical density values around 14th of September 2013, but we nevertheless consider the calculated downhole pressure as a first-order approximation of actual reservoir pressure.

Figure 9: (a) Temporal evolution of cumulative fluid pressure changes determined from repeated slips at the 4 locations in the reservoir indicated by corresponding colours in the small figure inset. A constant stress drop of 0.5 MPa is assumed. Contour line in the figure inset denotes the main region of seismic activity. (b) Temporal evolution of the calculated pressure changes in the wellbore at a depth of 3,818 m (open-hole section) during well control operations in July 2013 and well cleaning activities in September/October 2013.

Figure 9 shows the temporal evolution of inferred pressure changes at four different reservoir locations assuming a stress drop of 0.5 MPa (a) and the relative pressure changes in the open-hole section of the wellbore (b).

Locations 1 (red) and 2 (violet) exhibit repeated slip which occurred during the phase of gas kick and well control operations in July 2013. Deduced pressure changes at these two locations resulted in a systematic increase up to approximately 10 MPa.

Locations 3 (blue) and 4 (green) mainly show a pressure increase for the time period where well cleaning activities prior to the production test were performed. The pressure increase determined by SHPM is in the range of up to 11.5 MPa and consistent with basic physical principles, i.e., that the pressure changes in the reservoir do not exceed pressure changes applied at the wellbore.

We note, however, that the absolute level of inferred pressure changes is depending on the assumed stress drop, thus limiting quantitative interpretations. Although the 0.5 MPa stress drop assumption leads to a consistent picture of in situ pressure changes, we feel that we cannot answer whether or not an artesian gas pocket with higher overpressure has contributed to the induced seismicity sequence.

As mentioned before, it is also important to notice that Diehl et al. [7] speculate that some of the deepest earthquakes (red shaded area in Figure 1(a)) may actually be located inside the main cluster of seismic activity. This would have an implication for the inferred pressure changes, which would be smaller in the deeper section of the fault and larger at the main cluster of seismic activity.

3.2. Inferred Pressure Changes, Paralana

Following the same procedure as outlined above for inferring SHPM pressure changes, hypocentre locations were projected onto the best-fitting plane determined by least squares. Subsequently, projected hypocentres were rotated into the horizontal plane and aligned with the direction of the maximal horizontal stress () striking 97°N [17]. Given that corner frequencies and hence source radii could not be determined for this dataset (see Materials and Methods), mapping of the seismic moment was performed assuming a constant, event-independent stress drop.

Figure 10 shows the spatial distribution of the maximum fluid overpressure inferred from repeated slip for assumed stress drop values of 0.1, 0.5, and 1.0 MPa, respectively. A coefficient of friction of is assumed, which lies within the range deduced for porphyric rock from laboratory experiments [23].

(a)

(b)

(c)

Figure 10: Spatial distribution of the maximum cumulative fluid pressure changes inferred from repeated slip in top view of best-fitting plane. Pressure changes are stated in megapascals according to the colour map, which is saturated at 60 MPa. The arrow indicates northern direction. A constant stress drop of (a) 0.1 MPa, (b) 0.5 MPa, and (c) 1.0 MPa was assumed for determining fluid pressure changes. Contour line denotes region of main seismic activity as outlined by the hypocentre distribution. Coordinates are given with respect to the top of Paralana-2 (30.21287°S/139.72850°E).

Comparing inferred pressure for the different stress drop values tested here, we note that the absolute level of pressure changes scales with the assumed stress drop value, whereas the relative spatiotemporal distribution remains, to first order, similar.

Figure 11 shows the temporal evolution of inferred pressure changes at 8 reference points for the 0.1 MPa stress drop model. Except for the red curve, we note similar characteristics as observed by Baisch [5] in the sense that (i) maximum pressure changes are obtained near the injection well and tend to decrease with distance from the well, (ii) the magnitude of pressure changes is smaller, but in the same order of magnitude as changes measured at the wellhead, and (iii) the delay of the pressure signal increases with distance from the injection well, as could be expected for hydraulic pressure diffusion. The characteristic of the red curve, however, is slightly different. Here, repeated slip started early, and the inferred curve of pressure changes (red) intersects the blue curve. This can be explained by assuming that the pressure level at which seismicity started at the red location is ≥3.5 MPa higher compared to the respective level at the blue location. An increased pressure level at the red location might be a result from the prestimulation where seismicity has already occurred near the injection well. Due to the Kaiser effect, the pressure level at which subsequent seismicity is induced at this location is higher.

Figure 11: (a) Temporal evolution of cumulative fluid pressure changes determined from overlapping slips at the 8 locations in the reservoir indicated by corresponding colours in the small figure inset. A constant stress drop of 0.1 MPa is assumed. Contour line in the figure inset denotes the main region of seismic activity after projection and rotation operations. The black star denotes the injection point of the well Paralana-2. (b) Temporal evolution of the pressure changes measured at the wellhead of Paralana-2 during the main stimulation in July 2011. Pressure changes were determined by subtracting the initial overpressure of 27 MPa from pressure readings.

Even when accounting for an overpressure level in the order of 5 MPa at which seismicity is initially induced, inferred pressure changes are still smaller than the pressure changes measured at the wellbore, consistent with basic physical principles.

3.3. Inferred Pressure Changes, Habanero#1 Restimulation

We followed the same procedure for inferring SHPM pressure changes as outlined for the two other datasets. Hypocentre locations were projected onto best-fitting plane followed by rotation operations to align the projected hypocentres with the direction of the maximum horizontal stress () striking approximately westnorthwest–eastsoutheast [24].

Figure 12 shows the spatial distribution of maximum pressure changes assuming a constant stress drop of 0.1, 0.5, and 1.0 MPa, respectively. Consistent with our previous findings, the absolute level of inferred pressure changes scales with the assumed stress drop, while the spatiotemporal pattern of pressure changes remains similar.

Figure 12: Spatial distribution of the maximum cumulative fluid pressure changes inferred from repeated slip in top view of best-fitting plane. The time at which the maximum cumulative fluid pressure change is reached varies over the reservoir. Pressure changes are stated in megapascals according to the colour map, which is saturated at 50 MPa. The yellow star denotes the fault intersection of the well Habanero#1. The arrow indicates northern direction. A constant stress drop of (a) 0.1 MPa, (b) 0.5 MPa, and (c) 1.0 MPa was assumed for determining fluid pressure changes. Contour line denotes region of main seismic activity as outlined by the hypocentre distribution. Coordinates are given with respect to the top of Habanero#1.

Figure 13 shows the temporal evolution of inferred pressure changes at 9 reference points for the 0.1 MPa stress drop model. In contrast to SHPM pressure changes inferred for the Habanero#4 stimulation [5], the SHPM curves during the Habanero#1 restimulation (Figure 13) systematically intersect each other. This could be expected given that the critical pressure at which the first slip occurred is strongly heterogeneous due to the Kaiser effect resulting from the initial stimulation.

Figure 13: (a) Temporal evolution of cumulative fluid pressure changes determined from overlapping slips at the 9 locations in the reservoir indicated by corresponding colours in the small figure inset. A constant stress drop of 0.1 MPa is assumed. Contour line in the figure inset denotes the main region of seismic activity after projection and rotation operations. The black star denotes fault intersection of the well Habanero#1. (b) Temporal evolution of the pressure changes measured at the wellhead of Habanero#1 during the restimulation in September 2005. Pressure changes were determined by subtracting the initial overpressure of 33 MPa from pressure readings.

Based on an analytical solution for fluid injection into an infinite homogeneous fault [25], the maximum overpressure prevailing at the end of the initial Habanero#1 stimulation has been modelled (Figure 14(a)). In a simplified approximation, the resulting spatial distribution of maximum overpressure has been taken as an estimate of the pressure level at which further seismicity can occur.

Figure 14(b) shows absolute overpressure reconstructed by adding the initial pressure to . In a simplified model of a fault with homogeneous hydraulic properties, basic physical principles require that the overpressure curves at different locations do not intersect. Besides some data scattering, progression of the pressure curves in Figure 14(b) are indeed almost in line with this condition showing no major intersections anymore. Furthermore, the maximum reconstructed overpressure in the reservoir remains approximately 3 MPa below the maximum injection pressure, which is a plausible scenario for a diffusion process. These observations add confidence in the physical meaning of the values and demonstrate that the SHPM approach can be applied even if initial stress conditions are strongly heterogeneous.

Figure 14: (a) Analytical solution [25] for pressure resulting from injection into an infinite fault with homogeneous properties. Parameters were selected to approximate the initial 2003 stimulation. Homogeneous fault transmissibility of  Dm and a constant injection rate of 25 l/s are assumed. The spatial pressure distribution is shown after 11.5 days. Red squares denote the distance of the 9 reference points and the associated pressure level which is assumed to equal the level P₀(r) at which further seismicity can occur. (b) Same as Figure 13 while absolute overpressure is reconstructed by adding the initial pressure to .

3.4. Common Characteristics

For the three datasets investigated here, source parameter could not be determined from source spectra due to attenuation of high frequencies. Following Baisch [5], we have assumed an earthquake-independent, constant value for stress drop. This simplified assumption results in a mismapping of seismic slip of individual earthquakes, which tends to average out with an increasing number of overlapping slips.

In our datasets, the amplitude of inferred pressure changes is sensitive to the assumed stress drop, while the spatiotemporal distribution of pressure changes remains similar. In practice, this bears the possibility of deducing an average (constant) stress drop by matching SHPM pressure changes near an injection well with measured pressure changes in the well.

To further quantify the impact of the assumed stress drop on the amplitude of inferred pressure changes, we have determined the maximum number of repeated slips in the datasets investigated here. Figure 15 shows the percentage of maximum cumulative fluid pressure change for different stress drop models as a function of the maximum number of repeated slips.

Figure 15: Percentage of maximum cumulative fluid pressure change for different stress drop models as a function of the maximum number of repeated slips for the three geothermal projects mentioned here (St. Gallen, Paralana, Habanero#1-restimulation) and the Habanero#4-stimulation analysed in Baisch [5].

The impact of the assumed stress drop tends to decrease with increasing number of repeated slips, demonstrating that the individual source geometry tends to average out with increasing number of overlapping slips [5]. Nevertheless, the impact of the assumed stress drop is still pronounced even for the largest dataset of the Habanero#4 stimulation indicating that a considerable large number of repeated slips are necessary to get an efficient “averaging out”-effect.

4. Conclusions

A recently developed approach for inferring in situ pressure changes from induced seismicity observations (SHPM) is applied to data obtained from injection activities conducted in three geothermal reservoirs. We find that all three reservoirs meet the general requirements for applying the SHPM method in terms that the reservoirs are fault-dominated, where most of the induced seismicity aligns along a single plane with slip being driven by the regional stress field. For the Paralana (Australia) geothermal reservoir, we find that the previous reservoir interpretation of a complex fracture network was biased by instrument issues.

Our findings support the hypothesis that geothermal systems in crystalline rock could typically be fault-dominated.

The stress drop of induced earthquakes, which is an input parameter for SHPM, could not be determined for the current datasets due to attenuation of the high signal frequencies. Instead, a constant stress drop value is assumed. We speculate that other datasets from geothermal reservoirs exhibit similar limitations.

While the absolute value of inferred pressure changes scales with the assumed stress drop value, the spatiotemporal pattern of pressure changes remains similar even when varying stress drop by one order of magnitude. If event-specific stress drop cannot be determined from seismogram data, we suggest assuming a constant stress drop value which can be adjusted by matching measured pressure changes in an injection well. For large, densely spaced hypocentre distributions, the impact of the assumed stress drop on inferred pressure changes decreases.

Data Availability

Time-continuous seismogram recordings of the seismological stations near St. Gallen are available on the public website of the Swiss Seismological Service (SED) (http://arclink.ethz.ch). Station locations were provided by the SED [13, 14] within the “Science for Clean Energy (S4CE)”-project in the framework of the European Union’s Horizon 2020 research and innovation program. The induced earthquake catalogue (occurrence time, magnitude, and relative hypocentre locations) and hydraulic data were provided by S4CE-project partner St. Galler Stadtwerke in the framework of the European Union’s Horizon 2020 research and innovation program. All information is also available at IS EPOS [26], retrievable via doi:10.25171/InstGeoph_PAS_ISEPOS-2018-007. Cooper Basin and Paralana waveform data, seismic event catalogues, and hydraulic data are available to the scientific community through the Department of State Development (SDS) of South Australia (https://energymining.sa.gov.au/geothermal/data_And_publications/Reservoir_Stimulation_and_Microseismicity_Data).

Conflicts of Interest

The authors declare that there is no conflict of interest regarding the publication of this paper.

Acknowledgments

We would like to thank Elisabeth Brzoska who has performed the reprocessing of the Paralana dataset. This project has received funding from the European Union’s Horizon 2020 research and innovation program under grant agreement (764810).

References

M. K. Hubbert and W. W. Rubey, “Role of fluid pressure in mechanics of overthrust faulting,” GSA Bulletin, vol. 70, no. 2, pp. 115–206, 1959.

K. M. Keranen, M. Weingarten, G. A. Abers, B. A. Bekins, and S. Ge, “Sharp increase in central Oklahoma seismicity since 2008 induced by massive wastewater injection,” Science, vol. 345, no. 6195, pp. 448–451, 2014.

M. Schoenball, F. R. Walsh, M. Weingarten, and W. L. Ellsworth, “How faults wake up: the Guthrie-Langston, Oklahoma earthquakes,” The Leading Edge, vol. 37, no. 2, pp. 100–106, 2018.

K. Bosman, A. Baig, G. Viegas, and T. Urbancic, “Towards an improved understanding of induced seismicity associated with hydraulic fracturing,” First break, vol. 34, no. 7, 2016.

S. Baisch, “Inferring in situ hydraulic pressure from induced seismicity observations: an application to the Cooper Basin (Australia) geothermal reservoir,” Journal of Geophysical Research: Solid Earth, vol. 125, no. 8, 2020.

Z. Fan and R. Parashar, “Analytical solutions for a wellbore subjected to a non-isothermal fluid flux: implications for optimizing injection rates, fracture reactivation, and EGS hydraulic stimulation,” Rock Mechanics and Rock Engineering, vol. 52, no. 11, pp. 4715–4729, 2019.

T. Diehl, T. Kraft, E. Kissling, and S. Wiemer, “The induced earthquake sequence related to the St. Gallen deep geothermal project (Switzerland): fault reactivation and fluid interactions imaged by microseismicity,” Journal of Geophysical Research - Solid Earth, vol. 122, no. 9, pp. 7272–7290, 2017.

Y. Okada, “Internal deformation due to shear and tensile faults in a half-space,” Bulletin of the Seismological Society of America, vol. 82, no. 2, pp. 1018–1040, 1992.

S. Heuberger, P. Roth, O. Zingg, H. Naef, and B. P. Meier, “The St. Gallen fault zone: a long-lived multiphase structure in the North Alpine Foreland Basin revealed by 3D seismic data,” Swiss Journal of Geosciences, vol. 109, no. 1, pp. 83–102, 2016.

H. Naef, “Die Geothermie-Tiefbohrung St. Gallen GT-1,” Swiss Bulletin for Applied Geology, vol. 20, no. 1, 2015.

B. Edwards, T. Kraft, C. Cauzzi, P. Kastli, and S. Wiemer, “Seismic monitoring and analysis of deep geothermal projects in St Gallen and Basel, Switzerland,” Geophysical Journal International, vol. 201, no. 2, pp. 1022–1039, 2015.

A. McGarr, “Maximum magnitude earthquakes induced by fluid injection,” Journal of Geophysical Research - Solid Earth, vol. 119, no. 2, pp. 1008–1019, 2014.

Swiss Seismological Service, National Seismic Networks of Switzerland, 1983.

Swiss Seismological Service, “Temporary deployments in Switzerland associated with aftershocks and other seismic sequences,” 2005.

S. Baisch, R. Vörös, R. Weidler, and D. Wyborn, “Investigation of fault mechanisms during geothermal reservoir stimulation experiments in the Cooper Basin, Australia,” Bulletin of the Seismological Society of America, vol. 99, no. 1, pp. 148–158, 2009.

T. Diehl, J. Clinton, T. Kraft et al., “Earthquakes in Switzerland and surrounding regions during 2013,” Swiss Journal of Geosciences, vol. 107, no. 2–3, pp. 359–375, 2014.

J. Albaric, V. Oye, N. Langet et al., “Monitoring of induced seismicity during the first geothermal reservoir stimulation at Paralana, Australia,” Geothermics, vol. 52, pp. 120–131, 2014.

M. Hasting, J. Albaric, V. Oye et al., “Real-time induced seismicity monitoring during wellbore stimulation at Paralana-2 South Australia,” in New Zealand Geothermal Workshop 2011 Proceedings, p. 9, Auckland, New Zealand, 2011.

S. Baisch, R. Weidler, R. Vörös, D. Wyborn, and L. de Graaf, “Induced seismicity during the stimulation of a geothermal HFR reservoir in the Cooper Basin, Australia,” Bulletin of the Seismological Society of America, vol. 96, no. 6, pp. 2242–2256, 2006.

S. Baisch, E. Rothert, H. Stang, R. Vörös, C. Koch, and A. McMahon, “Continued geothermal reservoir stimulation experiments in the Cooper Basin (Australia),” Bulletin of the Seismological Society of America, vol. 105, no. 1, pp. 198–209, 2015.

S. Baisch and H.-P. Harjes, “A model for fluid-injection-induced seismicity at the KTB, Germany,” Geophysical Journal International, vol. 152, no. 1, pp. 160–170, 2003.

I. Moeck, T. Bloch, R. Graf et al., “The St. Gallen project: development of fault controlled geothermal systems in urban areas,” in Proceedings World Geothermal Congress 2015, p. 5, Melbourne, Australia, 2015.

M. Beblo, A. Berktold, U. Bleil et al., Physical Properties of Rocks, Subvolume b, vol. 1b, Springer-Verlag, Berlin/Heidelberg, 1982.

S. D. Reynolds, S. D. Mildren, R. R. Hillis, J. J. Meyer, and T. Flottmann, “Maximum horizontal stress orientations in the Cooper Basin, Australia: implications for plate-scale tectonics and local stress sources,” Geophysical Journal International, vol. 160, no. 1, pp. 332–344, 2005.

C. V. Theis, “The relation between the lowering of the piezometric surface and the rate and duration of discharge of a well using groundwater storage,” Transactions, American Geophysical Union, vol. 16, no. 2, pp. 519–524, 1935.

K. Leptokaropoulos, S. Cielesta, M. Staszek et al., “IS-EPOS: a platform for anthropogenic seismicity research,” Acta Geophysica, vol. 67, no. 1, pp. 299–310, 2019.