Induced Seismicity and Detailed Fracture Mapping as Tools for Evaluating HDR Reservoir Volume

Abstract

The main goal of this paper was to estimate the heat exchange rock mass volume of a hot dry rock (HDR) geothermal reservoir based on microseismicity location. There are two types of recorded microseismicity: induced by flowing fluid (wet microseismicity) and induced by stress mechanisms (dry microseismicity). In this paper, an attempt was made to extract events associated with the injected fluid flow. The authors rejected dry microseismic events with no hydraulic connection with the stimulated fracture network so as to avoid overestimating the reservoir volume. The proposed algorithm, which includes the collapsing method, automatic cluster detection, and spatiotemporal cluster evolution from the injection well, was applied to the microseismic dataset recorded during stimulation of the Soultz-sous-Forets HDR field in September 1993. The stimulated reservoir volume obtained from wet seismicity using convex hulls is approximately five times smaller than the volume obtained from the primary cloud of located events.

1. Introduction

Monitoring microseismicity in geothermal reservoirs is a useful tool for estimating the spatial distribution of fractures and evaluating fluid flow paths [1,2,3]. When fluid is injected, the fluid pore pressure changes causing a variation in effective stress, which can lead to rock failure [4]. Microseismic events identify the location of fractures generated or reactivated by injection [5]. When interpreting microseismic events recorded during hydraulic fracturing, it is important to understand what mechanism they were caused by [6]. During hydraulic fracturing, two mechanisms can result in the creation or reactivation of natural fractures. Microseismic events can be induced by fluid pressure (i.e., ‘wet’ microseismicity associated with a wet fracture network) or stress mechanisms (remote dry events associated with dry fractures) [6]. Wet microseismic events occur when the injected water reaches the preexisting fractures causing an increase in the fluid pressure. They are linked with a continuous slip as flow propagates. Dry events are induced by very quickly transferred mechanical stress change. They often occur in the vicinity of the hydraulic fracture tip [6]. The physical details of wet and dry events were studied by Maxwell [6].

One of the goals of hot dry rock (HDR) reservoir analysis is to estimate preferential fluid flow paths. It is important to distinguish dry volumes of a rock mass with a system of fractures that are not in any hydraulic connection with the stimulated fracture network. Excluding disconnected systems of fractures when interpreting fracture geometry can help avoid overestimating HDR reservoir performance [6]. Typically, the size of a microseismic cloud is used to estimate the stimulated reservoir volume (SRV). However, the cloud of microseismic locations is blurred and dispersed due to location errors, which in turn leads to an overestimation of microseismically active regions [6]. In the last few decades, significant improvements have been made to resolve structures within such clouds [7,8].

In this paper, an attempt has been made to evaluate SRV by investigating the volume of wet microseismicity, connected as a result of stimulation. The authors evaluated the percentage of wet events associated with fluid flow. It was assumed that the reservoir volume was directly proportional to the amount of wet microseismicity. The proposed algorithm includes the collapsing method [2], automatic cluster detection, and spatiotemporal cluster evolution from the injection well and was applied to the microseismic dataset recorded during stimulation of the Soultz-sous-Forets HDR field in September 1993. To increase the resolution of a seismic cloud, the collapsing method was applied. Based on relocated events, the seismic events were automatically divided into clusters using density-based spatial clustering of applications with noise algorithm (DBSCAN) [9]. Automatic cluster identification makes it possible to recognize parts of increased intensity of collapsed events—both wet and dry types. The spatiotemporal evolution of detected clusters was investigated to identify both events between which there was a possible hydraulic connection and events in which there was no connection. The latter usually appear as a distant, separated cluster. The convex hull algorithm [10] was used to estimate the SRV.

2. 1993 Soultz EGS Stimulation Injection Test and the Microseismic Data

2.1. The Soultz Geothermal Site

The Soultz-sous-Forets geothermal site is located in northeastern France, along the western edge of the Rhine Graben, a part of the West European Palaeogene Rift System [11]. This region is characterized by the occurrence of high-temperature gradients. The reservoir region in Soultz is situated between 2 to 4 km in granitic rock mass [12]. This granite is classified as an uplifted horst block [13]. It is estimated that the main granite body is a 331 ± 9 million years Hercynian monzogranite [14]. It is composed of phenocrystals of alkali feldspar in a matrix of quartz, plagioclase, biotite, and amphibole [15]. The crystalline basement is overlain by an approximate 1400 m layer of Triassic to Tertiary sediments [5]. More geological details are described in [11,16,17].

The site is located in a minor natural earthquake zone [18]. The region of Soultz has been investigated by many authors using microseismic, geophysical, geological, and hydraulic data [15,19,20,21]. Analysis of the cores and borehole geophysical logs indicated significant fracturing in the granitic massif with the main orientation of the fracture system NNW-SSE [22,23]. Although most of the fractures are naturally permeable, the permeability values are very low. This is due to sealing with hydrothermal deposits, such as calcite, silica, and clays [24].

Different parameters were used to study the stress field in the region of Soultz [25,26,27,28]. The orientation of the maximum horizontal principal stress ($S{H}_{max})$ up to a depth of 5 km is in the range of 170° N to 180° N [26], but is characterized by the presence of inhomogeneities in some depth intervals through major fracture zones [29].

2.2. Injection Test

During stimulation of the geothermal reservoir around the GPK-1 well in September 1993, 25,000 m³ of fresh water was injected at a rate of up to 40 L/s and pressure of up to 10 MPa (Figure 1 and Figure 2) [30]. In the September test, more than 10,000 microseismic events were recorded with three downhole 4-axis detectors (#4550, #4616, and #4601) and a hydrophone deployed in well EPS1 [31]. The frequency band of the detectors ranged from 10 Hz to 1 kHz, and the signals were digitized with a 5 kHz sampling frequency [31]. The borehole was open from the casing shoe at a depth of 2850 m and down to a depth of 3400 m [32]. Major intervals of water flow outlets from the well were located at measured depths of 2850–2900 m, 3090–3100 m, 3230–3240 m, and 3320–3325 m [33]. It was reported [34] that the most intensive seismic activity migrates from the top of the outlet point (2850–2900) toward deeper levels when the flow rate and the wellhead pressure were increased. The injection created a nearly vertical cloud of microseismicity with the following rough size: 0.5 km in width, 1.2 km in length, and 1.5 km in height [35]. The test revealed that the permeability and transmissivity increased [31,34].

2.3. Microseismicity

Data on the location of events and seismic stations are available on the IS-EPOS platform [36]. The company CSMA located a total of 10,743 events during the period of 2–22 September 1993. All the coordinates were converted in order to indicate distances from the GPK-1 wellhead. The distances shown on the X-axis indicate the distances from the well to the east and the Y distances to the north. Depths are expressed in true vertical depths. The locations of microseismic events and detectors are presented in Figure 3. The positions of the 4-component detectors are marked with red triangles and a hydrophone with a blue triangle. The direction of the minimum horizontal stress ($S{H}_{min})$ is shown in Figure 3B. The moment magnitude of these events ranged between −2.4 to 0.8.

The event rate was over 1000 per day (Figure 4). The link between the daily number of events and the injection flow rate is clearly demonstrated. Initially, the hypocenters of seismic events were located near the casing shoe at a depth range of 2800–3000 m [32]. At the end of the stimulation, the microseismic cloud reached a distance of 800 m from the injection well [37]. During the first days of the stimulation, seismicity extended downwards and was radially distributed. A few days after the stimulation began, the microseismic cloud expanded vertically and horizontally [32]. The microseismicity cloud grew gradually during the stimulation, with events distributed along subvertical planar structures [21,38,39]. The dominant structural trends inside the cloud range from NW-SE to N-S [32]. Jones et al. [32] mapped depth slices through the seismic data and showed the presence of several distinct strike directions for different sections of the seismic cloud, which suggest variation in the jointing patterns.

A lot of research has been conducted with the aim of finding structures inside this cloud of hypocenters. Most of these studies have been fragmentary in character and have only focused on selected structures. Phillips [3] performed a relative relocation of two chosen clusters and grouped events based on S-waveform shapes. Following their relocation, both clusters collapsed into well-defined patterns, consisting of two intersecting and terminating planar patches. The early activity inside these clusters indicates their link with permeable zones. Moriya et al. [35] analyzed the microseismicity cloud from the September 1993 test using high-resolution mapping methods, including the collapsing method, and the multiplet analysis method. Their analysis, based on microseismic multiplets, revealed that the fracture system was formed by the interconnection of fractures ranging in dimension from 4 to 26 m in the maximum principal direction of the existing tectonic stresses. They also claimed that the permeability of the system is enhanced by shear slippage and linking to more permeable fractures. The next set of results, provided by Moriya et al. [40], were drawn from an analysis of microseismic multiplets using the coherency function. They pointed out that multiplets were usually observed in the initial stage of the experiment. They concluded that seismic activation of fractures was caused by a small increase of pore pressure and extend to over several hundred meters. Evans et al. [21] applied the collapsing method to the cloud of microseismicity and claimed that it contains mainly linear and planar structures. They distinguished structure that appears from the proximity of a flow inlet at the beginning of the injection to 350 m downwards. A collapsed image of this structure suggests that it is tubular in shape due to a ±50 m horizontal location error. After the application of multiplets and clustering, the structure appears to be a subvertical, NNW-SSE striking fracture zone 10–20 m in width.

So far, research of the cloud of hypocenters recorded during the 1993 stimulation has focused on finding the dominant trend of such structures and their mechanism. In this article, the authors discussed which localized events can be related to the injected fluid flow, which is crucial for determining the performance of a geothermal reservoir.

3. Detailed Fracture Mapping

To accurately identify the extent of the geothermal reservoir and structures that may be related to the injection fluid flow, the authors proposed a 3-step procedure. During each step, some events were rejected from further analysis. Rejected events are treated as noise and authors assumed them to be related to dry fractures with no connection to the main reservoir. They were induced by changes in the pre-existing tectonic stresses.

The first step, i.e., the application of the collapsing procedure, was intended to increase the resolution of the seismic cloud. Following the collapsing method, those events that were too distant location-wise from the others were rejected. In the next step, the authors performed the automatic cluster detection algorithm, i.e., DBSCAN. The clustering result made it possible to reject from further analysis events that occurred in spatial isolation from existing clusters and consider them as noise. The spatiotemporal evolution of detected clusters was studied with the aim of distinguishing between events possibly connected hydraulically and those that occur suddenly and are remote. A flowchart describing wet microseismicity identification is shown in Figure 5. Details of each procedure will be discussed in the following sections. The final step in the analysis is a convex hull algorithm, the purpose of which is to estimate the SRV.

3.1. Collapsing Method

The reservoir structure was assessed by examining induced microseismicity. The first step increasing the resolution of the seismic cloud was the collapsing algorithm. The collapsing method has been used to reduce location error in the seismic cloud since 1997 [2]. This statistical optimization technique simplifies the seismic cloud structure originally blurred by location error. The key step in the collapsing procedure involves establishing the error ellipsoids for all seismic locations. The ellipsoid size is determined by partial errors in the measured parameters. Events are shifted towards local source densities inside their error ellipsoids iteratively until the translation vector takes the form of a χ² distribution with 3 degrees of freedom. The distribution of sources after collapsing may reflect the real distribution of the hypocenters not affected by location error. The collapsing procedure exhibited an ability to delineate seismic structures in volcanic seismicity [2], mining seismicity [41], and induced seismicity from geothermal reservoirs [21].

The results of applying the collapsing method to all of the hypocenters recorded during the September 1993 injection test is shown in Figure 6A. If an event did not move during the algorithm, it meant that it did not contain any other events inside their uncertainty ellipsoid. These 1110 events (Figure 6B) were excluded from further analysis. After their relocation, some structures inside the seismic cloud are more distinctive. In general, the collapsed cloud of hypocenters is more heterogeneous. This means that there is a tendency towards events within the cloud to concentrate. To resolve the microseismic structures inside the collapsed cloud, clustering was performed on displaced points (Figure 6C).

3.2. Density-Based Spatial Clustering of Applications with Noise

According to density-based clustering approaches, clusters are regions where the density of points is much higher than outside the clusters [42]. One of the most commonly used density-based clustering algorithms is density-based spatial clustering of applications with noise (DBSCAN). This algorithm was proposed in 1996 [9] and is still widely used today [43,44,45,46]. One of the advantages of this method is the ability to extract clusters of arbitrary shape in dense regions and recognize noise points. The undoubted advantage of the algorithm is that there is no need to specify the number of groups to be found. The algorithm works assuming that points are assigned to the same group if they are density-reachable from each other.

The DBSCAN algorithm requires estimating two parameters: the minimum number of points and the radius. After conducting a clustering series with different parameters, the optimal results were obtained with the assumption that a minimum of 30 events within a neighborhood radius of 50 m constitues a cluster. The neighborhood radius value was chosen based on the distance of the data set using k-nearest neighbour distance graph.

The result of applying the DBSCAN algorithm to the relocated microseismic cloud is presented in Figure 7. The DBSCAN algorithm was able to identify 8 clusters and 1048 noise points. The points marked in gray represent noise and were rejected for further analysis. The points in other colors indicate particular clusters. Detected clusters have different shapes, sizes, and population numbers. The population in each cluster is listed in

Table 1. Population in each cluster.

Cluster Number	1	2	3	4	5	6	7	8
Nr of points	7446	449	250	101	185	56	44	54

. These clusters represent high-density regions of event occurrence. It is assumed that they symbolize fractures created as a result of fluid flow or as disturbances of concentrated stress and pressure away from the fluid flow path. The method of distinguishing between clusters according to wet and dry seismicity is described in the next chapter.

4. Spatiotemporal Evolution of Clusters

Fracture network connectivity is a key factor shaping generated microseismicity, stimulated rock volume, and production [47]. To find which clusters detected by DBSCAN are possibly hydraulically connected, the spatiotemporal evolution of clusters was analyzed. The spatiotemporal evolution of a collapsed and clustered cumulative microseismic cloud is presented in daily intervals in Figure 8. Additionally, the minimum, maximum, and mean values for injection rate and wellhead pressure for each day are shown. The different colors of the seismic events indicate that they belong to a different cluster.

Various factors were investigated as part of the spatiotemporal analysis: the cluster growth rate measured by the number of events gathered in the cluster, the distance from the injection well, the size of the cluster, and the time of its appearance. The absolute number of daily events divided by clusters is shown in Figure 9. This histogram shows only those events that the DBSCAN algorithm classified into clusters. The clusters growth rate are presented in Figure 10A–H. These graphs show the daily percentage growth of each cluster in relation to the maximum number of events in a given cluster. The colors in Figure 8, Figure 9 and Figure 10 match the colors representing clusters in Figure 7. The authors analyzed whether events in clusters occur gradually or suddenly. It was assumed that if seismic events appeared gradually, they were related to the fluid flow. However, when the appearance of a cluster was sudden, it tended to be associated with a change in tectonic stress.

Each of the detected clusters was analyzed in order of their appearance. The highest number of points is contained in the light blue cluster—about 80% of all points subjected to cluster detection. It is not assumed to be a single fracture. This area consists of many smaller structures that are very close to one another, which may indicate a connection between them. Throughout the stimulation period, this structure grew uniformly on both sides of the well. The first events inside this cluster occurred near the casing shoe at a depth of 2800 m and migrated downwards. After 4 days, events were observed at shallower depths. This structure was approximately 1 km deep, 0.5 km wide, and 0.8 km long. An analysis of the depth of occurrence of the biggest detected cluster, its position in relation to the GPK-1 well (origin of the coordinate system), and the rate of growth, lead to the conclusion that the events inside this cloud occurred as a result of contact with fluid flows.

The second cluster that appears in the seismic cloud is the one marked in light green. It contains less than 5% of all the points that were subject to the automatic clustering algorithm. Its rate of growth is variable, but not rapid. The first events classified as part of this cluster appeared along with an increase in the injection rate. They took place when the largest cluster was growing, and were located very close to it. This may indicate a hydraulic connection between these clusters.

The third cluster, marked in red, containing approximately 2.5% of all clustered events, appeared as an expanded light green cluster after an increase in the injected fluid. This may indicate the presence of a fluid flow path between them.

The next chronological structures are marked in pink and heather. They consist of 1% and 2% of the entire population of points that underwent clustering, respectively. Events inside these structures began to occur after the next injection rate increased and as the range of the biggest cluster with a hydraulic connection to the well increased. This may suggest that these events represent a case of wet seismicity.

The sixth structure farthest from the borehole and which at the same time had one of the fastest growth rates is shown in brown. This cluster has no connection to the well or to the structures that can be permeable. This means that it could not have arisen as a result of flowing water. It can be concluded that events inside this cluster were triggered by tectonic stress changes.

The next cluster, marked in yellow, appeared in its entirety at a time when the injection rate was at its height. It is located in close proximity to the structure which the authors considered wet. The structure composed by these events may be permeable.

The fastest rate of point growth occurred in the dark green cluster. Less than 1% of all the points are located in this structure. The first observed events in this structure occurred when the injection flow rate increased. This fact may indicate that events within that cluster will be classified as wet.

In order to assign to each cluster to its presumed origin—tectonic or directly related to the fluid flow—a spatiotemporal analysis of the evolution of the clusters was performed. The cluster that had rapid growth was temporally and spatially separated and had no close connections with other neighboring clusters; it is probably of tectonic origin due to changes in stress in the rock mass. The clusters inside which the events appeared gradually along with the increasing injection rate are probably related to the fluid flow.

5. Discussion

Despite many studies conducted to find structures within the microseismic cloud recorded during injection in Soultz in 1993, the research presented in this paper was the first attempt to divide the entire cloud into clusters. The authors attempted to determine which of the localized seismic events may be involved in heat transfer. The methodology proposed by the authors does not explicitly identify wet microseismic events but provides a premise for the identification of events that may be associated with the flowing fluid.

It is important to note that the Soultz reservoir is naturally fractured. The big complex structure may be formed of smaller natural discontinuities, which as a result of fracturing combine into larger structures according to the percolation process. These natural discontinuities are potential paths for fluid flow [15]. The orientation of detected clusters inside microseismic cloud appears to be in agreement with the N-S to NNW-SSE direction of most natural fractures demonstrated by Genter et al. [22]. The largest detected cluster (cluster 1) contains approximately the same structure described by Evans et al. [21]. The structure identified as Evans Line [21] ranging from a flow zone to 300 m downwards is considered to be permeable, supporting significant flow.

The distinction between wet and dry microseismic events was based on the spatiotemporal evolution of clusters and the rate of point growth in them. The authors considered the possibility that events inside the brown cluster may have been triggered by stress mechanisms. These events will not contribute to the SRV due to the absence of any hydraulic connection. The remaining detected clusters may indicate the presence of permeable structures. The authors calculated that about 80% of all localized events can be characterized as wet microseismicity. A flowchart describing wet microseismicity identification and the number of events in each processing step is presented in Figure 5.

The last step in the analysis presented in this article involved estimating the stimulated reservoir volume using the convex hull algorithm. The volume engaged by seismic activity can be approximated as a volume of the microseismic cloud [48]. The convex hull of the original microseismic cloud—before the 3-step procedure was performed—is presented in Figure 11A. The volume of the primary cloud defined by the microseismic events is 7.6 × 10⁸ m³. The authors proposed only calculating the SRV of those clusters that can be in a hydraulic connection with the well. The convex hulls for each cluster are shown in Figure 11B. The authors assumed that the stimulated reservoir volume is the sum of the volumes of all convex hulls. The volumes of each structure are shown in

Table 2. Volumes of structures identified as possibly permeable.

Cluster Number	1	2	3	4	5	7	8
Volume [${\mathrm{m}}^3$]	$1.25\times {10}^8$	$5.11\times {10}^6$	$3.63\times {10}^6$	$6.92\times {10}^5$	$1.33\times {10}^6$	$1.66\times {10}^5$	$2.13\times {10}^5$

. The sum of the volumes of the structures identified as possibly permeable is equal to 1.4 × 10⁸ m³. It is important to note that the real volume engaged with the heat transfer is slightly greater because there is a connection between particular clusters. With this algorithm, it is impossible to calculate the volume of paths connecting individual clusters. The calculated volume is the smallest volume containing the wet fracture system.

The dominant source mechanism during hydraulic fracturing is shear slip, which is represented by a double-couple source model. This mechanism is triggered by an increase in the pore pressure reducing normal stress along existing fractures [49]. Some events occur near the tips of generated fractures and are characterized by significant non-shear (tensile opening and closing) components [6]. Although new fractures are inevitably created by a non-shear mechanism during stimulation of the geothermal reservoir and should be related to the volume of the stimulated rock mass, a large number of shear mechanism events are also connected to the fluid flow. Dry microseismic events resulting from the disturbance of stresses are characterized by a shear stress mechanism. However, when only the source mechanism is examined, it cannot be unequivocally determined whether shear or tensile fracturing are associated with wet microseismicity [6,50]. This analysis does not always produce good results. Hence, for this reason a different analysis was performed in this article. Instead of resorting to the source mechanism, only localized seismic events were used to extract wet microseismicity. These events were subject to relocation, clustering, and spatiotemporal evolution.

6. Conclusions

The microseismicity detected during hydraulic stimulation in the Soultz HDR field in September 1993 provides a basis for making a rough estimation of the number of events associated with fluid flow. To calculate the stimulated reservoir volume, the authors proposed an algorithm comprising the collapsing method, DBSCAN automatic cluster detection, and spatiotemporal cluster evolution from the injection well. The use of the collapsing method made it possible to simplify the hypocenter cloud and exclude points that were too distant from others. Clustering provided a platform for detecting areas with an above-average intensity of seismic events. However, this method could not be used to distinguish between wet and dry events. Therefore, the next step was a spatiotemporal analysis of detected clusters. Various factors were investigated in this analysis: the rate of cluster growth measured by the number of events, the distance from the injection well, the size of the cluster, and the time of the appearance. This approach allowed the authors to identify seven clusters containing events that can be related to the wet fracture system. These structures were deemed to be the preferential paths for the fluid flow.

When the volume of a geothermal reservoir is estimated according to the extent of the recorded microseismic events, this volume is overestimated. The microseismicity cloud contains events that are related to the flow of the injected fluid as well as events resulting from pre-existing tectonic stresses in the area. The stimulated reservoir volume obtained from wet systems of fractures, excluding paths connecting individual clusters using convex hulls, is approximately five times smaller than the volume obtained from the entire microseismic cloud.

The method proposed by the authors to identify events directly associated with the hydraulic fracture network can be applied to any geothermal reservoir. Accurately identifying dry fractures is crucial for correctly estimating the volume of a geothermal reservoir containing only hydraulically connected fracture systems. The results presented in the article are significant from the perspective of evaluating HDR reservoir performance.