# Modeling of induced seismicity and ground vibrations associated with geologic CO2 storage, and assessing their effects on surface structures and human perception,

J. Rutqvist, F. Cappa, A. Rinaldi, M. Godano

► To cite this version:

J. Rutqvist, F. Cappa, A. Rinaldi, M. Godano. Modeling of induced seismicity and ground vi- brations associated with geologic CO2 storage, and assessing their effects on surface structures and human perception,. International Journal of Greenhouse Gas Control, 2014, 24, pp.64-77. 10.1016/j.ijggc.2014.02.017 . hal-00996718

HAL Id: hal-00996718

https://hal.science/hal-00996718

Submitted on 16 Sep 2021

**HAL** is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>
<figure><img src="image_2.png" /></figure>

# Modeling of induced seismicity and ground vibrations associated with geologic CO₂ storage, and assessing their effects on surface structures and human perception

<figure><img src="image_3.png" /></figure>

CrossMark

Jonny Rutqvist<sup>a, *</sup>, Frédéric Cappa<sup>a, b</sup>, Antonio P. Rinaldi<sup>a</sup>, Maxime Godano<sup>b, 1</sup>

## ARTICLE INFO

**Article history:**

Received 19 July 2013

Received in revised form 21 January 2014

Accepted 28 February 2014

Available online 22 March 2014

**Keywords:**

CO₂ injection

Fault rupture

Induced seismicity

Ground acceleration

## ABSTRACT

In this paper, we present model simulations of ground motions caused by CO₂-injection-induced fault reactivation and analyze the results in terms of the potential for damage to ground surface structures and nuisance to the local human population. It is an integrated analysis from cause to consequence, including the whole chain of processes starting from earthquake inception in the subsurface, wave propagation toward the ground surface, and assessment of the consequences of ground vibration. For a small magnitude (*M*<sub>w</sub> = 3) event at a hypocenter depth of about 1000 m, we first used the simulated ground-motion wave train in an inverse analysis to estimate source parameters (moment magnitude, rupture dimensions and stress drop), achieving good agreement and thereby verifying the modeling of the chain of processes from earthquake inception to ground vibration. We then analyzed the ground vibration results in terms of peak ground acceleration (PGA), peak ground velocity (PGV) and frequency content, with comparison to U.S. Geological Survey's instrumental intensity scales for earthquakes and the U.S. Bureau of Mines' vibration criteria for cosmetic damage to buildings, as well as human-perception vibration limits. Our results confirm the appropriateness of using PGV (rather than PGA) and frequency for the evaluation of potential ground-vibration effects on structures and humans from shallow injection-induced seismic events. For the considered synthetic *M*<sub>w</sub> = 3 event, our analysis showed that the short duration, high frequency ground motion may not cause any significant damage to surface structures, but would certainly be felt by the local population.

Published by Elsevier Ltd. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/3.0/).

## 1. Introduction

Deep underground injection of carbon dioxide (CO₂), or geo-logic carbon storage (GCS) in deep sedimentary formation such as saline aquifers and depleted hydrocarbon reservoirs is one impor-tant option being considered for sequestering CO₂ (and thereby reducing the emission of greenhouse gases to the atmosphere). Coupled fluid flow and geomechanical processes can be significant associated with such injection (Orlic, 2009; Rutqvist, 2012). In par-ticular the potential for fault reactivation and induced seismicity have recently received increased attention from CO₂ sequestration stakeholders and media, especially after two recent high-profile publications. First, in anticipation of public concern about the

potential for energy development projects to induce seismicity, the U.S. Congress directed the U.S. Department of Energy to request that the National Research Council examine the scale, scope, and consequences of seismicity induced during fluid injection and with-drawal activities related to geothermal energy development, oil and gas development (including shale gas recovery), and carbon capture and storage (CCS). That study, published in 2012 (National Research Council, 2012), concluded that projects that involve injecting or extracting large net volumes of fluids over long periods of time (such as CCS) may have the potential to induce large seismic events, though insufficient information exists for understanding this potential, because no large-scale CCS projects are yet in opera-tion. In another article, published in the Proceedings of the National Academy of Sciences of the United States of America (PNAS), Zoback and Gorelick (2012) concluded that there is a high probability that earthquakes would be triggered by injection of large volumes of CO₂ into the brittle rocks commonly found in continental interiors, and because even small- to moderate-sized earthquakes threaten the seal integrity of CO₂ repositories, large-scale CCS would be

<footnote><sup>a</sup> Earth Sciences Division, Lawrence Berkeley National Laboratory, University of California, Berkeley, CA 94720, USA</footnote>

<footnote><sup>b</sup> Geoazur, University of Nice-Sophia Antipolis, Côte d'Azur Observatory, Sophia-Antipolis, France</footnote>

<footnote>* Corresponding author at: One Cyclotron Road, MS 90-1116, USA. Tel.: +1 510 486 5432; fax: +1 510 486 5686.</footnote>

<footnote>E-mail address: jrutqvist@lbl.gov (J. Rutqvist).</footnote>

<footnote n="1">¹ Now at Institut de Physique du Globe de Paris, Paris, France.</footnote>
a risky and likely unsuccessful strategy for significantly reducing greenhouse gas emissions to the atmosphere. Zoback and Gorelick (2012) were concerned about basin-scale fluid-pressure disturbances that would be expected from multiple large-scale injection operations associated with coal power plants in areas such as the midwestern U.S.—even though they still maintained that geologic CO₂ sequestration could be valuable and useful at a local level, in specific situations, such as in the case of the injection at the Utsira Formation, at the Sleipner Gas Field in the North Sea.

In their article, Zoback and Gorelick (2012) mentioned reservoir-induced earthquakes associated with dam construction and water-reservoir impoundment as a good analog for the seismicity potential that could be induced by large-scale CO₂ injection, because both activities cause pressure changes that act over large areas and are persistent for long periods. They also listed a number of recent seismic small-to-moderate events in the U.S. that seem to have been triggered by deep waste-water injection, referring to the critically stressed nature of the Earth's crust, which suggests that most everywhere in the continental interior, a subset of preexisting faults in the crust is potentially active in the current stress field (Zoback and Zoback, 1980). Because of the critically stressed nature of the crust, fluid injection in deep wells can trigger earthquakes when the injection increases the pore pressure in the vicinity of potentially active faults. In this regard, it is worth mentioning recent seismic events in Prague, Oklahoma in 2010 and 2011, including a damage-causing M<sub>w</sub> 5.7 earthquake, that according to some scientists were associated with deep wastewater injection (Keranen et al., 2013). However, it should also be mentioning that the subsurface at many locations is not necessarily critically stressed, especially not in sedimentary basins at shallow depths, (800 m to perhaps 2–3 km) which are of interest to CO₂ storage. Thus, the potential for induced seismicity will be dependent on the site specific conditions, including local stress field and rock properties.

With regard to shallow injection-induced seismicity of relatively small magnitudes (e.g. M<sub>w</sub> <3), it has been recognized that the frequency content might be too high to cause any structural damage, but they could still be felt and considered alarming by humans (Majer et al., 2012). There are no accepted specific guidelines on the control or definition of acceptable levels of induced seismic ground-shaking. In relation to enhanced geothermal system (EGS) activities, Majer et al. (2012), in developing a best-practice protocol for induced seismicity, suggested considering criteria developed for blasting and construction activities. The motivation is that seismic energy, frequency bandwidth, and range of events generated from these activities is similar to that from induced seismic events. Such an approach was also adopted by Bommer et al. (2006) and Häring et al. (2008), related to the development of EGS sites involving hydraulic stimulation of fractured rock masses. In these cases, a traffic-light system was used, including thresholds for human disturbance and for damage on vulnerable houses, defined in terms of peak ground velocity (PGV) or peak ground acceleration (PGA), with thresholds inferred from recommendations regarding tolerable vibration levels for blasting and pile driving. Similar guidelines as those developed and applied for EGS may also be applied for risk managements related to new GCS projects as proposed in a recent International Energy Agency Report (IEAGHG, 2013).

In this paper, we present model simulations of ground motions caused by CO₂-injection-induced fault reactivation and analyze the results in terms of the potential for damage to ground-surface structures and nuisance to the local population (Fig. 1). Starting from a previously developed numerical model by Cappa and Rutqvist (2012) we present simulation results of an integrated analysis from cause to consequence, including the whole chain of processes, from earthquake inception in the subsurface, wave propagation toward

<figure><img src="image_4.png" /><figcaption>**Fig. 1.** Schematic of injection-induced fault reactivation, wave propagation, and ground motions, and potential impact on surface structures and human perception analyzed in this study.</figcaption></figure>

the ground surface, and to assess the consequences of ground vibration (Fig. 1). In this paper, the main focus is on the analysis of ground vibrations in terms of the potential for inducing damage to building structures and human perception. To verify that our ground vibration results make sense related to the simulated earthquake source (fault reactivation), we analyzed the ground-motion wave train to independently invert for source parameters (moment magnitude, rupture dimensions, and stress drop). A good agreement between the synthetic and inverted source parameters verifies that we correctly model the chain of processes from earthquake inception to ground vibration. We present results of ground acceleration and velocity, and show how these are distributed over multiple locations on the ground surface. Finally, we compare calculated values with intensity scales and vibration criteria from the construction and mining industry, showing that for such a small and shallow fault reactivation event the short duration, high frequency ground motion may not cause any significant damage to surface structures, but would certainly be felt by the local population.

## 2. Ground vibrations and potential effects on structures and humans

In this study, we consider research and guidelines developed in seismology and earthquake engineering, as well as standards and guidelines related to ground vibrations associated with construction and mining activities, including blasting. In particular, blasting has been suggested as a good analog to injection-induces seismic events, because of their similarities in terms of duration and frequency content (Bommer et al., 2006; Dowding and Meissner, 2011; Majer et al., 2012).

In the context of seismology and earthquake engineering, ground acceleration and velocity are important input parameters
<figure><img src="image_5.png" /></figure>

<table><caption>**Fig. 2.** USGS ground-surface shake map around the November 6, 2011, $M_w$ = 5.7, Prague earthquake and a table with the color scale related to PGA, PGV, and Modified Mercalli Intensity. (http://earthquake.usgs.gov/earthquakes/shakemap/global/shake/b0006klz/).</caption><thead><tr><th>PERCEIVED SHAKING</th><th>Not felt</th><th>Weak</th><th>Light</th><th>Moderate</th><th>Strong</th><th>Very strong</th><th>Severe</th><th>Violent</th><th>Extreme</th></tr></thead><tbody><tr><th>POTENTIAL DAMAGE</th><td>none</td><td>none</td><td>none</td><td>Very light</td><td>Light</td><td>Moderate</td><td>Moderate/Heavy</td><td>Heavy</td><td>Very Heavy</td></tr><tr><th>PEAK ACC.(kg/s)</th><td>&lt;1.7</td><td>1.7-1.4</td><td>1.4-3.9</td><td>3.9-9.2</td><td>9.2-18</td><td>18-34</td><td>34-65</td><td>65-124</td><td>&gt;124</td></tr><tr><th>PEAK VEL.(cm/s)</th><td>&lt;0.1</td><td>0.1-1.1</td><td>1.1-3.4</td><td>3.4-8.1</td><td>8.1-16</td><td>16-31</td><td>31-69</td><td>60-116</td><td>&gt;116</td></tr><tr><th>INSTRUMENTAL INTENSITY</th><td>I</td><td>II-III</td><td>IV</td><td>V</td><td>VI</td><td>VII</td><td>VIII</td><td>IX</td><td>X</td></tr></tbody></table>

for earthquake engineering in developing response spectra used in the seismic design of buildings (Chopra, 2012). The acceleration and velocity from a wave train may be recorded by seismic instruments and the maximum values, i.e., PGA and PGV, provide measurements of instrumental intensity. Other intensity scales measure felt intensity, e.g., the Mercalli scale, and are based on eyewitness reports, felt shaking, and observed damage. The United States Geological Survey (USGS) has developed an instrumental intensity scale which maps PGA and PGV on a scale similar to the felt Mercalli scale (Wald et al., 1999). These values are used to provide near-real-time maps of ground motion and shaking intensity following significant earthquakes. For example, Fig. 2 shows such a map developed from the aforementioned November 6, 2011, $M_w$ = 5.7 Prague earthquake. Fig. 2, bottom shows the color code for the Modified Mercalli Intensity ($I_{mm}$) values; these have also been empirically correlated with PGA and PGV. However, for injection-induced events, the USGS scale may be very approximate, since it was developed based upon eight larger California earthquakes of $M \ge 5.8$ (Wald et al., 1999; Worden et al., 2012) that were tectonic events, which occur much deeper than shallow injection-induced seismic events.

In the case of shallow, injection-induced seismicity, the frequencies generated are generally too high to cause significant structural damage. As pointed out by Van Eck et al. (2006), large PGAs caused by small, shallow earthquakes are not unusual (Fletcher et al., 1983; McGarr, 1984; Wu et al., 2003; Sharma et al., 2013). However, in those cases, the strong accelerations are usually short, often only one cycle, and the PGV is in these cases often a more robust damage indicator (Van Eck et al., 2006). Both Van Eck et al. (2006) and Bommer et al. (2006) studied the response spectra from recordings of injection/production-induced and natural events of low to moderate magnitudes showing a pronounced peak at a frequency of

about 10 Hz. A frequency of 10 Hz is near the upper limit of the natural frequency range of a one-story building, and near the lower limit of the natural frequencies for walls and floor structures (Dowding, 1996; Dowding and Meissner, 2011).

Following suggestions by Majer et al. (2012), in their best-practice protocol for induced seismicity associated with EGS activities, here we consider criteria developed for blasting and construction activities. The seismic energy, frequency bandwidth, and range of events generated from these activities are similar to that from induced seismic events. Thresholds for human disturbance and for damage on vulnerable houses are defined in terms of PGV, with thresholds inferred from recommendations for tolerable vibration levels due to blasting and pile driving, and from correlations between PGV and seismic intensity. These PGV criteria are almost universally used by the construction and mining industry to assess the potential for threshold cracking due to blasting, and are also employed in many commercially available vibration-monitoring systems (Svinkin, 2004).

For building damage criteria associated with blasting and mining activities, the “Z-curve” or “Siskind curve” is the information most often cited. It was published by Siskind et al. (1980) as a result of an extensive study conducted in the late 1970s using data from numerous blasts, gathered from monitoring for different types of structures (Fig. 3a). The resulting graph shows the limits in peak particle velocity (PPV) in inch per second (ips), recommended by the U.S. Bureau of Mines (USBM) to preclude cosmetic damage to plaster and drywall, i.e., the most fragile building materials. As the frequency of the ground motion changes from low to high (1–40 Hz) the structure responds less and the limits increase.

Fig. 3b shows the USBM PGV limits for cosmetic damage (and modifications to these limits by the Office of Surface Mining–OSM), together with curves for human perception. The limits for human perceptions are those specified, as proposed guidelines in the U.S. Army Engineering Manual EM 1110-2-3800 (USACE, 1972) and also considered for ground-vibration-induced pile driving and associated with EGS developments (Athanasopoulos and Pelekis, 2000; Bommer et al., 2006). Fig. 3b shows that humans feel high frequency events, e.g., high acceleration. However, these human-perception curves, while applied for recommendations related to blasting, actually originate from studies of human exposure to steady-state sinusoidal vibrations, i.e., not short-duration events. For occasional injection-induced events, it is reasonable to use these criteria for a general description of human perception, but not as an absolute vibration limit. Thus, in this paper, we will use these human perception criteria merely for discussion in terms of potential nuisance to the local population.

### 3. Simulation of injection-induced fault reactivation and wave propagation

The simulation presented in this paper was conducted with the coupled thermo-hydro-mechanical code TOUGH-FLAC (Rutqvist et al., 2002; Rutqvist, 2011), which has previously been applied to study fluid-induced fault reactivation, in both generic and site specific studies (Rutqvist et al., 2007; Cappa et al., 2009; Cappa and Rutqvist, 2011a,b, 2012; Mazzoldi et al., 2012; Rinaldi and Rutqvist, 2013). TOUGH2 is a finite-volume code for the simulation of multiphase fluid flow (Pruess et al., 2011); FLAC³ is a finite-difference code for the simulation of geomechanics (Itasca, 2011), here applied for modeling of fault rupture as well as wave propagation. The approach adopted here involves coupled fluid flow and geomechanical numerical modeling to simulate CO₂ injection and fault rupture. Then, seismological theories are used to calculate the seismic source parameters. A strain-softening fault constitutive model enabled us to model sudden, dynamic fault rupture, as
<figure><img src="image_6.png" /><figcaption>**Fig. 3.** Ground-vibration damage and human-perception limits from mining and construction discussed in this study of CO₂-injection induced-seismicity: (a) ground vibration damage summary including USBM RI 8507 (Siskind et al., 1980) recommended safe limits, structural damage data means (lines), and damage observations (symbols) (modified from Svinkin, 2004); (b) USBM-recommended limits for cosmetic damage in plaster stucco and drywall and human-perception USAGE (1972) limits for blast vibration (after Bommer et al., 2006).</figcaption></figure>

well as provide a source for wave-propagation and ground-motion calculations (Cappa and Rutqvist, 2012). The wave-propagation was calculated using the FLAC³D dynamic analysis option based on an explicit finite difference scheme to solve the full equations of motion in a fully nonlinear analysis (see Itasca, 2011 for details).

The basic model geometry and input is the same as that used in Cappa and Rutqvist (2012) and involves a reservoir bounding fault with a dip angle of 80°, width of 2.5 m, and tectonic shear displacement offset (throw) of 125 m (Fig. 4). The model simulations

<figure><img src="image_7.png" /><figcaption>**Fig. 4.** Model geometry and boundary conditions and ground-motion recording stations (blue and red triangles) after Cappa and Rutqvist (2012). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of the article.)</figcaption></figure>

were conducted in a vertical cross section normal to the strike of the fault represented in a 2D plane-strain model (2 km × 2 km). The storage formation is 100 m thick and bounded at the top and bottom by low-permeability 150 m thick formations, which, in turn, is surrounded by two other permeable formations. Dynamic analysis required a very fine mesh with element sizes of 0.25 m along the rupture zone, in order to resolve the weakening process over the nucleation zone and to avoid mesh-induced artificial reflections (Cappa and Rutqvist, 2012).

Initial conditions for the model simulations were derived assuming (1) hydrostatic fluid pressure with the ground water table adjacent to the ground surface, (2) a temperature assuming a depth gradient of 25 °C/km with a ground-surface temperature of 10 °C, and (3) a vertical stress from the weight of the overburden rock for a bulk density of density $ρ=2260$ kg/m³, with the initial minimum horizontal stress being scaled by a factor 0.7 of the vertical stress, i.e. $σ_h = 0.7σ_v$. The minimum horizontal stress is assumed to be oriented normal to the strike of the fault plane, whereas the maximum horizontal stress is assumed to be oriented parallel to the strike of the fault plane and with a magnitude equal to the vertical stress, i.e., $σ_H = σ_v$. With these parameters, the initial fluid pressure and temperature at the depth of the CO₂ injection zone (1000 m) is about 10 MPa and 35 °C, respectively, whereas the vertical stress is 22.2 MPa and minimum horizontal stress is 15.5 MPa. The studied case in which a normal fault bounds a storage reservoir between two caprocks corresponds to a critical geometrical case for fault activation during CO₂ injection (Hawkes et al., 2005).

In the model simulations, CO₂ was injected as a point source within the storage formation, with a constant rate of 0.02 kg/(m s) (i.e. 630.72 tons/m/year). This is the injection rate for the half-symmetric model and per meter normal to the 2D model, and results in a reservoir pressure increase that would lead to fault reactivation within a few months of injection. Our simulation was conducted to intentionally induce fault reactivation, which occurs at a high reservoir pressure in an unfavorable stress regime. In fact, to reactivate the fault, we had to increase the reservoir pressure to exceed the initial minimum principal stress. This is an extreme case that in practice could be avoided by careful site characterizations, monitoring, and injection control. However, we
<table><caption>**Table 1**
Material parameters.</caption><thead><tr><th>Parameters</th><th>Storage formation</th><th>Caprock</th><th>Others formations</th><th>Fault</th></tr></thead><tbody><tr><td>Young's modulus (GPa)</td><td>10</td><td>10</td><td>10</td><td>5</td></tr><tr><td>Poisson' ratio</td><td>0.25</td><td>0.25</td><td>0.25</td><td>0.25</td></tr><tr><td>Porosity</td><td>0.1</td><td>0.01</td><td>0.1</td><td>0.1</td></tr><tr><td>Permeability (m<sup>2</sup>)</td><td>10<sup>-13</sup></td><td>10<sup>-19</sup></td><td>10<sup>-14</sup></td><td>10<sup>-16</sup></td></tr></tbody></table>

may compare this injection rate to a field situation. For example, assuming that injection would take place in a 1 km long horizontal well, the total injection rate for that well would be 630.72 × 2 × 1000 ≈ 1.26 × 10⁶ tons/year, i.e., about 1.3 million tons per year. We may also compare this to real injection rates at the In Salah CO₂ storage project, where injection rates have been 0.5–1.0 million tons per year distributed over 3 horizontal injection wells, each 1–1.5 km long (Rutqvist et al., 2010).

Properties for the permeable formations and the bounding low-permeable layers represent sandstone and shale, respectively (Table 1). The properties for the formation layers are not taken from a specific site but are realistic for sandstone and shale units (e.g. Zoback, 2007) and appropriate for a CO₂ sequestration site. That is, injection zone consists of a sufficiently porous and permeable sandstone and the caprock of a sufficiently tight shale to be an effective seal. The exact values of the elastic properties are not critical for the analysis of the potential for fault reactivation, but affect the wave propagation speed. In this case, based on the elastic constants in Table 1, the P-wave velocity, $C_p = 2360$ m/s, and the S-wave velocity, $C_s = 1402$ m/s. The fault properties, on the other hand, are very important for the potential for fault reactivation and induced seismicity, especially the frictional strength. In this study, the fault is modeled using an elasto-plastic anisotropic Mohr–Coulomb model, a built-in FLAC³D constitutive model (Cappa and Rutqvist, 2011a,b). The strain-softening fault geomechanical behavior was achieved by a reduction of the coefficient of friction during co-seismic shear slip from a peak static value, $\mu_s = 0.6$, to a dynamic value, $\mu_d = 0.2$, over a critical plastic strain of $10^{-3}$ (Cappa and Rutqvist, 2012). A coefficient of friction of 0.6 is commonly observed for laboratory samples as a lower limit value for the most common rocks and has also been inferred from field observations as the Earth's shallow crust being critically stressed for frictional failure (Byerlee, 1978; Zoback, 2007). However, clay-rich fault rock could have a much lower coefficient of friction, especially under wet conditions (Zoback, 2007; Samuelson and Spiers, 2012). A residual shear strength of 0.2 is not unusual for clay-rich fault gauge (Ikari et al., 2009), but we should relate our residual coefficient of friction to the dynamic friction coefficient in seismology, which depends on displacement rate, among other factors (Samuelson and Spiers, 2012). Thus, the values taken for the input parameters of the fault are uncertain, and we consider the selected parameters in this analysis as one set of possible parameters designed for generating a reasonable injection-induced seismic event for our study.

Fig. 5 presents a set of results related to the reactivation of the fault located at a depth of about 1000 m. For the assumed system and injection rate, the simulation showed that after more than 1 month of injection (~40 days), a dynamic fault rupture nucleates at the base of the CO₂ reservoir. At the initiation of the fault rupture, the fluid pressure had increased by about 7.5 MPa, i.e., achieving a total fluid pressure of 17.5 MPa, a few MPa higher than the local initial minimum principal stress at 1000 m depth. Moreover, as shown in Fig. 5b, the fluid pressurization of the reservoir caused a simultaneous increase in shear stress and reduction of effective normal stress, until the stress state reached the $\mu_s \cdot \sigma$ failure surface (quasi-static model phase in Fig. 5b). The shear stress on the fault increases at this monitoring point at the base of the reservoir, because of

pressure-induced vertical expansion of the reservoir to the left of the reservoir bounding fault.

After the stress-state reaches the $\mu_s \cdot \sigma$ failure surface, the fault reactivates with a sudden shear stress drop of about 1.4 MPa, as the friction drops to residual value and the shear stress equilibrates with the residual shear strength (dynamic phase in Fig. 5b). The fault reactivation induces localized plastic shear strain distributed over a length of about 290 m, with a maximum value of $4.5 \times 10^{-3}$, in a portion of the fault just below the reservoir (red zone in Fig. 5a). The full dynamic analysis showed that the rupture self-propagated as a result of fault-strength weakening, such that the slip propagated outside the pressurized zone (Fig. 5a). The reactivation occurred over about 0.4 seconds with a maximum co-seismic shear slip of about 4 cm (Fig. 5c). Based on the calculated co-seismic fault rupture length and average shear-slip displacement (290 m and 0.03 m), we quantified the overall size of the seismic event. A seismic moment, $M_0 = 3.48 \times 10^{13}$ Nm and a moment magnitude in the range of $M_w \approx 2.5–3$ were calculated using seismology relationships by Hanks and Kanamori (1979) and Kanamori and Anderson (1975) (Appendix A).

The simulated sudden fault rupture generated seismic waves that propagated and hit the ground surface at the speed corresponding to the P and S wave velocities for the medium of respectively $C_p = 2360$ m/s and $C_s = 1402$ m/s. In the dynamic simulation, the computational domain was limited laterally and at the base by absorbing conditions to avoid spurious reflections from boundaries. The mesh size and the time step were designed to accurately simulate the wave-propagation (Cappa and Rutqvist, 2012). The simulation includes energy dissipated by the plastic response within the fault as well as through interaction with the pore fluid, which reduces the ground accelerations compared to those calculated in a purely elastic medium. Under these conditions, no other artificial damping had to be imposed. Thus, the calculated wave attenuates by geometric spreading, pore fluid interaction, scattering, and by a substantial dissipation associated with plastic flow in the fault.

In order to verify our model simulation of ground motion, we used the simulated ground-motion wave train at several monitoring stations in an inverse analysis to estimate source parameters (moment magnitude, rupture dimensions and stress drop). The details of this analysis are presented in Appendix B showing a good agreement between simulated and inverted source parameters. In the source parameter analysis, the inverted Q-factor (inversely proportional to attenuation) was in the range of 20–30. These values are compatible with Q-factors observed in sedimentary rock reservoirs (Rutledge et al., 2004; Talebi and Bone, 1998; Godano et al., 2010). Although we do not rule out mesh effects on the very high frequency wave propagation results, the reasonable good agreement between the simulated and back-calculated source parameters show that the simulated ground motion is sensible. The resulting ground motion is in the next section used for a detailed analysis of potential damage of surface structures and nuisance to the local population.

## 4. Assessment of surface structural damage potential and human perception

In Figs. 6 and 7 we present acceleration and velocity recorded at 15 assumed stations along the ground surface (red and blue triangles in Fig. 4), for a complete picture of the ground motion and how it is distributed. We present both acceleration and velocity for the comparison with various criteria related to potential damage on structures and human perception. Moreover, we conducted a number of new simulations, including a somewhat softer soil layer at the ground surface that is either 50 or 100 m thick, to
<figure><img src="image_8.png" /><figcaption>**Fig. 5.** Simulations results after Cappa and Rutqvist (2012): (a) Plastic shear strain in the ruptured area and slip profile at the time of seismic rupture (white dots are control points), (b) shear stress-versus-effective normal stress path in the nucleation zone (control point 5 in (a)) of the rupture within the fault (the slow quasi-static phase is in blue and the dynamic phase is in red), (c) slip as function of time at different control points along the fault (white dots in (a)). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of the article.)</figcaption></figure>

investigate potential damping or amplification of the seismic response. The Young's modulus of this layer was set to be half that of the underlying rock formation. A more detailed view of acceleration and velocity records of two relevant stations (red triangles in Fig. 4) are presented in Figs. 8 and 9. Fig. 8 shows the results of a station located where the fault surface at the ground, a location where the highest acceleration and velocity were recorded, whereas Fig. 9 presents acceleration and velocity records for a station located at the approximate epicenter.

Fig. 6a shows that the high values of horizontal acceleration are concentrated around the fault (i.e. at about $x=673.3$ m and $x=675.8$ m), whereas the acceleration values quickly attenuate away from the fault. This is also clearly shown in Figs. 8a and 9a for a location just 300 m away from the fault where the horizontal acceleration is much smaller. We can also see that the high frequency accelerations are substantially damped in the case of a 50 or 100 m thick superficial soil layer (compare red and blue curves in Figs. 6, 8a and 9a). The same phenomenon of localized very high PGA values around faults has also been observed in the field (Wu et al., 2003; Ripperger et al., 2009). For example, Wu et al. (2003) reported, from their dense ground motion network in Taiwan, that large PGA values were usually characterized by a singular high values in space and time, with sharp peaks rich in high frequencies, and with the high values dropping off quickly outside certain neighborhoods. In our case, earthquake directivity, i.e., focusing or piling up of wave energy in the direction of rupture, may have contributed to the strong ground-motion response around the fault.

In general the results in Fig. 6 show a PGA less than 0.1 g, except for the horizontal component near the fault, where the

PGA is 0.6 g. Note that according to the USGS' empirically correlated instrumental intensity scale (see bottom Fig. 2), a PGA of 0.6 g would correspond to an instrumental intensity of VIII, which could give rise to perceived severe shaking and moderate-to-heavy damage. On the other hand, PGA's of about 0.1 g would correspond to an instrumental intensity of V, which could give rise to strong perceived shaking and light damage. However, as mentioned, those correlations were made for larger and deeper earthquakes involving lower frequencies of ground motion, whereas the dominant frequency in this case is 30–40 Hz for maximum acceleration magnitude (Cappa and Rutqvist, 2012). Such high-frequency ground motion is not expected to induce structural damage, although cosmetic damage might be possible and could certainly be felt by humans at least within a few kilometers from the epicenter.

Fig. 7 shows that the highest velocities are also associated with horizontal ground motion around the fault, with a maximum PGV of about 30 mm/s (Fig. 7a at $x=675.8$ m with PGV after about 0.7 s). As with PGA, the PGV magnitude also decreases for stations away from the fault being in the range of 2–10 mm/s, with the lowest values at the most distant stations. From the close-up view in Figs. 8b and 9b, we can see that the maximum PGV is associated with one main jolt having a period of about 0.1–0.2 s, i.e., a frequency of about 5–10 Hz. This main jolt is followed by smaller velocity peaks at somewhat higher frequencies. The soil layers have little effect on the main jolt, with some slight amplification in the case of thicker soil, whereas the subsequent higher frequency velocity peaks are substantially damped in the case of thicker soil. The peak velocity values of about 2–10 mm/s at stations away from the fault are consistent with observations of small and shallow earthquakes at field
a

<figure><img src="image_9.png" /></figure>

b

<figure><img src="image_10.png" /><figcaption>**Fig. 6.** Simulation results of (a) horizontal and (b) vertical component acceleration (expressed in $g=9.81 \text{ m/s}^2$) at 15 stations along the ground surface, with and without a softer top soil layer (red = no soil layer, green = 50 m thick soil layer, and blue = 100 m thick soil layer). Note the difference in y-axis scale in (a) and (b). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of the article.)</figcaption></figure>
<figure><img src="image_11.png" /></figure>

<figure><img src="image_12.png" /><figcaption>**Fig. 7.** Simulation results of (a) horizontal and (b) vertical component velocity at 15 stations along the ground surface with and without a softer top soil layer (red = no soil layer, green = 50 m thick soil layer, and blue = 100 m thick soil layer). Note the difference in y-axis scale in (a) and (b). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of the article.)</figcaption></figure>
<figure><img src="image_13.png" /><figcaption>**Fig. 8.** Close-up view of simulated horizontal (a) acceleration and (b) velocity component at x = 676 m (on top of fault) with and without a softer top soil layer (red = no soil layer, green = 50 m thick soil layer, and blue = 100 m thick soil layer). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of the article.)</figcaption></figure>

sites (Van Eck et al., 2006; Sharma et al., 2013). Comprehensive data sets for very shallow earthquake data are rare, but recently Sharma et al. (2013) derived PGV and PGA values as a function of distance for a large data set of earthquakes in the range of M<sub>w</sub> 0.5–5, at hypocenter distances of 0.5–15 km. Their results show average PGA value PGV values of about 0.1 g and 10 mm/s, respectively, for an M<sub>w</sub> = 3 event at a hypocenter distance of 1 km. The data in Sharma et al. (2013) also show that at hypocenter distances within a few km, local extreme values of PGA = 0.6 g and PGV = 0.05 m/s exist, i.e. values consistent with our simulated extreme near-fault ground

<figure><img src="image_14.png" /><figcaption>**Fig. 9.** Close-up view of simulated horizontal (a) acceleration and (b) velocity component at x = 348 m (approximate epicenter, about 300 m from fault at ground surface) with and without a softer top soil layer (red = no soil layer, green = 50 m thick soil layer, and blue = 100 m thick soil layer). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of the article.)</figcaption></figure>

motion. The frequency range of 5–10 Hz obtained in the simulations is also consistent with observations of shallow induced seismicity by Van Eck et al. (2006) and Bommer et al. (2006). Although, the resulting ground motion will always be site specific the good agreement with available field data shows that our simulated ground motion results are sensible (in terms of duration, magnitude and frequency) for a small, shallow event (i.e. M<sub>w</sub> = 3, at 1 km depth).

For a detailed comparison of the simulated PGV and its frequency content with the USBM building damage criterion and human perception, we conducted a frequency analysis of the

<figure><img src="image_15.png" /><figcaption>**Fig. 10.** Detailed frequency analysis of the velocity calculated at x = 673.324 m close to the fault for the case of no-soil layer: (a) simulated velocity, (b) filtered velocity between 5 and 50 Hz (Butterworth filter), (c) the spectrogram (frequency versus time with spectral amplitude in color, (d) a close-up view of the filtered velocity around the peak value, and (e) the frequency spectrum of the filtered velocity in the band between 5 and 50 Hz. The blue point in (4) and (5) corresponds to the peak value (around 1 s in the time domain) estimated at 0.028 m/s. This peak has a frequency of 6.85 Hz. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of the article.)</figcaption></figure>

<figure><img src="image_16.png" /></figure>
<figure><img src="image_17.png" /><figcaption>**Fig. 11.** Comparison of simulated PGV and frequency at different stations to the USBM criterion for cosmetic building damage and the human perception criterion.</figcaption></figure>

horizontal velocity record at the top of the fault. From the horizontal velocity records shown in Fig. 10a, we can see that relatively low frequency oscillations associated with the PGV are intermingled with higher frequency oscillations. This high-frequency oscillation gives rise to high accelerations that continue for several seconds (Fig. 10a), but with relatively small displacement, and with velocity values that are well below the cosmetic damage threshold. In order to obtain a good frequency analysis of the signal, one relevant for comparison to the Siskind curve, we filtered it in a frequency band near that of peak value. Using a Butterworth filter in signal processing, we filtered the data over a 5–50 Hz range and generated the filtered curve shown in Fig. 10b. Moreover, we constructed a spectrogram of frequency versus time and color-coded spectral amplitude in Fig. 10c. Finally, a close-up view of the filtered velocity in Fig. 10d and a frequency spectrum in Fig. 10e enabled us to make a precise identification of the PGV to 28 mm/s at a frequency of 6.85 Hz. However, Fig. 10d shows rather a plateau around 25 mm/s over a frequency 6–15 Hz.

In Fig. 11, we overlay the calculated frequency spectrum from Fig. 10e on the USBM criterion (Siskind curve) and the human perception criterion. We see that the simulated velocity of around 25 mm/s over a frequency of 6–15 Hz (green line in Fig. 11) is just above the blue Siskind curve and could therefore potentially cause cosmetic damage (e.g., hairline fracture in drywall or plaster). Moreover, the simulated velocity is above the red line for unpleasant steady-state vibration which implies that this short duration event would be clearly perceptible by humans. In Fig. 11 we also insert single peak velocity and frequency values at different stations along the ground surface. In this case, the frequency was evaluated directly from velocity record (Fig. 7) using the “zero-crossing” method, which is a method widely adopted by mining and construction industry for determining and reporting a single frequency value at the peak velocity of ground motions (Aimone-Martin et al., 2003). It is not as accurate as a full frequency analysis especially when the peak frequency occurs in a complex vibration time history containing a variety of frequencies and amplitudes. Nevertheless, the results in Fig. 11 show that for locations within a few hundred meters from the fault, the peak velocity is around 20–30 mm/s at a frequency of about 10 Hz. At

larger distances the peak velocity and frequency decreases and stabilizes at about 2.5 mm/s and 2.5 Hz for stations at distances from the fault exceeding 1 km ($x = 1500–2000$ m). According to Fig. 11, a PGV of 2.5 mm/s at 2.5 Hz might be barely perceptible by humans. Finally, we may also compare our simulated PGV against the USGS' empirically correlated instrumental intensity scale (bottom Fig. 2). Here a PGV of 30 mm/s (or 3 cm/s in Fig. 2 units) close to the fault, would correspond to perceived light shaking an no potential damage, whereas a PGV of 2–10 mm/s (or 0.2–1 cm/s) would be perceived as weak shaking. We can conclude that the USGS empirical PGV values are more consistent with the USBM damage criterion; whereas damage criterion based on USGS PGA values alone would substantial overestimate the potential damage and human perception.

Based on the exact location, the human perception of the event might be quite different. For example, to the left of the fault, there is a relatively strong z-component of velocity, associated the arrival of the P-wave, that gives rise to high frequency oscillations along with a downward movement occurring over a few seconds (see for example in Fig. 7b at $x = 486.8$ m where strong motions start already after about 0.2 s and in Fig. 7b where downward movement occurs from 0.2 s). Similar early P-wave response can be observed just to the right of the fault, but in this case with a systematic upward movement. The horizontal velocity components are very small during this early P-wave-dominated response. However, the highest velocity occurs close to the fault, with high values in both horizontal and vertical velocity—and at this location the ground motion would be distinctly perceptible by humans, perhaps perceived as unpleasant.

The simulation results are valid for the assumed model setting. Overall, the model setting adopted in this study is a critical case in which the fault can be reactivated to induce a seismic event. This include a set unfavorable conditions such as, (1) a reservoir bounding fault with pressurization on just one side of the fault, (2) a confined reservoir allowing for the pressure build-up, (3) a steeply dipping fault almost optimally oriented for shear reactivation relative to the in situ stress field, (4) a stress field with relatively low least principal oriented normal to the fault strike. Finally, (5) the fault was assumed brittle, meaning that the fault can reactivate in a seismic rather than aseismic manner. Other parameters such as stress ratio (Cappa and Rutqvist, 2012), the distance of the injection from the fault (Konstantinovskaya et al., 2014), fault reservoir thickness, caprock thickness and fault initial permeability (Rinaldi et al., 2014a) also impact the likelihood of reactivation and the magnitude and a seismic event. However, the reservoir bounding fault as a critical case was already mentioned by Hawkes et al. (2005), and have recently been shown quantitatively by Rinaldi et al. (2014b), who found higher likelihood for seismic response in a reservoir bounding fault compared to that of a reservoir crossing fault. This is not to say that a larger than $M_w$ 3 earthquake could never happen during geologic $CO_2$ sequestration operation. For example, the aforementioned $M_w$ 5.7 Prague earthquake (linked with nearby waste water injection) was associated with a reservoir bounding fault attached to the brittle crust (Keranen et al., 2013). Thus, the likelihood for induced seismicity will be site specific and a quantitative dynamic model analysis like the one demonstrated in this study could be applied to a real site, for example to study potential induced seismicity and ground motion associated with an identified critical fault.

Finally, we note that the seismic event calculated in our synthetic case would not cause any loss of caprock integrity or upward $CO_2$ leakage toward the ground surface as the fault reactivation took place below injection zone and not through the overlying caprock. This is in line with recent simulation results in Rinaldi et al. (2014b), who found poor correlation between the seismic events and $CO_2$ leakage, as relatively small-magnitude (between
2 and 3.5) events are not sufficient to substantially change the per- meability along the entire fault length. They found that a reservoir bounding fault with large shear offset was especially prone to cause injection-induced seismicity, but not to cause upward CO₂ leakage. However, we also note that the calculated ground motion and sur- face effects in this synthetic study are similar to that associated with the well-known 3.4 magnitude seismic event that occurred associated with a hydraulic stimulation operation at a geothermal project in Basel, Switzerland (Häring et al., 2008). The Basel event occurred in crystalline rock at a depth of 5 km, but the local ground- surface effects were similar, including slight nonstructural damage to buildings (e.g., hairline cracks in walls or minor pieces of plas- ter falling down—Ripperger et al., 2009). However, the magnitude 3.4 event was widely felt in the area of Basel, received nation- wide and international media attention and put the project on hold (Ripperger et al., 2009). One obvious lesson would be to place any injection site away from highly populated areas as even a small seis- mic event, even if not causing a loss of seal integrity, could likely unsettle the local population.

## 5. Concluding remarks

In this paper, we analyzed simulations of ground motions caused by CO₂-injection-induced fault reactivation in terms of the poten- tial for damage to ground surface structures and nuisance to the local human population. We conducted an integrated analysis from cause to consequence, including the whole chain of processes start- ing from earthquake inception in the subsurface, wave propagation toward the ground surface, to assess the consequences of ground vibration. From co-seismic average shear displacement and rupture area, we determined the moment magnitude to about M<sub>w</sub> = 3 for an injection-induced fault reactivation at a depth of about 1000 m. To verify our modeling approach, we used the resulting ground- motion wave train at several monitoring stations to conduct an inverse analysis by estimating source parameters (moment mag- nitude, rupture dimensions and stress drop), showing a reasonably good agreement. The good agreement between the synthetic and inverted source parameters verifies that we correctly model the chain of processes from earthquake inception to ground vibration. We then presented the results in terms of ground acceleration and velocity, and how these are distributed over multiple locations on the ground surface. In particular, we compared values of the peak ground acceleration (PGA) and peak ground velocity (PGV) to U.S. Geological Survey's instrumental intensity scales for earthquakes, and U.S. Bureau of Mines vibration criteria, as well as human per- ception vibration limits.

Our results showed (for the particular case simulated) locally strong PGA values up to 0.6 g close to the fault intersection with the ground surface, whereas PGA values away from the fault quickly attenuate to values below 0.1 g. A maximum PGV of about 30 mm/s also occurred close to the fault at the ground surface. This maximum PGV was associated with an initial pulse, one jolt, at a frequency of about 5–10 Hz, which (according to the USBM Siskind curve) could induce cosmetic damage in the form of hairline cracks in the plaster and drywall of buildings. The entire wave train from this 1000 m deep M<sub>w</sub> = 3 event lasted for about 1.5 s, but with apparent high-frequency free oscillations at the fault that continue for several seconds thereafter. Based on research into human tolerance of vibration, this event would likely be clearly felt and considered unpleasant by humans around the fault, where the PGV = 30 mm/s during one jolt, and veloci- ties oscillate above 10 mm/s for about 0.5 s. If one would look at the maximum PGA = 0.6 g, and the USGS' empirically correlated instrumental intensity scale (see bottom Fig. 2), a PGA of 0.6 g would correspond to an instrumental intensity of VIII, which would

correspond to perceived severe shaking and moderate-to-heavy damage. However, as previously discussed, the USGS' criteria were developed based on tectonic events, which occur much deeper than shallow injection-induced seismic events. Thus, our results confirm the appropriateness of using PGV (rather than PGA) and frequency for the evaluation of potential ground-vibration effects on structures and humans from shallow injection-induced seismic events.

Finally, in this study we demonstrated for a synthetic case that the integrated analysis involving the chain of processes from cause to consequences can be done quantitatively. At a future injection site, such an analysis will require site specific input parameters, including depth of injection zone, likely fault orientations and stress field, injection rates, and site specific material properties. For a dynamic wave propagation and ground vibration analysis it could also involve model calibration and testing of the velocity and atten- uation properties against seismic data (if available). One of the key properties for fault reactivation modeling is the coefficient of fric- tion of the fault and how it drops with shear. In our study, the coefficient of friction dropped from 0.6 to a residual value of 0.2 and we consider this as a conservative (worst case) scenario. A smaller difference between peak and residual shear strength sig- nificantly reduces the seismic response such that seismicity would easily drop to levels that would result in ground vibrations not dis- cernable by humans. This would be sites with more ductile (rather than brittle) rock properties that could accommodate the injection and yield to increasing reservoir pressure in an asseismic (rather than seismic) manner. Thus, in additional to placing a CO₂ injec- tion site away from dense populations, settings with more ductile rock properties would be preferable to minimize the risk of induced seismicity.

## Acknowledgments

This work was supported by the Assistant Secretary for Fossil Energy, Office of Natural Gas and Petroleum Technology, through the National Energy Technology Laboratory, under the U.S. Depart- ment of Energy Contract No. DE-AC02-05CH11231.

## Appendix A. Seismic moment and moment magnitude calculation

First, the scalar seismic moment M<sub>0</sub> defined for a ruptured patch on a fault by Hanks and Kanamori (1979):

$$M_0 = \mu Ad \quad (A1)$$

where μ is the shear modulus, A is the rupture area, and d is the average slip on the rupture area. Second, the moment magnitude (M<sub>w</sub>) of an earthquake is given, in terms of seismic moment, by Kanamori and Anderson (1975) as

$$M_w = \frac{2}{3}(\log_{10}(M_0) - 9.1) \quad (A2)$$

where the seismic moment, M<sub>0</sub>, is in Nm. Using Eqs. (A1) and (A2), we calculate a seismic moment of M<sub>0</sub> = 4 × 10<sup>9</sup> × 1000 × 290 × 0.03 = 3.48 × 10<sup>13</sup> Nm, and a moment magnitude of M<sub>w</sub> = 2.96. Note that our simulation result was achieved with a 2D plane-strain model, and in calculating the seismic moment and magnitude, we assumed that a 1 km section of the fault ruptured, leading to a rupture area of 1000 × 290 = 290,000 m<sup>2</sup>. This is a rupture geometry that could be related to (for example) a 1 km long horizontal injection well. Alter- natively, we may assume a circular rupture area with a diameter equal to the 290 m rupture length calculated from the plane strain analysis. Then, M<sub>0</sub> = 4 × 10<sup>9</sup> × π × 145<sup>2</sup> × 0.03 = 7.93 × 10<sup>12</sup> Nm, and the moment magnitude M<sub>w</sub> = 2.53.
## Appendix B. Inversion of earthquake source parameters from ground motion

The seismograms (velocity) recorded at 15 stations on the ground surface are shown in Fig. B1; these are the records we use for inverting source parameters. The analysis was done for the case without any shallow soil layer.

The frequency content of waves generated by a seismic event is empirically related to the magnitude of the event. First, the signal from three component sensors can be inverted for seismic moment and size of rupture by assuming a model for the source and relating that to an amount of slip. The Brune source model (Brune, 1970) stipulates that frequency response of the signal in displacement is flat until the corner frequency is reached, at which point the amplitude falls off as $f-2$ (although other models yield high-frequency asymptotes of $f-2.5$ or $f-3$). By fitting such a model to the displacement spectrum, we can estimate the low-frequency plateau $\Omega_0$ and the corner frequency $f_c$, and use them directly in the calculation of seismic moment and size of rupture area (source radius), respectively.

The general expression of the far-field displacement spectrum for P or S wave (Fig. B2) is:

$$U_{ij} = \frac{\Omega_0}{[1 + (f/f_c)^n]^l/n} \exp(-\pi f t_{ij}^*) \quad (B1)$$

where $\Omega_0$ is the amplitude of the spectrum low frequency part (plateau), $f_c$ is the corner frequency, $\gamma$ is the high frequency fall-off rate (on a log-log plot) and $n$ is a constant. Initially, Brune (1970) proposed a static rupture model for a circular fault with $\gamma = 2$ and $n = 1$. Boatwright (1980) proposed a modified version with $n = 2$, which was used in the present study. The Boatwright model produces a sharper corner frequency than the Brune model, which seems to generally better match observed data (Abercrombie, 1995). The parameter $t_{ij}^*$ is a measurement of the inelastic attenuation between an earthquake i and a seismological station j and is defined as:

$$t_{ij}^* = \frac{t_{ij}}{Q} \quad (B2)$$

<figure><img src="image_18.png" /><figcaption>Fig. B1. Seismograms (velocity) recorded at each of the 15 stations. The velocity records for the 15 stations are the same as shown in Fig. 15a, but normalized to the PGV to more clearly show the P wave arrival time, indicated with vertical bars. Inset: cross section displaying the position of the station (triangles) and the hypocenter. Station st05 and st06 have nearby location; only station st05 is represented on the cross section.</figcaption></figure>

<figure><img src="image_19.png" /><figcaption>Fig. B2. Observed (black) and theoretical (red) P-wave far-field displacement spectrum at station st01, st07, and st15. Corner frequency is indicated by the blue straight line. Plateau value is indicated by black dashed line. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of the article.)</figcaption></figure>

where $t_{ij}$ is the earthquake-station propagation time and $Q$ is the quality factor characterizing the inelastic attenuation of the medium.

In this study, we will consider the P-wave exclusively. The P-wave is preferable to the S-wave for spectral analysis, because P-waves are not affected by prior waves—contrary to the S-wave, which can be affected by converted waves or by the P-wave itself if the delay time between P- and S- is short (as in this case).

The far-field displacement spectrum $U_{ij}^{\text{obs}}$ is computed at each station for a time window (0.4 s) containing the P-wave. For 3-component recordings, the displacement spectrum is obtained by summing the components (in our case, we have only two components):

$$U_{ij}^{\text{obs}} = \sqrt{(U_{ij}^{X})^2 + (U_{ij}^{Y})^2 + (U_{ij}^{Z})^2} \quad (B3)$$

The parameters $\Omega_0, f_c$, and $t_{ij}^*$ are determined next, by inverting the observed spectrum. In practice, an iterative exploration is made
Eshelby, J.D., 1957. The determination of the elastic field of an ellipsoidal inclusion and related problems. Proceedings of the Royal Society, London, Series A 241, 376–396.

Fletcher, J.B., Boatwright, J., Joyner, W.B., 1983. Depth dependence of source parameters at Monticello, South Carolina. Bulletin of the Seismological Society of America 73, 1735–1751.

Godano, M., Gaucher, E., Bardainne, T., Regnier, M., Deschamps, A., Valette, M., 2010. Assessment of focal mechanisms of microseismic events computed from two three-component receivers: application to the Arkema-Vauvert field (France). Geophysical Prospecting 58, 775–790.

Hanks, T.C., Kanamori, H., 1979. A moment magnitude scale. Journal of Geophysical Research 84, 2348–2350.

Hanks, T.C., Wyss, M., 1972. The use of body-wave spectra in the determination of seismic source parameters. Bulletin of the Seismological Society of America 62, 561–589.

Hawkes, C.D., McLellan, P.J., Bachu, S., 2005. Geomechanical factors affecting geological storage of CO₂ in depleted oil and gas reservoirs. Journal of Canadian Petroleum Technology 44, 52–61.

Häring, M.O., Schanz, U., Ladner, F., Dyer, B.C., 2008. Characterisation of the Basel 1 enhanced geothermal system. Geothermics 37, 469–495.

IEAGHG, 2013. Induced seismicity and its implication for CO₂ storage risk. Report 2013/9. International Energy Agency Greenhouse Gas R&D Program (IAEGHG), Cheltenham, UK.

Ikari, M.J., Saffer, D.M., Marone, C., 2009. Frictional and hydrologic properties of clay-rich fault gouge. Journal of Geophysical Research 114, B05409.

Itasca, 2011. FLAC3D V5.0, Fast Lagrangian Analysis of Continua in 3 Dimensions, User's Guide. Itasca Consulting Group, Minneapolis, MN.

Kanamori, H., Anderson, D.L., 1975. Theoretical basis of some empirical relations in seismology. Bulletin of the Seismological Society of America 65, 1073–1095.

Keranen, K.M., Savage, H.M., Abers, G.A., Cochran, E.S., 2013. Potentially induced earthquakes in Oklahoma, USA: links between wastewater injection and the 2011 Mw 5.7 earthquake sequence. Geology, Published online on 26 March 2013 as doi:10.1130/G34045.1

Konstantinovskaya, E., Rutqvist, J., Malo, M., 2014. CO₂ storage and potential fault instability in the St. Lawrence Lowlands sedimentary basin (Quebec, Canada): insights from coupled reservoir-geomechanical modeling. International Journal of Greenhouse Gas Control (in press).

Madariaga, R., 1976. Dynamics of an expanding circular fault. Bulletin of the Seismological Society of America 66, 639–666.

Majer, E., Nelson, J., Robertson-Tait, A., Savy, J., Wong, I., 2012. Protocol for addressing induced seismicity associated with enhanced geothermal systems. U.S. Department of Energy, January 2012; DOE/EE-0662.

Mazzoldi, A., Rinaldi, A.P., Borgia, A., Rutqvist, J., 2012. Induced seismicity within geologic carbon sequestration projects: maximum earthquake magnitude and leakage potential. International Journal of Greenhouse Gas Control 10, 434–442.

McGarr, A., 1984. Scaling of ground motion parameters, state of stress, and focal depth. Journal of Geophysical Research 89, 6969–6979.

National Research Council, 2012. Induced Seismicity Potential in Energy Technologies. National Academies Press, Washington, DC, pp. 300.

Orlic, B., 2009. Some geomechanical aspects of geological CO₂ sequestration. KSCE Journal of Civil Engineering 13, 225–232.

Pruess, K., Oldenburg, C., Moridis, G., 2011. TOUGH2 User's Guide, Version 2.1, LBNL-43134 (revised). Lawrence Berkeley National Laboratory, Berkeley, CA.

Rinaldi, A.P., Rutqvist, J., 2013. Modeling of deep fracture zone opening and transient ground surface uplift at KB-502 CO₂ injection well, In Salah, Algeria. International Journal of Greenhouse Gas Control 12, 155–167.

Rinaldi, A.P., Jeanne, P., Rutqvist, J., Cappa, F., Guglielmi, Y., 2014a. Effects of fault-zone architecture on earthquake magnitude and gas leakage related to CO₂ injection in a multilayered sedimentary system. Greenhouse Gases: Science and Technology 4, 99–120.

Rinaldi, A.P., Rutqvist, J., Cappa, F., 2014b. Geomechanical effects on CO₂ leakage through fault zones during large-scale underground injection. International Journal of Greenhouse Gas Control 20, 17–31.

Ripperger, J., Kästli, P., Fäh, D., Giardini, D., 2009. Ground motion and macroseismic intensities of a seismic event related to geothermal reservoir stimulation below the city of Basel—observations and modelling. Geophysical Journal International 179, 1757–1771.

Rutledge, J.T., Phillips, W.S., Mayerhofer, M.J., 2004. Faulting induced by forced fluid injection and fluid flow forced by faulting: an interpretation of hydraulic-fracture microseismicity. Carthage Cotton Valley Gas Field, Texas. Bulletin of the Seismological Society of America 94, 1817–1830.

Rutqvist, J., 2011. Status of the TOUGH-FLAC simulator and recent applications related to coupled fluid flow and crustal deformations. Computers & Geosciences 37, 739–750.

Rutqvist, J., 2012. The geomechanics of CO₂ storage in deep sedimentary formation. Geotechnical and Geological Engineering 30 (3), 525–551.

Rutqvist, J., Wu, Y.-S., Tsang, C.-F., Bodvarsson, G.A., 2002. Modeling approach for analysis of coupled multiphase fluid flow, heat transfer, and deformation in fractured porous rock. International Journal of Rock Mechanics and Mining Sciences 39, 429–442.

Rutqvist, J., Birkholzer, J., Cappa, F., Tsang, C.-F., 2007. Estimating maximum sustainable injection pressure during geological sequestration of CO₂ using coupled fluid flow and geomechanical fault-slip analysis. Energy Conversion Management 48, 1798–1807.

Rutqvist, J., Vasco, D.W., Myer, L., 2010. Coupled reservoir-geomechanical analysis of CO₂ injection and ground deformations at In Salah, Algeria. International Journal of Greenhouse Gas Control 4, 225–230, http://dx.doi.org/10.1016/j.ijggc.2009.10.017.

Samuelson, J., Spiers, C.J., 2012. Fault friction and slip stability not affected by CO₂ storage: evidence from short-term laboratory experiments on North Sea reservoir sandstones and caprocks. International Journal of Greenhouse Gas Control 11S, S78–S90.

Siskind, D.E., Stagg, M.S., Kopp, J.W., Dowding, C.H., 1980. Structure response and damage produced by ground vibrations from surface blasting. Report of Investigations 8507. US Bureau of Mines, Washington, DC.

Sharma, N., Convertito, V., Maercklin, N., Zollo, A., 2013. Ground-motion prediction equations for the Geysers geothermal area based on induced seismicity records. Bulletin of the Seismological Society of America 103, 117–130.

Svinkin, M.R., 2004. Minimizing construction vibration effects. Practice Periodical on Structural Design and Construction 9, 108–115.

Talebi, S., Bone, T., 1998. Source parameters of injection-induced micro seismicity. Pure and Applied Geophysics 153, 113–130.

USACE, 1972. Systematic drilling and blasting for surface excavations. Engineering Manual EM 1110-2-3800. US Army Corps of Engineers.

Van Eck, T., Goutbeek, F., Haak, H., Dost, B., 2006. Seismic hazard due to small-magnitude, shallow-source, induced earthquakes in The Netherlands. Engineering Geology 87, 105–121.

Wald, D.J., Quitoriano, V., Heaton, T.H., Kanamori, H., 1999. Relationships between peak ground acceleration, peak ground velocity and modified Mercalli intensity in California. Earthquake Spectra 15 (3), 557–564.

Worden, C.B., Wald, D.J., Rhoades, D.A., 2012. Probabilistic relationships between ground-motion parameters and modified Mercalli intensity in California. Bulletin of the Seismological Society of America 102, 204–221.

Wu, Y.-M., Teng, T.-L., Shin, T.-C., Hsiao, N.-C., 2003. Relationship between peak ground acceleration, peak ground velocity, and intensity in Taiwan. Bulletin of the Seismological Society of America 93, 386–396.

Zoback, M.D., 2007. Reservoir Geomechanics. Cambridge University Press, Cambridge, UK.

Zoback, M.D., Gorelick, S.M., 2012. Earthquake triggering and large-scale geologic storage of carbon dioxide. Proceedings of the National Academy of Sciences, http://dx.doi.org/10.1073/pnas.1202473109.

Zoback, M.L., Zoback, M.D., 1980. State of stress in the conterminous United States. Journal of Geophysical Research 85, 6113–6156.