Groningen Gas Field Research Articles

Incorporating dwelling mounds into induced seismic risk analysis for the Groningen gas field in the Netherlands

In order to inform decision-making regarding measures to mitigate the impact of induced seismicity in the Groningen gas field in the Netherlands, a comprehensive seismic risk model has been developed. Starting with gas production scenarios and the consequent reservoir compaction, the model generates synthetic earthquake catalogues which are deployed in Monte Carlo analyses, predicting ground motions at a buried reference rock horizon that are combined with nonlinear amplification factors to estimate response spectral accelerations at the surface. These motions are combined with fragility functions defined for the exposed buildings throughout the region to estimate damage levels, which in turn are transformed to risk in terms of injury through consequence functions. Several older and potentially vulnerable buildings are located on dwelling mounds that were constructed from soils and organic material as a flood defence. These anthropogenic structures are not included in the soil profile models used to develop the amplification factors and hence their influence has not been included in the risk analyses to date. To address this gap in the model, concerted studies have been identified to characterize the dwelling mounds. These include new shear-wave velocity measurements that have enabled dynamic site response analyses to determine the modification of ground shaking due to the presence of the mound. A scheme has then been developed to incorporate the dwelling mounds into the risk calculations, which included an assessment of whether the soil-structure interaction effects for buildings founded on the mounds required modification of the seismic fragility functions.

Compaction of the Groningen Gas Reservoir Sandstone: Discrete Element Modeling Using Microphysically Based Grain‐Scale Interaction Laws

AbstractReservoir compaction, surface subsidence, and induced seismicity are often associated with prolonged hydrocarbon production. Recent experiments conducted on the Groningen gas field's Slochteren sandstone reservoir rock, at in‐situ conditions, have shown that compaction involves both poroelastic strain and time independent, permanent strain, caused by consolidation and shear of clay films coating the sandstone grains, with grain failure occurring at higher stresses. To model compaction of the reservoir in space and time, numerical approaches, such as the Discrete Element Method (DEM), populated with realistic grain‐scale mechanisms are needed. We developed a new particle‐interaction law (contact model) for classic DEM to explicitly account for the experimentally observed mechanisms of nonlinear elasticity, intergranular clay film deformation, and grain breakage. It was calibrated against both hydrostatic and conventional triaxial compression experiments and validated against an independent set of pore pressure depletion experiments conducted under uniaxial strain conditions, using a range of sample porosities, grain size distributions, and clay contents. The model obtained was used to predict compaction of the Groningen reservoir. These results were compared with field measurements of in‐situ compaction and matched favorably, within field measurement uncertainties. The new model allows systematic investigation of the effects of mineralogy, microstructure, boundary conditions, and loading path on compaction behavior of the reservoir. It also offers a means of generating a data bank suitable for developing generalized constitutive models and for predicting reservoir response to different scenarios of gas extraction, or of fluid injection for stabilization or storage purposes.

Full-Wavefield, Full-Domain Deterministic Modeling of Shallow Low-Magnitude Events for Improving Regional Ground-Motion Predictions

ABSTRACT At the fundamental level, seismic risk analysis relies on good modeling tools for predicting the ground motion resulting from hypothetical earthquake events, which is traditionally approximated using many variations of ground-motion prediction equations (GMPEs). The main benefit of these equations lies in their low computational cost, allowing one to run Monte Carlo simulations in which event probabilities are dictated by regional catalogs comprising historical observations. These equations, however, rely on approximations that are only accurate in a statistical sense. In this study, we consider cases in which regional high-resolution 3D earth models are available from exploration reflection seismology. These high-fidelity velocity models allow us to perform deterministic elastic ground-motion simulations at local distances, given a prescribed synthetic earthquake event, and compare the results with those predicted by GMPEs. This full-wavefield full-domain modeling approach is significantly more costly and particularly challenging due to the slow shear-wave velocity at the near surface, which requires fine spatial and temporal discretizations. With the aid of powerful computational resources, we use an adaptive mesh generator and an efficient wave solver to model the 3D elastic and anelastic wave propagation from the hypocenter all the way to the ground surface. This approach can simultaneously account for 3D subsurface structures, near-surface site effects, topographic relief, and the radiation pattern of the source. In areas where observations are sparse, the modeling results can fill the gap between stations and provide a test bed that can be used for improving the development and accuracy of GMPEs. This approach is well suited for areas where shallow low-magnitude-induced seismic events can occur. Lastly, to demonstrate our approach, we consider an observed seismic event at the Groningen gas field and compare the recorded ground motions with both—those predicted by our approach and those predicted by GMPEs.

In-Reservoir Waveform Retrieval for Monitoring at Groningen—Seismic Interferometry with Active and Passive Deep Borehole Data

The Groningen gas field in the Netherlands is an ideal test bed for in-situ reservoir monitoring techniques because of the availability of both active and passive in-reservoir seismic data. In this study, we use deconvolution interferometry to estimate the reflection and transmission response using active and passive borehole data within the reservoir at ∼3-km depth and separate up- and downgoing P- and S-wave fields by f-k filtering. We validate the results using synthetic data of a 1D elastic model built from sonic logs recorded in the well. The estimated full-waveform reflection response for a virtual source at the top geophone is consistent with the synthetic response. For the virtual source at the bottom geophone, the reflection response appears to be phase delayed, though its arrivals are consistent with the local subsurface geology. Similarly, the first-order estimated local transmission response successfully approximates that of the P-wave velocity in the reservoir. The study shows that reliable subsurface information can be obtained from borehole interferometry without detailed knowledge of the medium parameters. In addition, the method could be used for passive reservoir monitoring to detect velocity, attenuation, and/or interface time-lapse variations.

Quantification of damage evolution in masonry walls subjected to induced seismicity

This paper aims to quantify the evolution of damage in masonry walls under induced seismicity. A damage index equation, which is a function of the evolution of shear slippage and opening of the mortar joints, as well as of the drift ratio of masonry walls, was proposed herein. Initially, a dataset of experimental tests from in-plane quasi-static and cyclic tests on masonry walls was considered. The experimentally obtained crack patterns were investigated and their correlation with damage propagation was studied. Using a software based on the Distinct Element Method, a numerical model was developed and validated against full-scale experimental tests obtained from the literature. Wall panels representing common typologies of house façades of unreinforced masonry buildings in Northern Europe i.e. near the Groningen gas field in the Netherlands, were numerically investigated. The accumulated damage within the seismic response of the masonry walls was investigated by means of representative harmonic load excitations and an incremental dynamic analysis based on induced seismicity records from Groningen region. The ability of this index to capture different damage situations is demonstrated. The proposed methodology could also be applied to quantify damage and accumulation in masonry during strong earthquakes and aftershocks too.

Seismic Slip-Pulse Experiments Simulate Induced Earthquake Rupture in the Groningen Gas Field.

Rock materials show dramatic dynamic weakening in large‐displacement (m), high‐velocity (∼1 m/s) friction experiments, providing a mechanism for the generation of large, natural earthquakes. However, whether such weakening occurs during induced M3‐4 earthquakes (dm displacements) is unknown. We performed rotary‐shear experiments on simulated fault gouges prepared from the source‐, reservoir‐ and caprock formations present in the seismogenic Groningen gas field (Netherlands). Water‐saturated gouges were subjected to a slip pulse reaching a peak circumferential velocity of 1.2–1.7 m/s and total displacements of 13–20 cm, at 2.5–20 MPa normal stress. The results show 22%–81% dynamic weakening within 5–12 cm of slip, depending on normal stress and gouge composition. At 20 MPa normal stress, dynamic weakening from peak friction coefficients of 0.4–0.9 to 0.19–0.27 was observed, probably through thermal pressurization. We infer that similar effects play a key role during induced seismic slip on faults in the Groningen and other reservoir systems.

Defining the usable bandwidth of weak-motion records: application to induced seismicity in the Groningen Gas Field, the Netherlands

Seismic hazard and risk analyses are increasingly tapping into the previously underused resource of local weak-motion records. This is facilitating the development of local- or even application-specific models for the characterisation of earthquake ground motion. In turn, this offers the opportunity to derive non- or partially non-ergodic models and significantly reduce bias and uncertainty. However, weak-motion data, while carrying important information about local earthquake source, path and site effects, are susceptible to noise. We show that high-frequency noise has a record-, or region-specific, impact on pseudo-spectral acceleration (PSA). This impact depends on the shape of the records’ Fourier amplitude spectrum (FAS): PSA from moderately to highly damped ‘soil’ records (e.g. Groningen, the Netherlands) is much less susceptible to high-frequency noise than PSA from weakly damped ‘rock’ records (e.g. Eastern North America). We make use of simulated ground motion records to develop a parametric model for the lower usable period of PSA (Tmin). The model accounts for the impact of high-frequency noise on PSA, conditional on easily measured parameters characterising the shape of a record’s FAS. We then present a workflow, describing processing undertaken for records of induced seismicity from the Groningen gas field. The workflow includes the definition of maximum and minimum usable frequencies and periods of FAS and PSA, respectively. As part of the workflow, we present an approach that considers multiple estimates of Tmin. These include the parametric model and, additionally, record-specific hybrid simulations that artificially extend or modify time series’ FAS beyond the noise floor to assess subsequent impacts on PSA.

Characterization of the near-surface shear wave attenuation in the Groningen gas field using borehole recording

SUMMARYThe Groningen gas field located in the northern Netherlands is the largest gas field in Europe. Many induced earthquakes have occurred since the gas extraction commenced in the 1960s, especially in the last two decades, which have caused concerning social and economic problems. To better quantify the seismic hazards caused by the induced earthquakes, the near surface attenuation in the Groningen area needs to be characterized from actual earthquakes. The spectral ratio method based on multiwindow analysis has proven to be a powerful tool to estimate the S-wave attenuation parameters of shallow sediments. In this study, waveforms from earthquake events recorded by the dense network of 70 boreholes, each of which is equipped with four geophones are used to estimate the S-wave attenuation parameters of shallow sediments. There are strong lateral variations in the near-surface attenuation, which is structurally similar to the S-wave model from a previous study in the same area. Besides, the calculated $Q_S^{ - 1}$ values are inversely proportional to the frequency and depth. The averaged ${\rm{\ }}Q_S^{ - 1}$ values for depths between 0–150 m and 0–200 m are close to those of the soft soil Sendai basin in Japan. The $Q_S^{ - 1}$ values from this study have been used to model a surface peak ground acceleration (PGA) map. As expected, both near surface $Q_S^{ - 1}$ values and source mechanism have a profound effect on the modelled PGA.

Regional Subsidence Analysis Through a Multi‐Scale Modeling Framework Based on Breakage Mechanics

AbstractAlthough poroelastic models are often used to explain the delay between subsurface fluid depletion and ground subsidence, inelastic compaction involving permanent changes of rock microstructure may exacerbate hydro‐mechanical coupling, thus influencing the interpretation of measurements and long‐term forecasts. Here, a multi‐scale modeling approach is discussed, which accounts for the inherent connection between rock microstructure, hydraulic conductivity, and pore compaction. A constitutive model built within the framework of breakage mechanics is proposed to link the hydraulic conductivity of granular rocks with inelastic deformations and changes in grading caused by injection‐depletion cycles at stress levels far from yielding. The proposed model has been incorporated into large‐scale simulation frameworks, thus enabling the spatiotemporal mapping of regional subsidence through a hybrid, semi‐analytical approach. Numerical results based on this strategy show that the model allows isolating near‐field and far‐field effects into the computation of land subsidence and can generate forecasts for different modeling scenarios (e.g., elastic and inelastic compaction, constant permeability, and concurrent change of compressibility and permeability). In particular, examples of simulations for the case of the Groningen gas field are discussed, showing the model capabilities to use both field measurements and laboratory tests for the generation of reasonable subsidence maps, without expensive computational costs. Results indicate that ignoring coupled inelastic effects has major consequences on the predicted timescale of subsidence. Specifically, while all the model scenarios produced similar long‐term ground settlements, those ignoring breakage‐dependent permeability changes result in a variation of the temporal window of residual subsidence of the order of several decades.

Interevent-time distribution and aftershock frequency in non-stationary induced seismicity

The initial footprint of an earthquake can be extended considerably by triggering of clustered aftershocks. Such earthquake–earthquake interactions have been studied extensively for data-rich, stationary natural seismicity. Induced seismicity, however, is intrinsically inhomogeneous in time and space and may have a limited catalog of events; this may hamper the distinction between human-induced background events and triggered aftershocks. Here we introduce a novel Gamma Accelerated-Failure-Time model for efficiently analyzing interevent-time distributions in such cases. It addresses the spatiotemporal variation and quantifies, per event, the probability of each event to have been triggered. Distentangling the obscuring aftershocks from the background events is a crucial step to better understand the causal relationship between operational parameters and non-stationary induced seismicity. Applied to the Groningen gas field in the North of the Netherlands, our model elucidates geological and operational drivers of seismicity and has been used to test for aftershock triggering. We find that the hazard rate in Groningen is indeed enhanced after each event and conclude that aftershock triggering cannot be ignored. In particular we find that the non-stationary interevent-time distribution is well described by our Gamma model. This model suggests that 27.0(± 8.5)% of the recorded events in the Groningen field can be attributed to triggering.

Using machine learning for model benchmarking and forecasting of depletion-induced seismicity in the Groningen gas field

The Groningen gas field in the Netherlands is experiencing induced seismicity as a result of ongoing depletion. The physical mechanisms that control seismicity have been studied through rock mechanical experiments and combined physical-statistical models to support development of a framework to forecast induced-seismicity risks. To investigate whether machine learning techniques such as Random Forests and Support Vector Machines bring new insights into forecasts of induced seismicity rates in space and time, a pipeline is designed that extends time-series analysis methods to a spatiotemporal framework with a factorial setup, which allows probing a large parameter space of plausible modelling assumptions, followed by a statistical meta-analysis to account for the intrinsic uncertainties in subsurface data and to ensure statistical significance and robustness of results. The pipeline includes model validation using e.g. likelihood ratio tests against average depletion thickness and strain thickness baselines to establish whether the models have statistically significant forecasting power. The methodology is applied to forecast seismicity for two distinctly different gas production scenarios. Results show that seismicity forecasts generated using Support Vector Machines significantly outperform beforementioned baselines. Forecasts from the method hint at decreasing seismicity rates within the next 5 years, in a conservative production scenario, and no such decrease in a higher depletion scenario, although due to the small effective sample size no statistically solid statement of this kind can be made. The presented approach can be used to make forecasts beyond the investigated 5-years period, although this requires addition of limited physics-based constraints to avoid unphysical forecasts.

Variations between foundation-level recordings and free-field earthquake ground motions: numerical study at soft-soil sites

It is typically assumed that accelerograms obtained from ground motion recording instruments adequately installed inside a building can be considered as representative of free-field conditions if the structure is less than 2-3 storeys high and does not feature a basement or a massive foundation system. In this work, such an assumption is verified through nonlinear soil-block finite element analyses, considering the induced seismicity Groningen gas field as a case-study and both types of structures. The capability of the numerical model in adequately reproducing nonlinear soil response and capturing interaction between soil and structure is first verified through a number of analyses and cross-checks. Considering soil profiles at the locations of different recording stations in the Groningen field, together with accelerograms from recent events in the region, the results obtained with the numerical model confirm that recordings from instruments located at the base of lightweight structures are not affected by Soil-Structure Interaction (SSI) effects, neither kinematic nor inertial. Heavy structures with a basement, however, do exhibit recordings with variations with respect to the free-field; these are identified and possible corrections to the recordings are analysed. Another set of analyses then demonstrates that soil consolidation prior to the construction of lightweight structures resting on shallow foundations may instead lead to foundation-level recordings that feature high-frequency ground motion amplitudes slightly lower than their free-field counterparts. Additional simplified analyses show the dependency of the recordings on inertial SSI effects due to deformability of the ground supporting the structure.

Using horizontal-to-vertical spectral ratios to construct shear-wave velocity profiles

Abstract. For seismic hazard assessment and earthquake hypocentre localization, detailed shear-wave velocity profiles are an important input parameter. Here, we present a method to construct a shear-wave velocity profiles for a deep unconsolidated sedimentary layer by using strong teleseismic phases and the ambient noise field. Gas extraction in the Groningen field, in the northern part of the Netherlands, is causing low-magnitude, induced seismic events. This region forms an excellent case study due to the presence of a permanent borehole network and detailed subsurface knowledge. Instead of conventional horizontal-to-vertical spectral ratios (H∕V ratios) from amplitude spectra, we calculate power spectral densities and use those as input for H∕V calculations. The strong teleseisms provide resonance recordings at low frequencies, where the seismic noise field is too weak to be recorded well with the employed geophones and accelerometers. The H∕V ratios of the ambient noise field are compared with several forward modelling approaches to quality check the teleseism-based shear-wave velocity profiles. Using the well-constrained depth of the sedimentary basin, we invert the H∕V ratios for velocity profiles. A close relationship is observed between the H∕V spectral ratios from the ambient noise field, shear-wave resonance frequencies and Rayleigh-wave ellipticity. By processing only five teleseismic events, we are able to derive shear-wave velocities for the deeper sedimentary sequence with a 7 % bias in comparison with the existing detailed velocity model for the Cenozoic sediments overlying the Groningen gas field. Furthermore, a relation between resonance frequency and unconsolidated sediment thickness is derived, to be used in other areas in the Netherlands, where detailed depth maps are not available.

Stress‐Dependent Magnitudes of Induced Earthquakes in the Groningen Gas Field

AbstractGeological faults may produce earthquakes under induced stresses associated with hydrocarbon extraction, geothermal extraction, or CO2 storage. The associated risks depend on the frequency and magnitude of these earthquakes. Within seismic risk analysis, the exceedance probability of seismic moments, , is treated as a pure power law distribution, , where the power law exponent, β, may vary in time or space or with stress. Insights from statistical mechanics theories of brittle failure, statistical seismology, and acoustic emissions experiments all indicate this pure power law may contain an exponential taper, , where the taper strength, ζ, decreases with increasing stress. The role of this taper is to significantly reduce the probability of earthquakes larger than ζ −1 relative to the pure power law. We review the existing theoretical and observational evidence for a stress‐dependent exponential taper to motivate a range of magnitude models suitable for induced seismicity risk analysis. These include stress‐invariant models with and without a taper, stress‐dependent β models without a taper, and stress‐dependent ζ models. For each of these models, we evaluated their forecast performance within the Groningen gas field in the Netherlands using a combination of Bayesian inference, and simulations. Our results show that the stress‐dependent ζ model with constant β likely offers (75–85%) higher performance forecasts than the stress‐dependent β models with . This model also lowers the magnitudes with a 10% and 1% chance of exceedance over the next 5 years of gas production from 4.3 to 3.7 and from 5.5 to 4.3, respectively.

Earthquake ground motion modeling of induced seismicity in the Groningen gas field

AbstractA physics‐based numerical approach is used to characterize earthquake ground motion due to induced seismicity in the Groningen gas field and to improve empirical ground motion models for seismic hazard and risk assessment. To this end, a large‐scale (20 km × 20 km) heterogeneous 3D seismic wave propagation model for the Groningen area is constructed, based on the significant bulk of available geological, geophysical, geotechnical, and seismological data. Results of physics‐based numerical simulations are validated against the ground motion recordings of the January 8, 2018, ML 3.4 Zeerijp earthquake. Taking advantage of suitable models of slip time functions at the seismic source and of the detailed geophysical model, the numerical simulations are found to reproduce accurately the observed features of ground motions at epicentral distances less than 10 km, in a broad frequency range, up to about 8 Hz. A sensitivity analysis is also addressed to discuss the impact of 3D underground geological features, the stochastic variability of seismic velocities and the frequency dependence of the quality factor. Amongst others, results point out some key features related to 3D seismic wave propagation, such as the magnitude and distance dependence of site amplification functions, that may be relevant to the improvement of the empirical models for earthquake ground motion prediction.

Microphysics of Inelastic Deformation in Reservoir Sandstones from the Seismogenic Center of the Groningen Gas Field

Physics-based assessment of the effects of hydrocarbon production from sandstone reservoirs on induced subsidence and seismicity hinges on understanding the processes governing compaction of the reservoir. Compaction strains are typically small (ε < 1%) and may be elastic (recoverable), or partly inelastic (permanent), as implied by recent experiments. To describe the inelastic contribution in the seismogenic Groningen gas field, a Cam–clay-type plasticity model was recently developed, based on the triaxial test data obtained for sandstones from the Groningen reservoir (strain rate ~ 10−5 s−1). To underpin the applicability of this model at production-driven strain rates (10−12 s−1), we develop a simplified microphysical model, based on the deformation mechanisms observed in triaxial experiments at in situ conditions and compaction strains (ε < 1%). These mechanisms include consolidation of and slip on µm-thick clay films within sandstone grain contacts, plus intragranular cracking. The mechanical behavior implied by this model agrees favourably with the experimental data and Cam–clay description of the sandstone behavior. At reservoir-relevant strains, the observed behavior is largely accounted for by consolidation of and slip on the intergranular clay films. A simple analysis shows that such clay film deformation is virtually time insensitive at current stresses in the Groningen reservoir, so that reservoir compaction by these mechanisms is also expected to be time insensitive. The Cam–clay model is accordingly anticipated to describe the main trends in compaction behavior at the decade time scales relevant to the field, although compaction strains and lateral stresses may be slightly underestimated due to other, smaller creep effects seen in experiments.

Probabilistic Moment Tensor Inversion for Hydrocarbon-Induced Seismicity in the Groningen Gas Field, The Netherlands, Part 1: Testing

ABSTRACTSince 1991, induced earthquakes have been observed and linked to gas production in the Groningen field. Recorded waveforms are complex, resulting partly from a Zechstein salt layer overlying the reservoir and partly from free-surface reverberations, internal multiples, interface conversions, guided waves, and waves diving below the reservoir. Therefore, picking of polarities or amplitudes for use in moment tensor inversion is problematic, whereas phase identification may be circumvented employing full waveform techniques. Although moment tensors have become a basic tool to analyze earthquake sources, their uncertainties are rarely reported. We introduce a method for probabilistic moment tensor estimation and demonstrate its use on the basis of a single event within the Groningen field, concentrating on detailed tests of input data and inversion parameters to derive rules of good practice for moment tensor estimation of events recorded in the Groningen field. In addition to the moment tensor, event locations are provided. Hypocenters estimated simultaneously with moment tensors are often less sensitive to uncertainties in crustal structure, which is pertinent for the application to the Groningen field, because the task of relating earthquakes to specific faults hitherto suffers from a limited resolution of earthquake locations. Because of the probabilistic approach, parameter trade-offs, uncertainties, and ambiguities are mapped. In addition, the implemented bootstrap method implicitly accounts for modeling errors affecting every station and phase differently. A local 1D velocity model extracted from a full 3D velocity model yields more consistent results than other models applied previously. For all velocity models and combinations of input data tested, a shift in location of 1 km to the south is observed for the test event compared to the public catalog. A full moment tensor computed employing the local 1D velocity model features negative isotropic components and may be interpreted as normal fault and collapse at reservoir level.

Probabilistic Moment Tensor Inversion for Hydrocarbon-Induced Seismicity in the Groningen Gas Field, the Netherlands, Part 2: Application

ABSTRACTRecent developments in the densification of the seismic network covering the Groningen gas field allow a more detailed study of the connection between induced seismicity and reactivated faults around the gas reservoir at 3 km depth. With the reduction of the average station distance from 20 km to 4–5 km, a probabilistic full-waveform moment tensor inversion procedure could be applied, resulting in both improved hypocenter location accuracy and full moment tensor solutions for events of M≥2.0 recorded in the period 2016–2019. Hypocenter locations as output from the moment tensor inversion are compared to locations from the application of other methods and are found similar within 250 m distance. Moment tensor results show that the double-couple (DC) solutions are in accordance with the known structure, namely normal faulting along 50°–70° dipping faults. Comparison with reprocessed 3D seismic sections, extended to a depth of 6–7 km, demonstrate that (a) most events occur along faults with a small throw and (b) reactivated faults in the reservoir often continue downward in the Carboniferous underburden. From non-DC contributions, the isotropic (ISO) component is dominant and shows consistent negative values, which is expected in a compacting medium. There is some indication that events connected to faults with a large throw (&amp;gt;70 m) exhibit the largest ISO component (40%–50%).

Compaction of the Groningen gas reservoir investigated with train noise

SummaryInduced seismicity in the Groningen gas field in the Netherlands has been related to reservoir compaction caused by gas pressure depletion. In situ measurement of compaction is therefore relevant for seismic hazard assessment. In this study, we investigated the potential of passively recorded deep borehole noise data to detect temporal variations in the Groningen reservoir. Train signals recorded by an array of 10 geophones at reservoir depth were selected from the continuous noise data for two 5-month deployments in 2015. Interferometry by deconvolution was applied to the high-frequency train signals that acted as stable, repetitive noise sources. Direct intergeophone P and S wave traveltimes were then used to construct the P- and S-wave velocity structure along the geophone array. The resulting models agree with independently obtained velocity profiles and have very small errors. Most intergeophone P wave traveltimes showed decreasing traveltimes per deployment period, suggestive of compaction. However, the retrieved traveltime changes are very small, up to tens of microseconds per deployment period, with uncertainties that are of similar size, about 10 microseconds. An unambiguous interpretation in terms of compaction is therefore not warranted, although the 10 μs error per 5-month period is probably smaller than can be achieved from active time-lapse seismic surveys that are commonly used to measure reservoir compaction. The direct P-wave amplitudes of the train-signal deconvolutions were investigated for additional imprints of compaction. Whereas the P-wave amplitudes consistently increased during the second deployment, suggestive of compaction, no such trend was observed for the first deployment, rendering the interpretation of compaction inconclusive. Our results therefore present hints, but no obvious effects of compaction in the Groningen reservoir. Yet, this study demonstrates that the approach of deconvolution interferometry applied to deep borehole data allows monitoring of small temporal changes in the subsurface for stable repetitive noise sources such as trains.

Ground-Motion Attenuation, Stress Drop, and Directivity of Induced Events in the Groningen Gas Field by Spectral Inversion of Borehole Records

ABSTRACTProduction-induced earthquakes in the Groningen gas field caused damage to buildings and concerns for the population, the gas-field owner, and the local and national authorities and institutions. The largest event (ML=3.6) occurred in 2012 near Huizinge, and, despite the subsequent decision of the Dutch government to reduce the gas production in the following years, similar magnitude events occurred in 2018 and 2019 (ML=3.4). Thanks to the improvement of the local seismic networks in the last years, recent events provide a large number of recordings and an unprecedented opportunity to study the characteristics of induced earthquakes in the Groningen gas field and related ground motions.In this study, we exploit the S-wave Fourier amplitude spectra recorded by the 200 m depth borehole sensors of the G network from 2015 to 2019 to derive source and attenuation parameters for ML≥2 induced earthquakes. The borehole spectra are decomposed into source, attenuation, and site nonparametric functions, and parametric models are then adopted to determine moment magnitudes, corner frequencies, and stress drops of 21 events. Attenuation and source parameters are discussed and compared with previous estimates for the region. The impact of destructive interference of upgoing and downgoing waves at borehole depth on the derived parameters is also discussed and assessed to be minor.The analysis of the apparent source spectra reveals that several events show rupture directivity and provides clear observations of frequency-dependent directivity effects in induced earthquakes. The estimated rupture direction shows a good agreement with orientation of pre-existing faults within the reservoir. Our results confirm that rupture directivity is still an important factor for small-magnitude induced events, affecting the amplitude of recorded short-period response spectra and causing relevant spatial ground-motion variability.