Healing Behavior of Simulated Fault Gouges From the Groningen Gas Field and Implications for Induced Fault Reactivation.

We investigated the frictional strength recovery (healing) and subsequent reactivation and slip-weakening behavior of simulated fault gouges derived from key stratigraphic units in the seismogenic Groningen gas field (N. E. Netherlands). Direct-shear, slide-hold-slide (SHS) experiments were performed at in situ conditions of 100 °C, 40 MPa effective normal stress and 10-15 MPa pore fluid pressure (synthetic formation brine). Sheared gouges were allowed to heal for periods up to 100 days before subsequent reshearing. The initial coefficient of (steady) sliding friction μ was highest in the Basal Zechstein caprock (μ = 0.65 ± 0.02) and Slochteren sandstone reservoir (μ = 0.61 ± 0.02) gouges, and the lowest in the Ten Boer claystone at the reservoir top (μ = 0.38 ± 0.01) and in the Carboniferous shale substrate (μ ≈ 0.45). Healing and subsequent reactivation led to a marked increase (∆μ) in (static) friction coefficient of up to ~0.16 in Basal Zechstein and ~0.07 in Slochteren sandstone gouges for the longest hold periods investigated, followed by a sharp strength drop (up to ~25%) and slip-weakening trajectory. By contrast, the Ten Boer and Carboniferous gouges showed virtually no healing or strength drop. Healing rates in the Basal Zechstein and Slochteren sandstone gouges were significantly affected by the stiffness of different machines used, in line with the Ruina slip law, and with a microphysical model for gouge healing. Our results point to marked stratigraphic variation in healed frictional strength and healing rate of faults in the Groningen system, and high seismogenic potential of healed faults cutting the reservoir and Basal Zechstein caprock units, upon reactivation. ©2020. The Authors.

Introduction

The Groningen gas field, located in the northeast of the Netherlands, is one of the largest onshore gas fields in the world. It was discovered in the 1950s and has been producing continuously since 1963. Since the 1990s, induced seismic events have been recorded with more or less steadily increasing frequency and magnitude, notably in the center of the field (van Thienen‐Visser & Breunese, 2015). The largest event to date was the magnitude M  = 3.6 Huizinge earthquake (Dost & Kraaijpoel, 2013), which occurred in 2012 and raised major public concern. It is now becoming increasingly clear that the events recorded tend to occur within or close to the depth interval of the reservoir (Spetzler & Dost, 2017) and are the result of reactivation of the numerous preexisting (mainly normal) faults that cross‐cut the sandstone reservoir and the directly overlying and underlying anhydrite/carbonate and clay/siltstone formations (Dost & Haak, 2007; Van Eijs et al., 2006). In order to evaluate if and when faults in the Groningen field might be reactivated as a result of gas production, and with what effects, it is important to investigate the frictional behavior of the fault rock materials present, and to evaluate how fault strength evolves (or heals) following both tectonic and induced slip. Aside from induced seismicity, investigating fault healing behavior is central to understanding the seismic cycle exhibited by active faults in tectonically deforming terrains (Dieterich, 1972b, 1972a; Marone, 1998b; Scholz, 2002), as postslip strengthening is a prerequisite for repetitive earthquake occurrence on the same fault.

Faults encountered in boreholes, or exposed at the surface by uplift and erosion, typically display a fault core characterized by narrow principal slip zones embedded in and filled with fine‐grained wear material, called fault gouge, derived from the surrounding rocks (Chester & Chester, 1998; Collettini et al., 2009; Fisher & Knipe, 1998; Sibson, 1986). From experiments on both natural and simulated fault gouges, it is well established that gouge‐filled faults regain strength, or “heal,” during laboratory simulations of the interseismic period, that is, during periods of zero imposed slip achieved in so‐called “slide‐hold‐slide” (SHS) experiments (e.g., Beeler et al., 1994; Ikari et al., 2016; Karner et al., 1997; Marone, 1998a; Olsen & Scholz, 1998). However, most of these laboratory experiments are of short duration compared to natural interseismic intervals, having “hold” durations of at most hours to a few days. In addition, healing behavior is commonly found to vary strongly between gouges of different mineral composition (e.g., Chen et al., 2015a; Carpenter et al., 2016) and to depend on stress, temperature, and pore fluid conditions (e.g., Bos & Spiers, 2002; Karner & Marone, 1998; Nakatani & Scholz, 2004; Tenthorey & Cox, 2006; Yasuhara et al., 2005). All of these factors make it difficult to extrapolate laboratory healing data to natural faults and to timescales that are relevant to induced or natural seismicity. With respect to the Groningen gas field, only limited healing data exist for fault rocks derived from relevant sedimentary lithologies, and no healing data exist for the specific lithologies that occur there.

To bridge these gaps, we performed direct‐shear SHS experiments on simulated fault gouges derived from the complete stratigraphic sequence of the Groningen reservoir system. In these experiments, initial shearing was followed by healing for prescribed hold durations ranging from 10 s to as long as 100 days (i.e., timescales which are of the same order as seasonal fluctuations in gas production from the field), and by subsequent reactivation of slip in a multistep sequence. The experiments were conducted using synthetic formation brine as pore fluid at simulated in situ reservoir conditions, that is, at 100 °C, 40 MPa effective normal stress, and pore pressures of 10–15 MPa. Two types of triaxial deformation machine were used, which differed in stiffness by a factor ~30. We analyzed the results in terms of time‐dependent frictional healing and slip weakening associated with reactivation, focusing on the effects of gouge composition as well as the stiffness of the experimental systems used. We compare our results with simulations performed using the classical rate‐and‐state friction (RSF) laws (Dieterich, 1979; Ruina, 1983), and we attempt to explain the observed healing behavior in terms of the active mechanisms, utilizing a microphysical model. Finally, we evaluate the implications for fault strength recovery and induced seismicity in the Groningen field. Besides the Groningen case, our results are also directly relevant to both natural and induced seismicity occurring on faults that cut sedimentary units of similar composition, such as the Permian sediments that make up many of the hydrocarbon reservoirs present in the Southern North Sea basin area (Glennie & Provan, 1990), or the anhydrite/carbonate sequences present in the Italian Apennines (Collettini et al., 2009) and Khuff gas fields in Qatar (Alsharhan & Nairn, 1994).

Experimental Methods

Sample Materials

Drill core and cuttings derived from the main stratigraphic units in the Groningen reservoir were obtained from the Stedum‐1 (SDM‐1) and Zeerijp‐2 (ZRP‐2) wells, located in the central to northern part of the Groningen field, courtesy of the Nederlandse Aardolie Maatschappij (NAM). The formations sampled were, from bottom to top, the Carboniferous shale/siltstone substrate, the (Permian) Upper Rotliegend Slochteren reservoir sandstone, the overlying Ten Boer claystone, and the Basal Zechstein anhydrite/carbonate evaporate caprock (see Hunfeld et al., 2017 for details). Simulated fault gouges were prepared from these lithologies by crushing and sieving the samples to obtain a grain size fraction <50 μm, simulating the fine‐grained wear material or fault gouge typically present in the fault core and principal slip zones of natural fault zones at shallow crustal depths (Sibson, 1986).

The mineral assemblages present in the gouges were determined using conventional X‐Ray Diffraction (XRD) analysis (see also Hunfeld et al., 2017). The Basal Zechstein caprock gouge consists of anhydrite (48 wt%), carbonates (32 wt% calcite and dolomite), barite (14 wt%), and minor amounts (<10 wt%) of quartz. The Ten Boer claystone material is composed predominantly of quartz (48 wt%), phyllosilicates (39 wt% in total, with 33 wt% muscovite/phengite and 6 wt% chlorite) and minor feldspar, dolomite, and hematite. The Slochteren sandstone gouge used here was derived from above the gas‐water contact and contains quartz (67 wt%), feldspar (19 wt%), and minor amounts of carbonates and phyllosilicates (kaolinite and phengite). The Carboniferous shale/siltstone samples investigated consist of quartz (55 wt%), phyllosilicates (35 wt% in total, with 21 wt% illite/chlorite, and 14 wt% muscovite/phengite), and feldspar (10 wt%).

All experiments were performed on samples saturated with a synthetic brine mimicking that present in the Groningen reservoir (in connate and pore‐saturating form). This brine contains 3.4 Molar (M) NaCl, 1.0 M CaCl2, 0.1 M MgCl2, and minor amounts (<0.1 M) of KCl, SrCl2, NasSO4, and BaCl2 and has a total ionic strength of 6.9 M.

SHS Experiments

Direct‐shear SHS experiments were performed on ~1 mm thick layers of simulated fault gouge using two types of triaxial deformation machines equipped with specially designed direct‐shear assemblies, located respectively at the High Pressure and Temperature Laboratory at Utrecht University (referred to as UU) and at the Geomechanics Laboratory at Shell Technology Centre Amsterdam (referred to as Shell). The triaxial “Heard” machine used at UU, is described by Hangx et al. (2010). The “Quadruple” machines used at Shell are similar to the triaxial machines described by Hol et al. (2018). Brief descriptions of both machines, including schematic diagrams, are included in the supporting information, Text S1. The key difference between the two machines is the way axial load is applied to the sample, which in the case of the Heard apparatus is done mechanically via a piston/yoke assembly driven by a motor/gearbox/ball‐screw system, whereas in the Quadruple machines, this is done hydraulically via a pressure‐driven piston. Due to this difference in machine design, the effective stiffness (K) in terms of shear stress increase on the sample per unit load point displacement was considerably higher in the Heard machine than in the Quadruple machines (~7 × 1010 Pa/m vs. 2.4 × 109 Pa/m, see supporting information, Text S1). One experiment in the Heard system (UU‐SS‐04) was run with an alternative loading configuration, replacing one of the elements in the internal loading column by a calibrated load washer spring (see supporting information, Figure S1), which reduced the effective stiffness to ~2 × 1010 Pa/m. The direct shear assemblies used in the two machine types were similar in design, measuring 1 inch in diameter in the Quadruple machines at Shell (holding a 25 mm wide × 24 mm long gouge sample), and 1.5 inch in diameter in the Heard apparatus at UU (35 mm wide × 47 mm long sample). See the supporting information (Text S1) for more details.

In the present experiments, shear displacement was imposed on the gouge layer by advancing the loading piston in the Heard/Quadruple machines. In the Heard machine, this was done via the motor/gear box/ball‐screw system, driving the yoke and loading piston at a constant velocity of ~1 μm/s. In the Quadruple machines, this was achieved by advancing the hydraulic driving‐pump piston (i.e., the load point) at a fixed displacement rate of ~7 μm/s, which translates into a displacement rate at the sample of ~1 μm/s. In both systems, sliding of the gouge layers was imposed at this rate until a steady‐state shear strength was reached. To investigate the healing or strength recovery behavior of the gouges, we then initiated a series of SHS steps (Dieterich, 1972b), whereby active shearing of the gouge layer is stopped (load‐point displacement rate is set to zero) for a predetermined amount of time (“hold” interval) and active shearing subsequently resumed, while continuously monitoring shear strength. The hold times (t hold) in our experiments varied from 10 s up to 100 days (8.7 × 106 s). Most experiments included multiple hold intervals, with hold durations that increased in steps of approximately 0.5 orders of magnitude (here referred to as a “normal” SHS sequence). Two runs, on Basal Zechstein and Slochteren sandstone material (runs Shell‐19 and Shell‐15, respectively), were performed using a “reverse” SHS sequence, that is, where t hold was successively decreased. One experiment on Ten Boer claystone gouge (run Shell‐12) employed a random SHS sequence. In runs where multiple hold intervals were performed, sliding displacements of 0.3–0.6 mm were applied between the holds. A total of 19 experiments were performed (11 at Shell and 8 at UU) on samples of simulated Basal Zechstein (5), Ten Boer claystone (4), Slochteren sandstone (9), and Carboniferous shale/siltstone (1 run at UU only) gouges. All experiments were performed at 100 °C, using 50–55 MPa confining/radial pressure (P ) and 10–15 MPa pore pressure (P ), so that the effective normal stress, defined as = P − P , was equal to 40 MPa in all tests. Table 1 lists all experiments performed, detailing the experimental conditions and SHS sequences imposed.

Table 1

List of Experiments, Conditions, and Key Results

Experiment no.	Testing machine	T (°C)	P c/P rad (MPa)	P f (MPa)	Pore fluid	Hold intervals (t hold) (s)	μ 1mm (mm)	D tot (mm)	μfinal (mm)
Basal Zechstein
UU‐BZ‐01	Heard	100	55	15	brine	10–30–100–300–1000–3000	0.633	4.11	0.586
UU‐BZ‐02	Heard	100	55	15	brine	10–30–100–300–1000–3000–10000–60000	0.648	4.76	0.550
Shell‐14	Qdr.	100	50	10	brine	3600–10800–72000–255600–860160–2672640	0.636	4.00	0.444
Shell‐19	Qdr.	100	50	10	brine	588300–10200–3600	0.669	2.95	0.439
Shell‐21	Qdr.	100	50	10	brine	10–30–100–300–1000–3000–64800	0.636	3.77	0.428
Ten Boer claystone
UU‐TB‐01	Heard	100	55	15	brine	10–30–100–300–1000–3000	0.379	4.23	0.366
Shell‐11	Qdr.	100	50	10	brine	7862400	0.376	2.09	0.390
Shell‐12	Qdr.	100	50	10	brine	3600–72780–336420–2678400–935220–10800	0.386	3.64	0.348
Shell‐23	Qdr.	100	50	10	brine	10–30–100–300–1000–3000–10800–76320	0.392	3.81	0.315
Slochteren sandstone
UU‐SS‐01	Heard	100	55	15	brine	10–30–100–300–1000–3000	0.594	4.17	0.598
UU‐SS‐02	Heard	100	50	10	brine	10–30–100–300–1000–3000–10000–45120–100800–270000–956160–3002400	0.603	6.10	0.573
UU‐SS‐03	Heard	100	50	10	brine	10–30–100–300–1000–3000–10000–60000–100000–335800	0.644	4.77	0.592
UU‐SS‐04	Heard b	100	50	10	brine	10–30–100–300–1000–3150–10200–54420–176400–3286800	0.619	7.00	0.561
Shell‐08	Qdr.	100	50	10	brine	3600–10800–82800–334800–864000–2419200	0.622	3.72	0.546
Shell‐10	Qdr.	100	50	10	brine	8712000	0.591	2.08	0.487
Shell‐15	Qdr.	100	50	10	brine	1037640–429000–85920–10800–3780	0.610	3.57	0.431
Shell‐20	Qdr.	100	50	10	brine	10–30–100–300–1000–3000–65580	0.590	3.29	0.527
Shell‐25	Qdr.	100	50	10	brine	10–30–100–300–1000–3000–10800–78060	0.603	3.87	0.361
Carboniferous shale/silts
UU‐C‐01	Heard	100	55	15	brine	10–30–100–300–1000–3000	0.445	4.18	0.450

Symbols used: Qdr. denotes use of the quadruple cells at Shell (low stiffness); Heard denotes use of the heard apparatus at UU (high stiffness); T = temperature; P /P rad = confining/radial pressure; P  = pore fluid pressure; t hold = hold time; μ 1mm and μ final = steady‐state friction coefficient, given at 1–1.5 mm and at maximum shear displacement (D tot), respectively.

Experiment UU‐SS‐04 was performed using the alternative, intermediate‐stiffness heard configuration

After axial unloading, depressurization/cooling and extraction of the sample assembly from the vessel, intact sample fragments were recovered for microstructural analyses. These were flushed with demineralized water, dried at 50 °C and then impregnated with epoxy resin. Sections of the impregnated fragments were cut in an orientation parallel to the shear direction and perpendicular to the shear plane and were examined by means of scanning electron microscopy (SEM), performed using an FEI Helios 650UC FIB‐SEM. Imaging was done in backscattered electron (BSE) mode, using an acceleration voltage of 5–15 kV.

Data Acquisition and Analysis

During the experiments, axial load, axial piston displacement, confining pressure, temperature, and pore fluid pressure signals were recorded at rates of 1–1,000 Hz using either a National Instruments, 16‐bit A/D data acquisition system and Signal Express VI logger software (Heard system) or custom‐built acquisition systems and recording software developed at Shell (Quadruple systems). The raw data were processed to obtain shear stress τ (MPa), effective normal stress σ eff (MPa) and apparent friction coefficient μ = τ/σ eff (−), versus displacement and time. Calibration runs using a steel dummy of known elastic properties were used to correct the recorded displacement data for elastic distortion of the testing machines. Note that in the Quadruple machines, axial displacement is recorded directly across the sample assembly, in addition to the load point displacement measured in the axial drive pump. The former is used in the stiffness calibration to correct the axial displacement data for elastic distortion. The maximum absolute measurement error due purely to transducer resolution was ±1.5 μm (Heard) and ±0.1 μm (Quadruple) in shear displacement, and ±0.003 in friction coefficient in both machines. However, due to the lack of an internal load cell in the Quadruple machines, shear stress could not be determined independently of seal friction acting on the axial loading piston. Under the current experimental conditions, seal friction amounts to ~0.7–0.8 MPa in axial stress on the sample or an equivalent of 0.02 in friction coefficient, for which our data were corrected. During SHS sequences, however, small variations in seal friction can occur, producing uncertainties in healing and strength evolution data of the order of 10% of the total seal friction. This amounts to ±0.002 in friction coefficient, which is similar to the uncertainty due to transducer resolution. Fluctuations in temperature and confining/pore fluid pressure (hence effective normal stress) were minor (<0.2 °C and <0.05 bar) in all experiments performed, including runs imposing long holds (up to 101 days, Table 1).

A number of parameters are used here to characterize the observed SHS behavior, following Karner et al. (1997). These are defined in Figure 1. During hold periods, the measured shear stress relaxes nonlinearly from a steady‐state value, due to sample creep and elastic unloading of the testing machine and sample assembly. Stress relaxation and accompanying shear displacement occurring during holds are denoted as Δμ hold and d hold, as shown in Figure 1. When shearing is resumed, shear stress increases (near) linearly due to (largely) elastic loading until the sample begins to slip, and a peak in friction is reached, at which turning point the slip velocity is momentarily equal to the load point velocity. After the peak, slip accelerates transiently, as friction decreases rapidly and the machine unloads, until a new steady state is established. The stress drop (Δτ) is defined as the difference in peak stress and the minimum shear stress supported following reactivation. It occurs over the slip‐weakening distance (d ), as shown in Figure 1. Frictional healing (Δμ) is defined as the difference between the peak friction value and the steady‐state value before the onset of the hold (see Figure 1), and is often found to increase log‐linearly with hold duration (e.g., Dieterich, 1972b). This type of log‐linear healing behavior is commonly expressed as where β is a constant known as the healing rate (typically of the order of 0.01/decade; Marone, 1998b), and t is a cutoff time of the order of 1–100 s, beyond which Δμ conforms to a log‐linear relation with t hold. We quantified the frictional healing rates in our experiments by fitting equation 1 to Δμ versus log t hold data, using least‐squares regression to obtain the best‐fit healing rate (β) and cutoff time (t ). In obtaining Δμ, we did not correct for any background (linear) displacement weakening or hardening trends, as these have virtually no influence on healing as defined in Figure 1.

Figure 1

Example of a slide‐hold‐slide (SHS) test and parameter definition (run Shell‐08, Slochteren sandstone gouge), showing the evolution of shear stress and friction coefficient μ (blue line, left‐hand axes) and shear displacement (orange line, right‐hand axis) versus time. Frictional relaxation and displacement occurring during a hold (t hold), are denoted as and d hold, respectively. Frictional healing is defined as the difference in peak strength upon reactivation and the prehold value, and is denoted . After the peak, friction or shear strength decreases rapidly (strength drop denoted as ) over a slip‐weakening distance d , before attaining a new steady‐state value.

Results

Friction Coefficient Versus Displacement Curves

Representative curves of μ versus shear displacement for all gouge materials tested using the Heard machine at UU and the Quadruple machines at Shell are shown in Figures 2a and 2b, respectively. Figure 2c shows the Quadruple runs (Shell‐15 and Shell‐19) imposing a reversed SHS sequence on Slochteren and Basal Zechstein gouges. At displacements up to 0.2–0.5 mm, all experiments showed a rapid, near‐linear increase in μ. In most gouges, this was followed by hardening at decreasing rate toward a steady‐state frictional strength, attained at 1–1.5 mm of shear displacement, that persisted through subsequent SHS steps with only minor background slip weakening or slip hardening (Table 1). In the case of Basal Zechstein samples, however, a peak frictional strength was usually attained at 0.5–1 mm displacement, followed by near‐linear slip‐weakening superimposed on the SHS steps, as illustrated by the difference in μ values (0.10 to 0.15) obtained at ~1 mm shear displacement versus maximum shear displacement (Table 1). Total shear displacements of 2.1–4.6 mm and 3.8–7.0 mm were reached in the Shell and UU experiments, respectively (Table 1). Ten Boer and Carboniferous gouges sometimes showed minor step changes (<0.01) in the steady‐state friction coefficient between SHS's (see Figure 2a). In the case of Basal Zechstein gouges in Shell experiments, we observed a more systematic series of steps in the background slip‐weakening trend, with a notable increase in friction level per step (up to ~0.02) following the 10‐day and 31‐day hold intervals (see Figure 2b).

Figure 2

Representative friction coefficient μ versus shear displacement curves obtained in slide‐hold‐slide experiments performed at (a) Utrecht University (UU) using the high‐stiffness heard machine, and (b) and (c) Shell, using the low‐stiffness quadruple machines. All curves are corrected for elastic machine distortion effects. Hold durations are indicated in seconds in (a) and in hours or days in (b) and (c). The experiments shown in (c) imposed a reversed SHS sequence (i.e., decreasing hold durations in successive holds).

The (near‐) steady‐state friction coefficient values (μ 1mm) attained at 1–1.5 mm shear displacement, that is, before any hold intervals were applied, are plotted for all experiments in Figure 3 (see also Table 1). Basal Zechstein gouges exhibited the highest friction coefficients (μ 1mm = 0.65 ± 0.02), followed by Slochteren sandstone gouges (μ 1mm = 0.61 ± 0.02). Ten Boer claystone and Carboniferous shale gouges were significantly weaker (μ 1mm = 0.39 ± 0.02 and ~0.45, respectively). For any given sample type, the values obtained in the experiments performed at UU (including using the reduced stiffness configuration ‐ UU‐SS‐04), were very similar to those obtained at Shell, using the lower‐stiffness Quadruple machines. Indeed, no differences in initial steady‐state friction coefficient could be detected, beyond the spread of ±0.02 obtained in repeat experiments performed using the same material and experimental setup (typically).

Figure 3

The (near‐) steady‐state friction coefficients obtained at 1–1.5 mm displacement for all experiments on simulated fault gouge derived from the four principal stratigraphic units present in the Groningen gas field (see also Table 1).

SHS Data

General Features

Regardless of machine configuration, we observed two types of SHS behavior in the present experiments superimposed on the background slip‐weakening trends, as seen in Figure 2. The first type, exhibited by Slochteren sandstone and Basal Zechstein gouges, is characterized by a distinct peak in frictional strength upon reshear after a hold interval, followed by slip weakening until a new steady‐state value. The amount of frictional healing (Δμ) generally increased with increasing hold duration for these materials, resembling the classical behavior observed by Dieterich (1972b). The second type is exhibited by Ten Boer and Carboniferous gouges. In these materials, no significant peak friction or slip weakening occurred after the hold intervals, and frictional healing (Δμ) was negligible.

Aside from these main types of behavior, systematic differences in SHS results were obtained using the two different testing machines, notably in the Slochteren sandstone and Basal Zechstein gouges. These are illustrated in Figure 4, where we compare data obtained in short holds of ~1 and ~3 hr, and in subsequent resliding, for Slochteren and Basal Zechstein gouges tested in both machines. In the low‐stiffness experiments (Shell), samples typically showed less frictional relaxation (Δμ hold) than in the high‐stiffness experiments (UU). At the same time, the shear displacement accumulated during the hold (d hold) was consistently larger in the Shell experiments. In addition, during the postpeak weakening stage in experiments on these gouges, the Shell experiments typically exhibited a transient low in frictional strength before evolving to a new steady‐state value, which did not occur in the UU experiments (Figure 4). Similar differences in Δμ hold and d hold were observed between machine types for the Ten Boer and Carboniferous samples, but no transient low in strength occurred during the reslides, as no (dynamic) slip‐weakening occurred in these gouges.

Figure 4

Representative examples of individual slide‐hold‐slides performed in Shell and UU experiments on Slochteren sandstone and basal Zechstein gouges. (a) Friction coefficient versus time. (b) Shear displacement versus time. (c) Friction coefficient versus load‐point displacement for the reslide stages only. Note that friction, shear displacement, and time are all normalized to the onset of the hold. Note that the load point velocity during sliding stages in the Shell setups is ~7 μm/s (displacement rate of the moving piston within the hydraulic source pump generating axial pressure) versus 1 μm/s in the UU setup (displacement rate of the yoke, see supporting information, Text S1). Therefore, the difference in slope during reslides when plotting friction versus time is only ~4, instead of ~30 as is expected based on the difference in effective stiffness between the two machines. Plotting friction versus load point displacement, as is done in (c), shows this stiffness contrast between the machines. During steady‐state sliding, the shear velocity is ~1 μm/s in all cases, as indicated in (b).

Frictional Healing

In Figure 5 (refer also Table 1), we plot healing magnitude Δμ versus t hold for all SHS runs performed on Basal Zechstein (Figure 5a), Slochteren sandstone (Figure 5b), and Ten Boer claystone plus Carboniferous shale samples (Figure 5c). For descriptive purposes, the data for each material and each experimental machine/configuration are fitted using equation 1, yielding the best‐fit values for healing rate (β) and cutoff time (t ) indicated.

Figure 5

Frictional healing versus hold time t hold for (a) basal Zechstein, (b) Slochteren sandstone, and (c) Ten Boer claystone plus carboniferous shale gouges. Note the different y scale in (a) versus (b and c). The data sets, including both normal and reverse data in the case of Shell experiments, are fitted with equation 1 (solid lines), and the best‐fit healing rate (ß) and cutoff time (t ) are shown. In the case of basal Zechstein gouges, fits to standard and reverse SHS sequence results are shown separately (dotted and dashed lines, respectively).

Both the Basal Zechstein and Slochteren sandstone samples showed a more or less log‐linear increase in Δμ with t hold, for hold durations longer than 100–1000 s (Basal Zechstein) and ~100 s (Slochteren sandstone). This behavior is well described by equation 1, as shown by the best‐fit solid lines in Figures 5a and 5b, except for the Basal Zechstein in the low‐stiffness (Shell) experiments, which showed exceptionally large scatter in the data and more complex healing behavior. In both materials, we observed clear effects of machine stiffness on healing rate and cutoff time. The healing rates in Slochteren sandstone gouges were 0.009 and 0.012 in the high‐stiffness and intermediate‐stiffness (UU) experiments, respectively, and 0.015 in the low‐stiffness (Shell) experiments. β values thus increased by a factor of ~1.5 for an approximately thirtyfold (i.e., 1 order of magnitude) decrease in effective stiffness. The cutoff time, below which no log‐linear relation between Δμ and t hold was observed, correspondingly increased from 54 s (high stiffness) to 85 s (intermediate stiffness) and then 136 s (low stiffness). Basal Zechstein gouges showed a similar inverse dependence of healing rate and cutoff time on stiffness, with β ranging from 0.02 in the high‐stiffness case to 0.033 in the low‐stiffness case (when considering all data combined), and t ranging from 96 to 830 s, respectively. For hold durations beyond ~1 hr imposed in any given machine type, Basal Zechstein gouges thus healed roughly twice as fast as Slochteren sandstone gouges.

In sharp contrast to Basal Zechstein and Slochteren sandstone gouges, the Ten Boer claystone and Carboniferous shale/siltstone samples showed negligible healing, with Δμ ≈ 0 for all hold durations investigated (Figure 5c). For completeness, β was nonetheless determined for these materials, although a cutoff time was not considered here because of the lack of log‐linear healing in these materials. The relation β = d(Δμ)/d(log (t hold)) was accordingly fitted to the data, instead of equation 1. For Ten Boer samples, this yielded β = 0 in the Shell experiments (maximum hold duration 91 days) and β = −0.001 in the UU experiments (maximum hold 3000 s). The single Carboniferous sample tested at UU showed a just‐positive healing rate (β = 0.0004; maximum hold 3000 s).

Note that in Figure 5, healing rates are determined from undifferentiated data sets that include experiments in which both normal and reverse SHS sequences were applied (see Table 1). However, in the case of the Basal Zechstein gouges, the data show additional trends related to the type of sequence imposed. In experiment Shell‐14, in which we imposed a normal SHS sequence with holds ranging from 1 hr up to 1 month (blue, solid data points), Δμ remained more or less constant from step to step, even decreasing slightly for hold durations longer than ~3 days (2.6 × 105 s), suggesting a saturated healing magnitude of Δμ = 0.08 to 0.09. By contrast, the first and longest hold (~7 days or 5.9 × 105 s) applied in reverse SHS experiment Shell‐19 on the Basal Zechstein gouge (Figure 5a, blue triangles, see also Figure 2c), showed a much higher Δμ value (~0.16). When considering the results obtained in run Shell‐19 only, the fitted healing rate β was 0.104 (dashed blue line), which is much higher than the value obtained from all Shell data for this material combined (solid blue line). Specifically, the healing rate for Basal Zechstein gouge obtained using only the normal SHS sequence results of experiments Shell‐14 and Shell‐21 was 0.03 (dotted blue line fitted to the square data points, Figure 5a). These results indicate a strain or hold history dependence of healing in Basal Zechstein gouges, whereby holds imposed early during the experiments, that is, at low shear strains or when no previous holds have been applied, produce more healing than holds of the same duration applied at larger shear strains or after a series of earlier holds (e.g., compare the ~1 hr and ~20 hr hold results in experiments Shell‐14 and Shell‐21). In the tests on Slochteren sandstone and Ten Boer claystone gouges, on the other hand, no significant differences were observed between data obtained in normal and reverse or random SHS sequences. Data obtained in these sequences are therefore not differentiated in Figure 5.

Slip Weakening in Basal Zechstein and Slochteren Sandstone Samples

We now consider the behavior observed in samples that showed measurable time‐dependent frictional healing and slip weakening upon resliding, that is, the Slochteren sandstone and Basal Zechstein gouges. Specifically, we focus on the relation between stress drop (Δτ) and weakening distance (d ) during the postpeak weakening phase, that is, when friction evolves from the peak value to a new steady‐state value via slip‐weakening (refer Figure 1).

The Δτ versus d data obtained for the Basal Zechstein and Slochteren sandstone gouges, using the low‐ (Shell) and high‐stiffness (UU) machines, are plotted in Figure 6a. The majority of the data show stress drops in the range 0–5 MPa, with corresponding weakening distances up to ~0.3 mm. Note that in the high‐stiffness UU results, Δτ and d do not exceed values of ~2.3 MPa and 0.16 mm, respectively, and are approximately half those seen in the low‐stiffness Shell results. The most extreme values were obtained in experiment Shell‐19 on Basal Zechstein gouge following a 5.9 × 105 s (~7 days) hold, which resulted in a stress drop of ~8 MPa and weakening distance of ~0.45 mm. Nonetheless, all data represented in Figure 6a fall approximately on a straight line, demonstrating that stress drop tends to be linearly related to weakening distance, independently of gouge composition (quartz‐dominated Slochteren sandstone versus anhydrite/carbonate‐dominated Basal Zechstein) or of the stiffness of the experimental fault plus testing machine. This relation is characterized by a slope of 17.4 MPa/mm, as indicated by the linear fit to the data shown in Figure 6a.

Figure 6

(a) Stress drop versus weakening distance d obtained in experiments on Slochteren sandstone and basal Zechstein gouges, using the high‐stiffness (UU) and low‐stiffness (Shell) setups. (b, c, and d) Examples of the near‐linear slip‐weakening behavior observed during reslides in basal Zechstein and Slochteren sandstone gouges, plotting shear stress (black line, left‐hand axis) and time normalized to the onset of the reslide (blue line, right‐hand axis) versus shear displacement. Subscripts p and r denote peak and residual shear stress, respectively.

Representative examples of the evolution of shear stress with shear displacement during postheal resliding steps imposed on the Slochteren sandstone and Basal Zechstein gouges are shown in Figures 6b–6d. These data demonstrate near‐linear, rapid slip weakening in all cases, although the stress drop and weakening distance in the high‐stiffness UU example (Figure 6c) is much less than seen in the Shell results obtained after similar hold times (Figure 6b, also Figure 6d). From the corresponding displacement versus time data, the peak velocities reached during the dynamic weakening stage were of the order of ~1 cm/s in the Shell experiments, and ~1 mm/s in the UU experiments. Closely similar behavior was seen in all Slochteren sandstone and Basal Zechstein samples per machine type.

Microstructures

After the experiments, the samples broke into cohesive, platy fragments that split on planes running parallel to the shear direction and close to the sample‐piston interface, which are interpreted as boundary (B) shears, and on R 1 Riedel shears (terminology after Logan et al., 1979). Representative microstructures from the deformed Slochteren sandstone, Basal Zechstein, and Ten Boer gouges are shown in Figure 7. The microstructures obtained in experiments conducted at Shell versus UU were very similar for each of these materials but distinctly different between materials. Sample material was not recovered from the single experiment on Carboniferous shale/siltstone.

Figure 7

BSE images showing representative microstructures of sheared Slochteren sandstone (a and b), basal Zechstein (c and d), and Ten Boer claystone (e and f). Shear sense is right lateral in all images. (a) Part of the microstructure of run Shell‐10, showing localization into an R 1 Riedel shear band. (b) Enlargement of the bulk gouge shown in (a). (c) Part of the microstructure of run Shell‐14, showing sharp localization into a boundary (b) shear at the bottom of the sample, as well as an R 1 Riedel shear band. (d) High‐magnification image taken in the bulk gouge of (c), which shows indented grain contacts (white arrows) between anhydrite and carbonate grains, and overgrowths (black arrows) in adjacent pore space. (e) Typical microstructure of Ten Boer claystone gouge (Shell‐12). (f) High‐magnification image showing anastomosing phyllosilicate minerals curving around quartz grains, defining a P and R 1 foliation.

The Slochteren sandstone samples (Figures 7a and 7b) show a gouge bulk with relatively large (10–40 μm), subrounded to angular quartz and feldspar grains surrounded by a very fine (<1 μm) grained matrix. Microfractures are common in the larger quartz and feldspar grains, whereas the finer grains show fewer fractures. Where visible, grain‐to‐grain contacts are sharp, but most are separated by the very fine, submicron‐sized material. Grain indentations or truncations were not observed in the larger grains. The resolution of the SEM images was insufficient to inspect such features in the fine‐grained material. The shear bands inferred from sample splitting are visible as boundary and R 1 shears in the SEM, characterized by a reduced internal grain size (<1 μm up to ~10 μm). The width of boundary shears was hard to determine, since the samples split on these bands so they were not fully preserved. R 1 Riedel shears were found to be spaced roughly 3–5 mm apart, typically 30–100 μm in thickness and characterized by irregular, convoluted boundaries.

In the Basal Zechstein samples (Figures 7c and 7d), the gouge bulk consists of subrounded to angular grains varying in size from ~5 to 50 μm, surrounded by a matrix of grains smaller than 1 μm. Anhydrite, barite, and carbonate grains are often fractured, whereas fracturing is less common in the quartz grains. Barite grains (white contrast in Figure 7c) in particular show marked comminution and often form fine‐grained trains in the P, R 1, and Y shear band orientations, both throughout the gouge bulk and within or close to shear bands. At higher magnifications, anhydrite and carbonate grains in the gouge bulk (Figure 7d) often show indented contacts (white arrows) associated with fine idiomorphic crystallites suggestive of anhydrite or carbonate overgrowths (black arrows) in adjacent pores. The dominant shear bands in Basal Zechstein samples are R 1 and B shears, with occasional Y shears that link up with R 1 shears. All are characterized by a strongly reduced grain size (below 1 μm) with respect to the bulk gouge and are typically 10–50 μm in thickness. The B shears are especially straight and are developed along the entire length of the samples, at either one or both boundaries.

The Ten Boer claystone samples (Figures 7e and 7f)) show a matrix of aligned, very fine grained phyllosilicate platelets of submicrometre to micrometer length scale, that anastomose and wrap around larger, angular to subrounded grains of mostly quartz. These quartz grains vary in size from ~5 to 50 μm and show little evidence of grain‐scale fracturing or comminution. The phyllosilicates define foliations in both the P and R 1 orientation (Figure 7f). Numerous fractures in the P, R 1, and Y orientations are visible throughout the samples, suggesting the presence of shear planes that opened following unloading. Other than these decompression features, there is no evidence for localization features or grain size reduction indicating the presence of shear bands. Instead, shear deformation seems to have been accommodated on the phyllosilicate foliation, in a more or less distributed fashion throughout the full width of the gouge layer.

Discussion

The present results demonstrate clear differences in frictional strength, healing and reactivation behavior of simulated fault gouges derived from the main stratigraphic units of the Groningen gas field. In terms of steady‐state friction coefficient (Figure 3), the highest values were found in the Basal Zechstein caprock (μ 1mm = 0.65 ± 0.02), and the lowest in the Ten Boer claystone (μ 1mm = 0.39 ± 0.02). The Slochteren sandstone reservoir and Carboniferous substrate showed intermediate values (μ 1mm = 0.61 ± 0.02 and ~0.45, respectively). These results are very similar to the values obtained by Hunfeld et al. (2017) for the same gouge compositions and highlight the marked mechanical stratigraphy that the key formations in the field exhibit, with respect to fault gouge frictional strength (Figure 3).

Our SHS experiments further revealed that the quartz‐rich Slochteren sandstone and the anhydrite‐carbonate‐dominated Basal Zechstein show considerable time‐dependent frictional healing (Δμ up to ~0.07 in Slochteren and ~0.16 in Basal Zechstein samples, see Figure 5), leading to significant stress drops (up to 5 MPa in Slochteren sandstone, and 8 MPa in Basal Zechstein) and associated slip weakening upon reactivation (see Figure 6). Additionally, we observed an inverse dependence of healing rate (), in these materials, on the stiffness of the experimental set up. By contrast, the phyllosilicate‐rich Ten Boer claystone and Carboniferous shale/siltstone showed negligible healing and negligible stress drop.

In the following, we further analyze and aim to explain the effects observed in our experiments. We first compare our healing data with previous experimental work. We go on to compare our healing data for the Basal Zechstein and Slochteren sandstone gouges with rate‐and‐state friction simulations, in an attempt to elucidate the observed effects of machine stiffness. A microphysical explanation is then proposed for the observed healing and stiffness effects. Finally, we discuss the implications of our results in the context of induced seismicity in the Groningen field.

Comparison With Previous Data

Influence of Gouge Composition on Frictional Healing

Carpenter et al. (2016) investigated frictional healing behavior of a wide range of natural and synthetic fault gouges in SHS experiments performed at room temperature and 100% relative humidity. Fitting equation 1 to their data, they found that of all materials investigated, phyllosilicate‐rich fault gouges exhibit the lowest healing rates (β = 0.001–0.003). Intermediate healing rates were found in quartz‐ and feldspar‐rich gouges (β ~ 0.007), and the highest healing rates were found in carbonate‐rich gouges (β ~ 0.012). Other experimental studies also showed that the presence of phyllosilicates in halite and silicate gouges can significantly inhibit frictional healing, especially in the presence of water (Bos & Spiers, 2000; Carpenter et al., 2011, 2012; Ikari et al., 2012), due to prevention of solution transfer effects and contact cementation/healing as a result of “saturation” of the contact area (i.e., cessation of further growth in contact area) (Tesei et al., 2012).

Under hydrothermal conditions, healing rates in quartz‐ or carbonate‐dominated gouges are generally found to be faster than at room temperature (Karner et al., 1997; Nakatani & Scholz, 2004), whereas for phyllosilicate‐rich gouges, they remain low (Chen et al., 2015a). Healing rates in water‐saturated carbonate fault gouges, at 80–150 °C, was found to be significantly faster than in clay‐bearing gouges (β = 0.013–0.017 versus 0.004–0.007; Chen et al., 2015a). Pluymakers and Niemeijer (2015) reported healing rates in water‐saturated anhydrite gouges that were similar to those reported by Chen et al. for carbonate gouges, under similar P–T conditions. In both studies, pressure solution most likely dominated compaction in the gouge layers during hold intervals, reducing gouge porosity and thus introducing a transient strengthening effect during resliding (see also Chen & Spiers, 2016). Nakatani and Scholz (2004) investigated frictional healing in quartz gouge at temperatures of 100–200 °C and found healing rates that were considerably higher than at room temperature due to the operation of solution transfer processes (pressure solution). Other studies on pure quartz gouges under hydrothermal conditions (up to 927 °C) have shown similar results (Chester & Higgs, 1992; Karner et al., 1997; Tenthorey et al., 2003; Tenthorey & Cox, 2006; Yasuhara et al., 2005). The consensus from these studies is that thermally activated pressure solution significantly enhances healing rates in experimental faults and fault gouges, with an increase in friction coefficient per decade increase in hold time of up to ~0.1.

Our data agree well with all of the above studies, exhibiting virtually no healing in the phyllosilicate‐rich Ten Boer and Carboniferous samples, intermediate healing rates in the quartz‐rich Slochteren sandstone, and the highest healing rates in the anhydrite‐carbonate rich Basal Zechstein samples (Figure 5). The healing rates obtained in the Slochteren and Zechstein gouges are similar to those reported by Nakatani and Scholz (2004) for quartz gouges, and by Chen et al. (2015a) and Pluymakers and Niemeijer (2015) for carbonate/anhydrite gouges, under hydrothermal conditions where pressure solution is likely active. The microstructures in our sheared Basal Zechstein samples also show clear evidence for solution transfer processes in the form of grain indentations and overgrowths. The fact that gouges made up of carbonate and sulfate minerals heal faster than quartz‐rich gouges can be attributed to the faster kinetics of pressure solution in these relatively soluble ionic minerals (see Spiers et al., 2004, and references therein).

The steps observed in background slip‐weakening trends and the superimposed increase in μ following long holds in our Basal Zechstein data (Figure 2b) are similar to the “non‐Dieterich” type healing behavior observed by Chen et al. (2015b) in carbonate gouges (see also Karner et al., 1997; Muhuri et al., 2003). Such changes in (slip‐weakening) friction level after versus before hold intervals have been attributed to effects of healing on slip localization (Sleep et al., 2000) or to permanent changes in sample porosity and/or grain size (Chen et al., 2015a; Chen & Spiers, 2016). Background slip weakening in general was more pronounced in the 1 inch Shell setup (Figure 2). It is associated with ongoing shear localization (Kanagawa et al., 2000; Niemeijer, Spiers, & Peach, 2008; see also Pluymakers et al., 2014), but can also be due to a decrease in load‐bearing area with increasing displacement in a direct shear setup, which affects the smaller setup used at Shell more at any given displacement (i.e., larger relative change in area). This is considered in more detail in the supporting information (Text S2). However, these strain‐dependent effects do not significantly influence healing data, as shear displacements between the onset of a hold and the peak strength upon reactivation are small.

Influence of Stiffness on Frictional Healing

Several previous SHS studies have demonstrated that frictional healing in gouges is faster when the shear stress on the sample is allowed to relax or is kept constant during hold periods, as opposed to conditions where the shear stress is reduced to zero in the beginning of the hold period (Bos & Spiers, 2002; Karner & Marone, 1998; Nakatani & Mochizuki, 1996). In experiments at room temperature and humidity, performed by Nakatani and Mochizuki (1996) on granitic gouge, at fixed servo‐controlled hold stresses (shear stress), healing rate was found to increase with increasing hold stress—which is inconsistent with constitutive equations, such as equation 1 (Dieterich, 1978) for purely time‐dependent healing. Marone and Saffer (2015) modeled SHS cycles with the classical RSF laws, using both the Dieterich aging law and the Ruina slip law to describe the evolution of the RSF state parameter . They found that, while neither model showed complete agreement between the experimental observations and simulation results, the Ruina law generally matched the observations much better. All of these results suggest that the accrual of sliding (or shear) displacement during hold periods promotes frictional strengthening upon reshear. Since more shear displacement occurs during hold periods executed in a low‐stiffness setup (e.g., Figure 4), healing should therefore increase with decreasing machine stiffness. Our results for Basal Zechstein and Slochteren sandstone are consistent with this effect (Figures 5a and 5b), regardless of the underlying physical mechanism.

Interestingly, Beeler et al. (1994) found no effect of stiffness on healing behavior in SHS experiments performed on initially bare granite and quartzite blocks using a rotary shear apparatus where stiffness was varied. Their results imply that healing was independent of fault slip during hold periods and is better described by the time‐dependent Dieterich law. However, Bhattacharya et al. (2017) recently reevaluated the Beeler et al. data, pointing out that while no effect of stiffness on healing was observed, high‐stiffness experiments showed systematically more relaxation than low‐stiffness experiments. These authors modeled the experiments reported by Beeler et al. (1994) using the classical RSF models and concluded that the SHS behavior observed was much better captured by the slip‐dependent Ruina law rather than the time‐dependent Dieterich law and that the lack of an effect of stiffness on healing could be reproduced just as well by the slip law for a specific set RSF parameters.

Present Healing Data Versus Rate‐and‐State Friction Simulations

We now follow the previous studies mentioned above by comparing our healing results for Basal Zechstein and Slochteren sandstone gouges with model simulations conducted using the classical RSF equations. We use the most common form of the RSF equations (Dieterich, 1979; Ruina, 1983), given which describes frictional strength as a function of two variables: (1) slip rate (V) and (2) “state” (), recalling that the latter is generally viewed as the average lifespan of a population of grain‐to‐grain contacts in the case of gouge friction (Marone, 1998b). In equation 2, a is the “direct effect” parameter accounting for the variations in frictional strength μ due to changes in slip rate, and b is the “evolution” parameter which determines the change in friction due to evolution of fault state over the slip distance d . The parameters and represent the values of friction coefficient and slip rate at an arbitrary reference steady‐state. The evolution of the state parameter is formulated using the “aging law” (Dieterich, 1979) or the “slip law” (Ruina, 1983) referred to in the preceding section. These are respectively given by In the aging law 3a, evolves linearly with time for stationary contacts. In the slip law 3b, state evolution occurs only for nonzero slip velocities.

We performed forward simulations using equations 2, (3a) and (3b), coupled with the relation d which describes the interaction with the elastic surroundings, where k is the stiffness of the loading frame, normalized by the normal stress, and V lp is the load point velocity. In the simulations, we used fixed values for the RSF parameters a, b, and d , based on velocity‐stepping friction experiments published by Hunfeld et al. (2017) (using average values obtained from velocity steps of 0.1–1 μm/s), which were performed on the same gouge compositions, and at similar pressure, temperature, and pore fluid conditions as the present experiments. From the simulations, Δμ is calculated as a function of hold time. In Figures 8 and 9, we plot the experimentally obtained healing data for the Slochteren sandstone and Basal Zechstein gouges, alongside the simulations performed for these materials using a stiffness of 2.4, 20, and 70 × 109 Pa/m—corresponding to the effective stiffnesses of the apparatus configurations used in the present experiments.

Figure 8

Comparison between healing data and RSF model simulations for Slochteren sandstone gouges in (a) low‐stiffness Shell experiments, (b) intermediate‐stiffness, and (c) high‐stiffness UU experiments. The RSF model simulations are shown separately in (d), where “low K” equals 2.4 × 109 Pa/m, “intermediate K” equals 2.0 × 1010 Pa/m, and “high K” equals 7.0 × 1010 Pa/m. the values for a, b, and d used are based on Hunfeld et al. (2017).

Figure 9

Comparison between healing data and RSF model simulations for basal Zechstein gouge in (a) low‐stiffness Shell experiments and (b) high‐stiffness UU experiments. The RSF model simulations are shown separately in (c), where “low K” equals 2.4 × 109 Pa/m and “high K” equals 7.0 × 1010 Pa/m. The values for a, b, and d used are based on Hunfeld et al. (2017).

From Figures 8 and 9, it is evident that the simulations performed with the Dieterich aging law show a systematic increase in healing with increased stiffness for a given hold time, with no change in healing rate (slope). Aside from the increase in cutoff time seen with decreasing stiffness, this trend is opposite to that seen in our experiments, in which decreasing stiffness produced increased healing (rates) beyond the cutoff time in both Slochteren sandstone and Basal Zechstein gouges. In contrast to the Dieterich simulations, the results obtained using the Ruina slip law not only show the same trends as our experimental data but also show a good fit to the data for both materials, especially for the high‐ and intermediate stiffness runs at UU experiments. Despite the scatter in the experimental data obtained in the low‐stiffness runs (Shell), the Dieterich simulations also show a reasonable match to these data, at least for hold durations <105 s. For longer hold durations, the Dieterich simulations overestimate healing.

From this comparison between experiments and RSF simulations (Figures 8 and 9), it is clear that only the Ruina slip law captures all trends in the healing data observed in our experiments on Basal Zechstein and Slochteren sandstone gouges. The observed stiffness effect thus points to a slip‐dependence of healing in these gouges, as captured phenomenologically by the Ruina/slip law. In the following, we offer a possible explanation for the physical origin of this effect.

Microphysical Basis for Frictional Healing and the Effect of Stiffness

The increase in static rock friction with time, seen in SHS experiments, is classically interpreted to originate from an increase in real contact area as a result of thermally activated creep processes operating on the interface contact asperities (e.g., Brechet & Estrin, 1994; Dieterich & Kilgore, 1994). Specifically in gouges, the effect has been attributed to grain contact “quality” or area growth through the operation of pressure solution (Bos & Spiers, 2002; Nakatani & Scholz, 2004; Niemeijer, Marone, & Elsworth, 2008; Yasuhara et al., 2005), or due to compaction by pressure solution, resulting in a decrease in gouge porosity and an increase in dilatation angle and hence frictional strength (Chen & Spiers, 2016; Niemeijer & Spiers, 2007). Alongside compaction, cohesion development during quasi‐static holds has been proposed to contribute to gouge restrengthening, based on experiments performed on halite gouges (Bos & Spiers, 2002; Hickman & Evans, 1992; Muhuri et al., 2003; van den Ende & Niemeijer, 2019). Current models for healing by pressure solution (Chen & Spiers, 2016; Niemeijer, Marone, & Elsworth, 2008; Yasuhara et al., 2005) assume that the gouge deforms as an isotropic medium with compaction by pressure solution during hold periods being driven solely by the effective normal stress acting on the gouge layer. In these models, healing is accordingly unaffected by shear stress and by different shear stress relaxation paths (and shear strain development paths) occurring during hold periods implemented in systems (machines) with different stiffness. For this reason, current models do not capture the effect of stiffness on healing (rates) seen in our experiments.

Building upon the existing microphysical model for gouge friction and healing developed by Chen, Niemeijer, and Spiers (from here on referred to as the Chen‐Niemeijer‐Spiers or CNS model, see Niemeijer & Spiers, 2007; Chen & Spiers, 2016; Chen et al., 2017), and drawing inspiration from the pressure solution microstructures observed in our Basal Zechstein experiments (see Figure 7d), we propose that healing in those experiments is the result of (1) compaction by pressure solution, causing a reduction in gouge porosity and an associated increase in dilatation angle, compare Chen and Spiers (2016), and (2) an increase in grain contact strength, or cohesion, due to overgrowths that develop as a result of mass removal from grain contacts by pressure solution and subsequent precipitation/cementation in the adjacent pore space. This is illustrated schematically in Figure 10, where we adopt the same idealized gouge microstructure as Chen and Spiers (2016), including the same geometric relation between porosity and grain‐contact/dilatation angle. Since the contact‐normal stress that drives mass removal (pressure dissolution) at grain contacts is a function of the externally supported normal and shear stresses, and since shear stress relaxation is affected by stiffness, the extent of mass transport and grain contact overgrowth/cementation occurring during hold periods, and the corresponding cohesion developed, should depend on stiffness as well as hold time. By incorporating this grain‐contact cohesion development into the CNS model, we test if the observed stiffness effect can be reproduced. We restrict ourselves to the Basal Zechstein gouges here, as these show the largest effects, and the pressure solution parameters are well constrained for the dominant anhydrite and calcite phases at the conditions of our experiments (unlike quartz, for which they are only constrained at much higher temperatures). The aim here is to qualitatively assess the impact of stiffness on healing using time‐dependent grain‐contact cohesion in the CNS model, thereby gaining insight into the potential physical cause(s) for this effect, as opposed to a detailed quantitative comparison with our experimental data.

Figure 10

Schematic illustration of the microstructure/geometry assumed in the CNS model and geometric relations (refer Chen & Spiers, 2016). The inset shows indented grain contacts and overgrowths. A is the total contact area including overgrowths and A the indented grain contact area where grains are dissolved. The area of the overgrowth annulus is denoted A . The total volume of overgrowths adding to A and thus cohesion is less than the total volume of dissolved material (dark gray) at the real grain contact, as some material is still in solution in the grain boundary fluid film or is precipitated elsewhere in open pore space. This is represented by weighting factor f < 1 in equation 6.

In the CNS model, the frictional strength of a gouge is assumed to be controlled by competition between (1) dilatant granular shear flow, with velocity‐strengthening atomic scale frictional interactions at grain contacts, and (2) time‐dependent shear and compaction by a thermally activated creep mechanism (such as linear viscous pressure solution) acting at the grain scale. The response of the gouge to this creep mechanism is further assumed to be isotropic, so that time‐dependent compaction and shear can be treated independently. The model has been shown to successfully reproduce the behavior seen in velocity‐stepping and SHS friction tests performed on fluid‐saturated carbonate or anhydrite‐dominated gouges using a single testing machine with fixed stiffness (Chen & Spiers, 2016; Hunfeld et al., 2019).

To capture the effects of varying stiffness described above, in the simplest possible way, we extend the CNS model (as formulated by Chen & Spiers, 2016) by assuming that, during granular shear flow of a deforming gouge layer with partial accommodation by pressure solution, individual grain contacts satisfy a Coulomb‐type slip criterion given Here, is a term representing cohesion developing at grain contacts due to overgrowths produced by pressure solution, is the logarithmically rate‐dependent grain boundary friction coefficient employed by Chen and Spiers (2016), and is the (intensified) normal stress acting at grain contacts. Following Niemeijer and Spiers (2007), shear resistance to pure granular flow is then written as where the first term on the right‐hand side represents the contribution of cohesion to shear resistance, and the second the contribution of both friction and dilatant slip on inclined contacts. In this equation, k and H are dimensionless constants related to the geometry of the grain pack (see Niemeijer & Spiers, 2007), ψ is the dilatation angle that is linked in the assumed microstructure to gouge porosity, φ, via the relation tanψ = H(q − 2φ) (Figure 10), and is the effective normal stress acting on the gouge layer. During rapid granular flow, grain contact cohesion will be negligible, as grain rolling and neighbor swapping will largely prevent overgrowths enveloping and cementing grain contacts. During hold periods, on the other hand, the rate of granular shear flow and of grain neighbor swapping in the gouge layer will quickly diminish due to arrest of the loading piston. At that point, gouge deformation will fast become dominated by pressure solution creep operating at grain contacts in both the shear and normal (compaction) directions, causing ongoing shear stress relaxation to levels that preclude grain boundary frictional slip due to its much higher (exponential) sensitivity to stress (see Chen & Spiers, 2016). Cohesion at grain contacts will then increase, with hold time, in proportion to the amount of mass removed from grain contacts and deposited around them as load‐supporting cement. To represent this in the CNS model, we assume that, during an individual hold period, cohesion increases with growing cement‐rimmed contact area or annulus A according to the relation where and A 0 are reference values of cohesion and small overgrowth area, corresponding to the prehold sliding period. The growth in A with increasing pressure solution strain , measured from the onset of the hold in the direction normal to the grain contact, is then approximated via which gives In these equations, f is a dimensionless factor prescribing the contribution of mass dissolved at grain contacts to load‐supporting overgrowth area (as opposed to mass precipitated along free grain boundaries away from grain contacts, producing little contribution to cohesion). Based on our microstructures (Figure 7d), we assume f to be of the order of 0.5 (as opposed to of the order of 0.01, as expected for uniform overgrowths on pore walls, e.g., Yasuhara et al., 2005), thus assuming that half of the mass is precipitated around the “neck” of the contact (see Figure 10).

We implemented the relations outlined above in the CNS model. Specifically, equation 5 replaces the cohesion‐free “friction law” in the model (see equation (32c) in the paper by Chen & Spiers, 2016). The amount of pressure solution strain () occurring normal to the inclined grain contacts during an individual hold period, and hence the corresponding change in grain contact area and cohesion via equations 6–8, was obtained using the geometrical relation where and are the total pressure solution strains accumulated during an individual hold period in the direction normal and parallel to the shear plane, respectively (see Chen & Spiers, 2016).

We simulated SHS behavior using the CNS model with this simple grain‐contact cohesion model, using the input parameters listed in Table S1 presented in the supporting information. The values for the parameters we used are based on the conditions applied in our experiments, the microstructures observed in our Basal Zechstein gouges, and the closely similar pressure solution parameters for anhydrite and calcite from the literature, which together constitute ~80 wt% of the total mineral composition in the Basal Zechstein gouge material. For the model implementation and solution method used to simulate healing versus log (t hold) in SHS simulations, see Chen and Spiers (2016).

In Figure 11a, we show versus hold duration as calculated from the modified CNS model, varying system stiffness in the range 0.2 × 109 to 7 × 1011 Pa/m (i.e., an order of magnitude above and below the range of machine stiffness in our experiments). The trends are very similar to the healing behavior observed in our experiments, showing both an increase in cutoff time and healing rate with decreasing stiffness. In the model, this is the result of a higher shear stress being maintained during the hold periods in the low‐stiffness experiments (due to reduced shear stress relaxation compared with high‐stiffness tests). These higher shear stresses (i) prolong the influence of stress sensitive intergranular friction and hence granular shear flow to longer hold durations, which counteracts compaction and thus increases the apparent cutoff time for healing, and (ii) increase the grain contact‐normal stress, speeding up pressure solution and cohesion development, and thus healing rate for longer hold durations. These two competing effects result in a healing versus log (t hold) curve that initially shows sluggish healing for the low stiffness case, and then accelerates to cross the high‐stiffness healing curve (Figure 11a). Note that no crossover occurs for the lowest stiffness simulated in Figure 11a, for which pressure solution plus cohesion development outweigh granular flow from the onset of healing. The position of the crossover, the magnitude of healing, and the contribution of cohesion strengthening strongly depend on the kinetics of pressure solution (varying DCS for a diffusion‐controlled case, see Table S1), as shown in Figure 11b. In Figure 11c, we show simulation results in which the product of S 0 and f are varied or set to zero, also producing significant effects, with higher values of S 0 or f speeding up healing. Note that when the reference cohesion S 0 is set to zero (inhibiting any cohesion effect) the model is identical to the standard CNS model formulation, and no stiffness effect occurs apart from the change in cutoff time related to shear stress relaxation affecting granular flow and thus porosity evolution during the early stages of a hold.

Figure 11

CNS with grain‐contact cohesion simulations representing basal Zechstein gouges (see Table S1 in the supporting information for the input parameters). (a) Healing versus hold time for varying stiffness, including values corresponding to the present experiments (2.4 × 109 Pa/m for the Shell setup and 70 × 109 Pa/m for the UU setup). Decreasing stiffness results in (i) longer cutoff time due to prolonged influence of granular flow, and (ii) steeper healing rate due to faster development of cohesion (increased mass transfer) plus compaction and associated increase in dilatation angle. For the lowest stiffness investigated (2.4 × 108 Pa/m), cohesion‐strengthening dominates over granular flow from the onset of healing. (band c) Sensitivity of healing to the rate of pressure solution (via DCS) and the product of reference cohesion S times weighting factor f. increased pressure solution rates result in faster healing and a decrease in cutoff time. Increasing S or f has similar effects. The end‐member case where there is no cohesion effect at all ( = 0) is shown in (c). In all simulations, DCS = 1 × 10−24 m3/s and = 2.5 MPa, unless otherwise specified.

Overall, the behavior observed in our experiments on the Basal Zechstein gouge is qualitatively captured by the extended CNS model. In the case of materials that deform by fluid‐assisted creep mechanisms like pressure solution, as is likely the case for sulfates, carbonates and quartz under conditions pertaining to the shallow upper crust, the model thus offers a physical explanation for the apparent “slip dependence” of healing as phenomenologically captured by the Ruina slip law.

One aspect of the experimental healing data that we report for the Basal Zechstein gouge that is not captured by the extended CNS model is the plateau in healing observed in the Basal Zechstein gouge experiments when subjected to a normal SHS sequence reaching hold periods of 3 days and beyond (Shell‐14, Figure 5a). The behavior suggests a limit to fault healing at a saturation level of Δμ = 0.08–0.09 for this material. However, our reverse SHS experiment (Shell‐19) showed almost twice as much healing as the normal sequence (Shell‐14) for a ~7 day hold. This implies that the observed plateau does not reflect the maximum strength gain in Basal Zechstein gouge. No such healing plateau has been observed in previous SHS experimental work, although it should be noted that the present study is the first to investigate hold durations of up to 3 months (8.7 × 106 s) under hydrothermal conditions. One exception is the study of Chen et al. (2015a) on carbonate‐rich gouges healed under hydrothermal conditions at 80 °C. This did show a saturation or plateau in healing for hold durations beyond 104 s up to the maximum investigated hold duration of 105 s.

We now consider the possible reasons for the observed plateau in our Basal Zechstein healing data. While only the study of Chen et al. (2015a) on calcite has shown similar behavior, a few other studies investigating healing (e.g., in halite gouges; Niemeijer, Marone, & Elsworth, 2008) have reported higher healing rates in reverse SHS sequences than in normal sequences (cf. Figures 2c and 5). This was interpreted to be consistent with compaction by pressure solution as the primary strengthening mechanism, as explained by Niemeijer, Marone, and Elsworth (2008) (see also Chen et al., 2015a). In short, in a reversed SHS sequence, the longest holds are applied when the gouge still has relatively high porosity and small grain‐contact area, so the gouge compacts and heals more than in a similar hold duration in a normal SHS sequence, where the gouge has already compacted during earlier, shorter holds. This mechanism can also explain the large scatter observed in Shell data for Basal Zechstein gouges, specifically the difference in healing observed for ~1 hr hold durations imposed in normal SHS sequence experiments Shell‐14 and Shell‐21(see Figure 5a). In the case of Shell‐21, the 1 hr hold was preceded by shorter holds and produced less healing compared to the 1 hr hold in Shell‐14, which was the first hold applied and produced much more healing.

Following this interpretation, the plateau observed in our Basal Zechstein healing data (Shell‐14, Figure 5a) may thus be the result of cumulative gouge compaction, producing low starting porosity at the onset of successive holds. According to the CNS model, such a decrease in porosity between successive sliding intervals should lead to systematically increasing frictional strength (Chen & Spiers, 2016), as is observed for long holds in our Basal Zechstein experiments in the Shell setup (see Figure 2b). Regarding the limit of healing of gouge‐filled faults, in the CNS model, a reduction in porosity from 30% to a near‐zero level leads to an increase in frictional strength of the order of Δμ ≈ 0.3 (Chen et al., 2020), which is much higher than the plateau observed in Figure 5a. Development of cohesion alongside compaction, as proposed in the previous section, would raise this upper limit even further, potentially even up to the level of intact rock strength (Tenthorey & Cox, 2006). Another potential mechanism capable of producing a plateau in healing is related to a local increase in pore fluid pressure (i.e., undrained behavior) during long hold periods, leading to a decrease in effective normal stress and hence slower compaction, as well as lower peak strength during the reactivation stage.

In light of the above, and recalling that the reverse SHS sequence produced significantly more healing than the plateau values seen in the normal SHS sequence results, we conclude that it is unlikely that fault restrengthening in Basal Zechstein gouges is limited to the increase in (static) friction of 0.08–0.09 suggested by the data in Figure 5a.

Implications

The present results demonstrate that the marked mechanical stratigraphy of faults cutting the Groningen reservoir system, previously recognized in terms of a contrast in (near) steady‐state gouge friction coefficients for different lithologies (Figure 3, see also Hunfeld et al., 2017), is enhanced by gouge healing (see Figure 5). The faults in the Groningen reservoir were tectonically active during Late Jurassic to Early Cretaceous rifting, followed by reactivation during the Late Cretaceous to Early Tertiary inversion that accompanied the Alpine orogeny (De Jager & Visser, 2017). Since then, these faults have likely undergone considerable strength recovery via healing, prior to reactivation in response to gas production in the past few decades. Our experiments show that fault segments cutting through and incorporating gouges derived from the Slochteren sandstone reservoir rocks and Basal Zechstein caprock are expected to have restrengthened the most, whereas segments dominated by Ten Boer‐ and Carboniferous‐derived gouges likely remain relatively weak. Extrapolating the healing data for Slochteren sandstone shown in Figure 5b from 3 months to 1, 10, and 100 M years, and assuming an initial friction coefficient of 0.61, yields a total healed (static) strength of up to 0.78, 0.795, and 0.81, respectively. Similarly, extrapolating the Basal Zechstein results (Figure 5a) to these healing times, with 0.65 as initial friction, gives values up to 1.0, 1.035, and 1.07.

Since reactivation of healed Basal Zechstein and Slochteren sandstone gouges is accompanied by a stress drop and approximately linear, rapid slip‐weakening behavior (see Figure 6), we infer significant potential for seismogenesis by slip weakening upon reactivation of faults cutting these formations (e.g., Andrews, 1976; Ida, 1972; Ohnaka, 2013). By contrast, the Ten Boer and Carboniferous gouges showed no healing or strength drop. Together with the velocity‐strengthening behavior of the Ten Boer and Carboniferous gouges (Hunfeld et al., 2017), the data suggest that these phyllosilicate‐rich gouges exhibit stable sliding behavior with no tendency for accelerating slip, whereas the Slochteren sandstone and Basal Zechstein gouges are much more prone to nucleating unstable slip via slip weakening, potentially aided by velocity‐weakening in the case of Basal Zechstein material or its mixtures with sandstone gouge [Hunfeld et al., 2017, 2019. Note here that the seismogenesis‐prone nature of the Basal Zechstein and Slochteren sandstone applies not only to first reactivation of long inactive faults but also to repeated reactivation due to gas production, as substantial healing already occurs within 3 months (see Figure 5)—which is rapid on production timescales.

The most accurate earthquake locations in the Groningen field (Dost et al., 2017; Spetzler & Dost, 2017) show that the majority of the events occur on faults within the reservoir‐caprock interval, which is consistent with our findings. However, the volume of data acquired since installation of a new (deep and shallow) borehole seismometer network in the center of the Groningen field is as yet insufficient to establish if there are correlations with a specific lithology or combination of juxtaposed lithologies (i.e., whether events occur preferentially in the Basal Zechstein, the Ten Boer claystone or in the Slochteren sandstone).

Besides the Groningen case, the present results are also applicable to faults cutting sedimentary sequences of similar composition elsewhere. The high friction coefficients and rapid healing behavior observed for anhydrite/carbonate rich gouges in particular imply that faults cutting such formations can produce large stress drops as well as repeated events nucleating on the same fault (Scholz, 2002). Indeed, moderate to large natural earthquakes commonly occur on faults cutting evaporite and carbonate sequences, for example, the 2009 M  = 6.1 L'Aquila earthquake (Chiarabba et al., 2009), the 2008 M  = 7.9 Wenchuan earthquake (Zhang et al., 2010), and the 2009 M 5.0–5.4 earthquakes in Albania (Ormeni et al., 2013).

In addition to the above aspects related to the mineralogical composition of faults, the present SHS experiments have revealed that the elastic properties of the fault surroundings, that is, the machine in our case, play a significant role in post‐slip strength recovery, with healing being enhanced by low stiffness and associated shearing during load relaxation. In modeling frictional strength recovery after natural and induced fault activation, the implication of our data for gouge types that heal is that the effective shear stiffness of the host rocks, which is controlled by their bulk elastic properties and the volume of the elastically deformed zone, should be considered. Similarly, stress evolution on faults, either during tectonic reloading of seismically active faults or production‐induced stressing in reservoir systems, is expected to affect healing rate via gouge compaction and cohesion development. Finally, the effective stiffness may also influence dynamic fault strength during reactivation by causing a larger dynamic strength drop than expected from healing alone (as occasionally observed in our experiments, e.g., Figures 2b and 4; see also Bhattacharya et al., 2017).

Conclusions

We performed direct shear, SHS friction experiments on simulated fault gouges prepared from core and drill cuttings recovered from the main stratigraphic units present in the reservoir system of the seismogenic Groningen gas field. We investigated the frictional strength, postslip strength recovery, reactivation, and subsequent slip behavior of these materials employing hold durations in the range of 10 s up to 100 days. The experiments were done at simulated reservoir conditions, that is, at temperatures of 100 °C, effective normal stress of 40 MPa, and using a 6.9 M ionic strength mixed chloride brine, mimicking the formation water present in the Groningen reservoir. Our main findings are summarized as follows:

In reverse stratigraphic order, that is, from caprock to substrate, our results in terms of (near‐) steady‐state friction coefficient in initial stages of sliding (displacements of 1–1.5 mm) show the highest values in the Basal Zechstein caprock gouges (μ = 0.65 ± 0.02), and the lowest in the underlying Ten Boer claystone (μ = 0.39 ± 0.02). These are followed by the Slochteren sandstone reservoir gouges (μ = 0.61 ± 0.02) and the Carboniferous shale substrate (μ ≈ 0.45).

Slochteren sandstone and Basal Zechstein gouges showed significant time‐dependent healing (Δμ up to ~11% and ~25%) and stress drop following reactivation. By contrast, no healing or strength drop was observed in Ten Boer claystone and Carboniferous shale gouges.

Healing rates in Slochteren sandstone and Basal Zechstein gouges increased by a factor ~2 for a thirtyfold decrease in effective shear stiffness of the (experimental) fault surroundings, along with a systematic increase in cutoff time for log‐linear healing. These trends are well captured by the Ruina/slip formulation of rate‐and‐state friction, implying that healing is purely slip dependent.

Based on microstructural observations showing evidence for solution transfer processes operating during healing of Basal Zechstein gouges, and making use of the Chen‐Niemeijer‐Spiers model for gouge friction (Chen et al., 2017; Chen & Spiers, 2016), we present a simple physical explanation for the observed healing behavior and stiffness effect seen in the Basal Zechstein data. In our model, healing is the result of combined gouge layer compaction and cohesion development at grain contacts due to mass transport and cementation. We anticipate similar processes in the Slochteren sandstone, but insufficient data are available to model the relevant mass transport processes under the temperature conditions of the Groningen reservoir.

With respect to induced seismicity in the Groningen gas field, our results point to (i) an enhancement, via healing, of the mechanical stratigraphy seen in steady‐state friction values, and (ii) significant potential for seismogenesis via slip‐weakening in healed gouges derived from the Basal Zechstein caprock and Slochteren sandstone reservoir units.

Supporting Information S1

Click here for additional data file.