Induced seismicity in EGS reservoir: the creep route

Observations in enhanced geothermal system (EGS) reservoirs of induced seismicity and slow aseismic slip ruptures on related faults suggest a close link between the two phenomena. We base our approach on the case study of the EGS site of Soultz-sous-Forêts where seismicity has been shown in particular during the 1993 stimulation to be induced not only by fluid pressure increase during stimulation but also by aseismic creeping effects. We propose an interpretation of the field observations of induced seismicity using a laboratory experiment that explores, in great detail, the deformation processes of heterogeneous interfaces in the brittle-creep regime. We track the evolution of an interfacial crack over 7 orders of magnitude in time and 5 orders of magnitude in space using optical and acoustic sensors. We show that a creep route for induced seismicity is possible when heterogeneities exist along the fault. Indeed, seismic event occurrences in time and space are in strong relation with the development of the aseismic motion recorded during the experiments. We also infer the statistical properties of the organization of the seismicity that shows strong space-time clustering. We conclude that aseismic processes might drive seismicity besides the classical effects related to fluid pressure and show that a creep route for induced seismicity is possible.


Background
The exploitation of a geothermal reservoir can lead to the occurrence of an abundant seismicity, particularly during phases of hydraulic stimulation. This seismicity has a strong societal impact, as it can be felt by the population, especially in densely populated area as in Europe, e.g., M L = 3.4, Basel, 2006(Häring et al. 2008. Reducing this impact is thus one of the main challenges in the development of geothermal energy production. In order to understand the mechanical processes responsible for seismic activity, all possible mechanisms for induced seismicity have to be studied (Majer et al. 2007). Among them the link between earthquakes and aseismic deformations is only at a starting state (Cornet et al. 2007;Bourouis and Bernard 2007). Deciphering this link between seismicity and aseismic motion can not only help to mitigate the risk posed by seismicity but also www.geothermal-energy-journal.com/content/2/1/14 help to monitor and model the evolution of the geothermal reservoir. The mechanism relating earthquakes and aseismic processes is still elusive due to the difficulty of imaging these phenomena with large spatiotemporal variability at depth. However, a good example of joint seismic and aseismic deformation has been obtained at the enhanced geothermal system (EGS) site of Soultz-sous-Forêts (France)  in particular during the 1993 water injection experiment for stimulation of the naturally fractured granite reservoir below 2,500 m.

Fluid-induced seismicity at Soultz-sous-Forêts
Seven large-scale fluid injections have been performed at Soultz-sous-Forêts on the four deep wells: GPK1 was stimulated in 1993, GPK2 in 1994 and 1995 for the upper reservoir and 2000 for the deep reservoir, GPK3 in 2003, and GPK4 in 2004and 2005(Dorbath et al. 2009Huenges and Ledru 2010). During all these stimulations, induced seismicity has been extensively studied (Charléty et al. 2007;Cuenot et al. 2008). For all of them, fluid pressure increase has been considered as the main source of seismicity initiation (Shapiro et al. 1999). The 1993 stimulation is certainly the first stimulation for which a precise influence of the fluid pressure has been demonstrated (Cornet et al. 1997;Cornet et al. 2007). Figure 1 reproduces the main results that have been obtained. Interpretation in the Mohr space relies on the difference between Coulomb's failure criterion of an initially intact material (i.e., a yield criterion for the initiation) and the friction law (see Figure 1a) (Burov 2009;Jaeger et al. 2009) on existing interfaces. Typically, when the fluid pressure p is increased, the effective normal stress σ e is reduced and the material reaches first the friction criterion (in green in Figure 1a) if pre-existing fractures exist and are favorably oriented within the local stress field. In this case, a shear slip, or hydro-shear, is initiated related to radiated waves well characterized by an earthquake double-couple solution. This process explains very well the onset of seismicity at the bottom of the stimulated zone (below 3000 msee Figure 1b,c). The orientation of the seismic cloud with respect to the maximum horizontal stress orientation gives the friction angle φ ≈ 30°, i.e., a static friction coefficient μ ≈ 0.60. Because of dilatancy of the fracture during shear offset, permeability is significantly and permanently increased during hydro-shear. At higher flow rates (≥ 18 l/s), pressure stabilizes to a plateau for which seismicity migrates upwards along the maximum horizontal stress direction (in red in Figure 1b,c). During that regime, the effective stress is supposed to be so low that Coulomb's failure criterion has been passed (in red in Figure 1a) and hydro-frac is initiated. Fractures are then expected to be perpendicular to the minimal principal horizontal stress with little shear offset and a poor permeability enhancement without the use of proppant (Huenges and Ledru 2010). An alternative scenario would occur if no pre-existing fractures are optimally oriented. Then, a fresh new shear fracture might be initiated following the failure of Coulomb's criterion (in blue in Figure 1a, in the domain above Byerlee's criterion in blue) along the τ = c + σ tan φ line. In that case, the orientation of the cloud should be closer to the maximum principal stress direction. Indeed, the internal friction coefficient φ is typically larger than the friction coefficient of Byerlee's law. Such orientation of the seismic cloud has been observed for the very basal part of the stimulated reservoir suggesting that the initiation of a new fracture might have existed there. www.geothermal-energy-journal.com/content/2/1/14  (Vilarrasa et al. 2013)) Induced failure criteria related to a pore fluid increase: 1) Friction criterion for hydro-shear initiation along a pre-existing fracture (Byerlee criterion -green); 1') Fracture criterion for hydro-shear initiation in an intact rock (Mohr-Coulomb criterion -blue); and 2) Fracture criterion for hydro-frac intiation in an intact rock (Mohr-Coulomb criterion -red). (b) (modified after (Cornet et al. 2007)) Top view of the seismicity cloud during the September 1993 stimulation (25,300 m 3 injected) for two different flow rates (below and above 18 l/s at two different depths: below 3,000 m where the seismic cloud is oriented N146°E consistently with a hydro-shear regime and between 2,850 and 3000 m where hydro-fractures develop with a cloud oriented N179°E along the direction of the regional maximal horizontal stress σ H ). The angle between the orientation of the cloud during hydro-shear (1) and hydro-frac (2) defines the friction angle (φ ≈ 33°) (c) (modified after (Cornet et al. 2007)) Depth-time history of the seismicity and hydraulic flow during the 1993 hydraulic stimulation at Soultz-sous-Forêts. www.geothermal-energy-journal.com/content/2/1/14

Aseismic slip at Soultz-sous-Forêts
Aseismic slips in a deep geothermal reservoir were first observed at Soultz-sous-Forêts by (Cornet et al. 1997). They showed that after the two stimulations of GPK1 performed in September and October 1993, the borehole has been shifted at several depths just below the casing shoe at 2,857 m. From ultrasonic borehole images between 2,853 and 3,104 m before and after the two consecutive hydraulic injections of fresh water, cumulating up to 44,000 m 3 , they evidenced fresh displacements along existing fractures between 2,867 and 2,976 m that were of the order of several millimeters to several centimeters. One of the largest offset was measured at 2,925 m (fault F): 4.3 cm and considered as part of a significant fault zone of the reservoir with a quasi-vertical fault plane (dip of 86°and a dip direction of N48° (Cornet et al. 1997)) (see Figure 2). Interestingly, this fault zone is not part of the major fault zones observed during borehole logging (Sausse et al. 2010). Borehole logging shows that it is largely filled with secondary minerals like illite which might influence the friction properties of the fault (Genter et al. 1997). An important complementary observation is the seismicity recorded during the injection period corresponding to the occurrence of the borehole offset.  (Cornet et al. 1997) and (Bourouis and Bernard 2007) is added in red. It is located on the top of the reservoir. It intersects GPK1 at 2,925 m with a strike of N138°E and a dip of 86°NE. Seismicity is observed along this fault down to the major fault FZ-4770 (N234°E, 71°W). www.geothermal-energy-journal.com/content/2/1/14 Indeed, 165 events were recorded from the surface network which was far much less than the number of events recorded by the borehole down-hole network (20,000 events) (Helm 1996). They showed however that the magnitude of the largest recorded event was M L = 1.9 corresponding to a seismic moment M 0 of 5.4 × 10 11 Nm and a rupture area of typically d = 50 m for a stress drop of 9 MPa (assuming a circular crack embedded in an elastic isotropic medium and a shear modulus of μ = 20 GPa) consistently with observations in similar contexts (Abercrombie and Leary 1993). From general scaling laws (Kanamori and Anderson 1975), the expected slip D from this largest event would be of the order of 1 cm: D = 4M 0 /π μd 2 which is significantly smaller than the measured offset in the borehole. (Cornet et al. 1997) concluded that the large slips evidenced in the borehole could not be explained by the recorded earthquakes. Subsequently, they proposed that these large slips are aseismic which was the first indirect observation of aseismic slip induced by a fluid injection at Soultz-sous-Forêts EGS. (Bourouis and Bernard 2007) re-explored 10 years after these outstanding observations trying to conduct a fine analysis of the micro-seismicity during the injection period. They used a multiplet approach and relocated up to 400 events within the fault zone where aseismic slip was observed. The location accuracy was of the order of 1 m and they obtained 30 multiplets or families of similar events within the fault zone. An interesting observation is that the rupture size for all these events is of the order of d = 10 m. From their observations, they inferred three important conclusions. First, events within a multiplet were located within the same rupture zone which showed that the same asperity was reloaded and broken several times during the injection period. The second conclusion is that several asperities along the fault were ruptured at the same time during the injection. Third, the cumulative slip at each asperity through several ruptures was consistent with the borehole offset measurement. The conclusion is a clear image of the fault behavior during loading: the fault undergoes a large aseismic slip which triggers multiple local asperity failures.

Micro-seismicity during aseismic slip at Soultz-sous-Forêts
The goal of the present paper is to propose a mechanical model of this fault behavior. To do so, we developed an experimental approach to mimic the response of a single fault when submitted to a global slow rupture propagation but locally unstable at asperities where micro-seismic activity is triggered. The model is analogous and incorporates a large space and time dynamics: 7 orders of magnitude in time and 5 orders of magnitude in space. A numerical approach would be difficult with such a large range of timescales and wavelengths (Kaneko et al. 2010).

An analogous fault model with random asperities
Samples are made of transparent Plexiglas which provides optical access to the rupture propagation. The analogous fault model is obtained by annealing two plates of 20 × 10 × 1 cm 3 and 23 × 2.8 × 0.5 cm 3 at 190°C for 45 min which is significantly above the glass transition of the material but below the melting point. Under normal load, the two plates get in close contact and stick together along a relatively weak plane (weaker than the bulk) (Schmittbuhl and Måløy 1997). The goal is to study the collective behavior of multiple www.geothermal-energy-journal.com/content/2/1/14 asperities when a slow crack propagates and locally pins. For this, we sandblast one of the plates before the annealing procedure to induce random local toughness fluctuations (Lengliné et al. 2011b). Sandblasting is obtained using (180 to 300 μ) particles, an air pressure of 3 bars for 4 min at a blowing distance of 20 to 40 cm (Grob et al. 2009).

A subcritical rupture propagation
Samples with a weak plane along which the toughness is fluctuating but with an average toughness lower than that of the bulk are submitted to a slow cantilever loading (mode I) (see Figure 3) (Lengliné et al. 2012). This configuration differs from a shear mode (mode II or III) of a fault configuration, but similarities in the stress intensity formulation exist and analogies can be established to extend the results to other loading conditions (Schmittbuhl et al. 2003).
The loading includes different phases: firstly, the plate is loaded by the slowly imposed displacement of the loading rod (see left side of Figure 3a) at a low constant speed up to the initiation of the fracture. At some point, the displacement rate is suddenly increased (see Samples are made of two plates that are annealed together (20 × 10 × 1 cm 3 for the upper one and 23 × 2.8 × 0.5 cm 3 for the lower one). One of the plates is longer and undergoes a cantilever load from the vertical displacement of a rod (white circle on the left) that bends the plate. This bending with respect to the other plate induces a mode I fracture that propagates along the weak annealing interface (green zone). A camera (blue) is sitting over the sample and takes an optical image of the front propagation using a lighting by transmission located below the sample and reflected through a mirror (blue). Two acoustic arrays of sensors (red) are also attached to the fixed plate and record acoustic emissions during propagation. The sensors have a central frequency of 500 kHz. The array along the propagation direction is made of 64 sensors. The array perpendicular to the front is made of six sensors (see also Figure 7). (b) Bottom view of the sample. The rod is loading the longest plate of the sample. The green area corresponds to the broken zone behind the fracture front. The blue zone shows the zone where pictures are taken. The red zones are transmission views of the acoustic arrays. www.geothermal-energy-journal.com/content/2/1/14 Figure 4b), and the crack propagates over a few millimeters at quasi-constant velocity. The loading rate is of the order of 0.4 mm/s during this acceleration phase. The mechanical response, i.e., the force undergone by the plate during the experiment, is then recorded (Figure 4a). The setup being at imposed displacement, the crack propagation is stable and the average rupture velocity is controlled by the loading rate in the subcritical regime. Lastly, the loading rate is reduced, and the fracture propagates owing to the relaxation of the rupture at a decreasing speed (Lengliné et al. 2011a).

Optical monitoring of the fracture front
Taking images at high resolution (from 1 Mpixels with the fast camera up to 12 Mpixels with the SLR digital camera), we can obtain a precise spatial description of the fracture front (see Figure 5). It is observed that Fourier analysis of the front shows a power law behavior of the power spectrum over more than two decades with a roughness exponent H = 0.6 (Lengliné et al. 2011b;Schmittbuhl and Måløy 1997;Schmittbuhl et al. 2003).

Acoustic emissions during fracture propagation
From the acoustic sensors, acoustic emissions are recorded at 5 M samples/s (sampling time is 0.2 μs) for up to 30 s (cf. the duration of each experiment). Figure 6 shows an example of such an acoustic event. During the acoustic acquisition, the fast camera is shooting at a frame rate of 1,000 images/s (i.e., time delay between images is 1 ms) which is much lower than the sampling rate of the acoustic acquisition. As a consequence, recorded acoustic events have a much higher temporal resolution than optical events. In contrast, spatial resolution of optical events is much higher than that of acoustic events. Moreover, duration of acoustic events is also very different from those of optical events. The former never exceeded 100 μs. On the contrary, the latter could last up to several seconds at the largest scale. The overall observation is that numerous acoustic events are triggered during a single large-scale optical event. It suggests that we directly observe a simultaneous aseismic and seismic deformation along the crack front.
Locations of acoustic events are obtained using a beam-forming procedure. Because of the large number of waveforms (up to 64 sensors) and the continuous recording at very high frequency (5 MHz), we introduced an array-processing scheme for locating the source of the acoustic emissions. We define a grid domain over which we search for the possible location of the events. For each position of the grid we compute the expected arrival times of the event at each recorder assuming a homogeneous medium and straight ray path. We then shift each waveform according to these computed arrival times. Signals from all shifted waveform are stacked, and we finally keep the stacked signal around the first arrival. If the computed delay corresponds to the one for the true position of the source, we expect to have a very good match of all waveforms. We finally use the amplitude of this stack as a measure of the probability of the event location. An example of the location of one event is shown in Figure 7.  Figure 8 shows the evolution of the locations of the acoustic emissions during one experiment where the fracture front advances by 1 cm. We are able to follow the front advance during its propagation even if the spatial resolution of the acoustic location is significantly larger than the optical resolution. During each experiment, i.e., the propagation of the fracture front over several millimeters, we recorded several hundreds of acoustic emissions. Figure 9 shows a distribution of the energy recorded for each event. A first simple calibration provides an estimate of the absolute energy of the order of 10 −11 to 10 −8 J per events. Magnitude estimates of these events are then in the range −10 to −8 using the classical energy-magnitude relation from (Kanamori and Anderson 1975): M = 2/3Log(E) − 3.2. Figure 9 shows that a power law behavior exists over two decades of energy consistently with Gutenberg-Richter distribution.

Combining aseismic and seismic events
An imposed loading perturbation To go one step further in understanding the link between acoustic emission and local creep events, we took advantage of the experimental configuration to introduce a perturbation in the loading and look for the response of the system. Figure 4 shows how this perturbation is introduced in the loading rate. We observed that the average position of the fracture front follows exactly the perturbation of the loading frame after a small shift in time. The force decreases as the front velocity increases which is a velocity weakening effect. www.geothermal-energy-journal.com/content/2/1/14  Figure 4).
Following the crack front advance with the digital camera as described in the 'Methods' section, we are able to monitor the local speed of the front (Grob et al. 2009;Måløy et al. 2006). We see in Figure 10 that the overall trend of the experiment is imaged by the general evolution of the color: slightly blue at the beginning (long waiting time corresponding to low speed), becoming very red (shorter waiting time corresponding to higher speed), and back to mostly blue during the speed reduction of the front. The interesting part is the details of the local velocity. Clearly patches of increased speeds (red patches) exist. By thresholding at ten times the average velocity of the front, the image of the local velocity field, we obtained a set of 'optical' events that are significantly faster than the loading speed but also significantly slower than the dynamical rupture velocity (i.e., Rayleigh speed). These zones correspond to aseismic events which are accelerated advances compared to the average speed. An interesting observation is that local velocities are distributed as a unique power law on a broad range of velocities (more than 3 orders of magnitude): P(v) ∝ v −2.55 (Måløy et al. 2006), suggesting that most of the energy is dissipated by creep events consistently with observations from (Gross et al. 1993) who showed that only 3% of the fracture energy is radiated and from (Kanamori and Anderson 1975) for earthquake dynamics where the radiated energy is always negligible in front of the seismic moment (M 0 /E ≈ 1/20, 000, where M 0 is the seismic moment describing the total dissipated energy and E the radiated energy).
The last input from the experiment is the comparison of the optical events and acoustic emissions during the load perturbation. Figure 11 shows how both families of events are Comparison of the rate of optical events (i.e., slow events) and acoustic events (i.e. fast events) during the average front velocity perturbation shown in Figure 4, bottom panel. www.geothermal-energy-journal.com/content/2/1/14 synchronous which clearly suggests that acoustic emissions are directly controlled by the large-scale creeping process.

Conclusions
Observations at the Soultz-sous-Forêts EGS site of aseismic slips that are synchronous with micro-seismic events located in multiplets (Bourouis and Bernard 2007;Cornet et al. 1997) suggest that both processes co-exist within the same fault zones during a fluid injection period. Here we proposed an experimental model that directly supports this observation. Our model is built on the interfacial failure along a heterogeneous weak plane and combines creep failure and brittle rupture without fluid. Experimental observations are numerous and consistent with large-scale measurements. They provide clear hints on the processes involved at the asperity scale. Our conclusion is therefore that fluids are not necessarily the driving force of the fault activity in terms of pore pressure. They might have rather a role on local creep acceleration because of their effect on subcritical crack growth. An interesting perspective would be to estimate the evolution of the permeability owing to aseismic slip.