(432a) Managing a Gigaton Carbon Sequestration Future: A Framework for Prediction of Fault Activation Caused By Basin-Scale CO2 Storage

After decades of research on geologic carbon storage, the world is finally moving from pilot tests and demonstration experiments to industrial-scale implementation. Many large-scale CO₂ storage hubs are currently in development across the world, notably in the United States, Europe and China. Even further expansion of carbon capture and sequestration to perhaps a gigaton scale is needed to achieve a net zero clean energy future. Sedimentary basins with interconnected reservoirs will likely host multiple CO₂ storage sites between which pressure interference can be expected. The cumulative effect of multiple CO₂ projects in a given region may limit its overall storage capacity, may enlarge the project area that needs to be characterized and monitored, and may cause unwanted geomechanical effects such as generating seismic events and caprock failure per reactivation of critically stressed faults. Such risks need to be carefully managed.

Here, we focus on regional-scale CO₂ storage management with particular emphasis on geomechanical risks and discuss our ongoing development of a versatile framework for simulation-driven storage optimization at a regional scale. First, we use coarse poroelastic flow simulations across an entire basin (with multiple CO₂ storage projects) to identify potentially critically-stressed faults. Second, we simulate rupture on these faults using coupled hydromechanical simulations of fault slip at the scale of individual CO₂ storage projects. Finally, we use results from the project-scale simulations to provide parameters for sophisticated simulations of dynamic rupture propagation using rate and state friction and Cam-Clay models.

The state-of-the-art for assessment of induced seismicity and caprock failure risks relies on overpressure as the main control for the fault stability via the Coulomb failure criterion. Such a simplified view disregards several key mechanisms. First, poroelastic stresses triggered by waste water injection proved to be an important factor driving the increase of seismicity in Oklahoma. Second, a few recent field tests, including the controlled injection tests at Mont Terri (Switzerland), showed a clear relationship between the seismic/aseismic fault slip and permeability increase in the caprock faults. Third, multiphase effects that accompany CO₂ filtration through reservoir faults may boost their susceptibility to slip. Our simulation-driven framework aims to incorporate these physical mechanisms explicitly across a broad range of scales. In the following, we present a hypothetical basin model along with several injection scenarios, which highlight the importance of project-to-project interaction as well as accurate upscaling of the fault physics.

Simplified Poroelastic Basin-scale Simulations

The key challenge of regional-scale optimization arises from the fact that pressure/stress changes in a given basin occur over several years/decades and have characteristic dimensions of tens of kilometers, while fault reactivation is rapid and local, and thus requires frequent time steps as well as accurate representation of the fine-scale fault geometry and rheological properties. We implement a number of practical strategies to bridge the gap between the basin-scale computations of pressure and stress changes and sophisticated coupled hydromechanical simulations of localized fault rupture. The simplest approach locates the segments of faults with sufficiently high Coulomb stress on the faults in the coarse-scale poroelastic simulations. These segments are then passed through a suite of project-scale coupled hydromechanical simulations with the boundary and initial conditions extracted from the coarse-scale simulations.

We illustrate the coarse-scale poroelastic simulations using a synthetic benchmark model shown in Figure 1a-c. The model is based on the local geology in the San Francisco Bay Area and some scaling relationships between the fault parameters, such as throw, length, spacing. This purely hypothetical study allows us to test optimizing the computational parameters as well as quantifying the impact of various physical effects, such as 2-phase flow versus single-phase injections, boundary conditions, hydraulic connectivity to the basement, and the injection schedules for different projects. We find for example that only one fault becomes critically stressed for a scenario with constant injection volumes and simultaneous injection start of five CO₂ storage projects (Figure 1d). We analyze the potential for the failure of this fault in the next section.

Advanced Simulations of a High-Resolution Rupture Patch for Complex Faults

On the next level of detail, we use the TOUGH-FLAC simulator to capture the reactivation of the fault with complex geometry highlighted in Figure 1d. Figure 2a shows the computational mesh for TOUGH-FLAC simulations with a complex fault geometry. The size of the mesh is 5 km x 5 km x 1.4 km. Three layers (caprock, reservoir, and basement ordered from the top to bottom) are used. We consider injection of 0.5 Mt of supercritical CO₂ injection at 2,100 m depth into the hanging wall of a normal fault, which is offset from the injection by 140 m. After 2 years, the rupture nucleates at the boundary between the reservoir bottom and the underlying low-permeability rocks (Figure 2b), while red indicates past reactivation after 20 years of injection (Figure 2c). A majority of the slipped elements gets activated within the first 2 years of the injection. The rupture nucleates at the top of the basement, below the bottom of the storage reservoir. As more CO₂ is injected, the rupture propagated upwards a little bit, but it mainly spreads horizontally: about 1,400 m along the fault strike and about 140 m up-dip. Caprock integrity is thus not affected. Interestingly, in the case of an injection in the foot wall of the fault, the rupture occurs at the top of the reservoir and tends to extend horizontally and down-dip.

Estimating Rupture Acceleration to Seismicity

Thus far, we have considered a Coulomb-Mohr failure criterion on the faults, which precludes estimation of the magnitude of the triggered events. Here, we derive the range of seismic magnitudes using highly complex models for rupture acceleration to seismicity on the relatively small fault zone modeled in Figure 2. We compare two theoretical frameworks for the fault stability (Figure 3): (1) the Rate and State Friction model (RSF) and the Modified Cam-Clay model (MCC). In both models, the injection into the hanging wall of the fault creates a 5 MPa pore pressure increase in the fault zone after 2 years. Poroelastic stress changes induce a small 8 mm/year tangential displacement on the fault and 3 mm/year normal compression (Figure 3a). Preliminary results are: (1) the RSF and MCC model predictions provide similar estimates of fault slip acceleration and amplitude (Figure 3b-e: interval 3-4 of plastic instability, 15cm slip, 30-40cm/s); (2) the RSF model underestimates the preceding slow slip period which could have strong implications on fault permeability variation associated to dilatant rupture (Figure 3d-e: interval 2-3 of plastic stability, 30 cm slip with MCC and 8mm slip with RST); (3) the MCC is better suited to consider the effects of fault bulk plastic deformation on slip and leakage.

Conclusions

We present a workflow for simulation-driven assessment of risks of induced seismicity and fault leakage potentially caused by basin-scale CO₂ storage. The workflow starts with simplified multiphase poroelastic flow simulations to detect potentially critically-stressed faults, then zooms in on the detailed hydromechanical rupture behavior of these faults with advanced coupled hydromechanical simulations, and finally estimates seismic magnitudes of relevant fault slip scenarios with various constitutive laws for the dynamic rupture propagation. We propose this workflow as a practical strategy for multi-scale simulation of both slow and seismogenic slip on faults as well as estimation of the resultant magnitude.