#### Abstract

The viability of carbon capture sequestration (CCS) is dependent on the secure storage of CO_{2} in subsurface geologic formations. Geomechanical failure of caprock is one of the main reasons of CO_{2} leakage from the storage formations. Through comprehensive assessment on the petrophysical and geomechanical heterogeneities of caprock, it is possible to predict the risk of unexpected caprock failure. To describe the fracture reactivation, the modified Barton–Bandis model is applied. In order to generate hydro-geomechanically heterogeneous fields, the negative correlation between porosity and Young’s modulus/Poisson’s ratio is applied. In comparison with the homogeneous model, effects of heterogeneity are examined in terms of vertical deformation and the amount of leaked CO_{2}. To compare the effects of heterogeneity, heterogeneous models for both geomechanical and petrophysical properties in coupled simulation are designed. After 10-year injection with petrophysically heterogeneous and geomechanically homogeneous caprock, CO_{2} leakage is larger than that of the homogeneous model. In contrast, heterogeneity of geomechanical properties is shown to mitigate additional escape of CO_{2}. Vertical displacement of every heterogeneous model is larger than homogeneous model. The model with compressive tectonic stress shows much more stable trapping with heterogeneous caprock, but there is possibility of rapid leakage after homogeneous caprock failure.

#### 1. Introduction

CO_{2} sequestration in aquifer is an effective and a verified method to reduce the atmospheric CO_{2} airborne fraction. The world’s first commercial approach of CO_{2} storage project is attempted at the Sleipner gas field in the North Sea, and CO_{2} sequestration in a saline aquifer has been regarded as the feasible technology [1]. For the geological storage of CO_{2}, there should be an impermeable caprock which acts as a physical barrier to prevent the leakage and migration of the buoyant CO_{2} [2]. One of the considerable risks is the leak of CO_{2} incurred by a mechanical failure of caprock. Injected CO_{2} can leak out and flow throughout the opened fractures and faults of caprock; this is an important point in the caprock’s sealing mechanism. Thus, it is important to monitor the geomechanical change of storage formation to evaluate the stability of CO_{2} injection. In the aftermath of large volume of CO_{2} injection, fluid pressure in the storage formation will rise, and increase of mechanical stresses will also occur. If injected CO_{2} induce too high fluid pressure, the induced stresses will exceed the summation of least principal effective stress and tensile strength and affect irreversible mechanical displacements of storage rock, reactivating natural fractures or creating additional fractures. Geological displacements of CO_{2} in storage formation have been observed at Weyburn CO_{2} monitoring and storage project [3], and surface uplift caused by increased subsurface pressure has been detected at the In Salah project [4]. These changes could generate new high-permeable zone as leakage path, thereby substantially decreasing the efficiency of CO_{2} storage.

To ensure the feasibility of potential storage site, intensive estimation of formation properties should be preceded [5]. Hurter et al. [6] emphasized the importance of numerical simulation as preinjection studies. Several numerical simulation studies have been conducted to assess the stability of CO_{2} injection with geomechanical module. Lucier et al. [7] applied a geomechanical analysis of a potential injection zone in the Ohio Valley region to evaluate the potential risks for injection-induced stresses. Rutqvist et al. [8] conducted coupled reservoir-geomechanical simulations to evaluate the shear and tensile failure in order to accurately estimate the potential of fracture reactivation for a geomechanical CO_{2} storage formation, and they suggested that long-lasting observation is needed to geomechanical changes resulting from the injected fluid pressure. Park et al. [9] evaluated the effects of geomechanical properties on the fracture reactivation by using two-dimensional simulation model, and they also examined the effects of the heterogeneity of geomechanical properties during CO_{2} storage.

For a much more accurate evaluation of geomechanical change, the coupled mechanism of geomechanical and petrophysical properties was considered. Cappa and Rutqvist [10] tried to reveal the interaction between fluid flow and mechanical deformation during CCS and demonstrated how to analyze coupled hydromechanical processes. Khajeh et al. [11] investigated the impact of coupled geomechanical and petrophysical properties including heterogeneity and range of uncertainties. Previous researches have considered this coupling, but they did not focus on the correlation among these properties that could affect the accuracy of numerical modeling.

In this study, we constructed a CCS model with heterogeneous caprock and compared the flow and geomechanical responses to CO_{2} injection with particular focus on the risks to storage security posed by geomechanical deformation. The correlations between geomechanical and petrophysical properties were regarded as a function of porosity, and several numbers of numerical simulations were conducted to evaluate the effects of the heterogeneous caprock during CCS. The amount of stored CO_{2} and vertical deformation of caprock within a two-dimensional formation models were compared for several rock properties including parameters such as Poisson’s ratio, Young’s modulus, porosity, and permeability. These properties were distributed in randomly generated heterogeneity conditions. The effects of heterogeneity on caprock were revealed by comparing the average values of heterogeneous model simulation results with homogeneous caprock model. Considering correlation between geomechanical and petrophysical properties can help assessing formation stability and mitigating geological risk. Additionally, insitu tectonic stress is known to play an important role in fracture propagation [12, 13]. Horizontal tectonic stress also influences the porosity and permeability of formation [14]. However, there appears to be a lack of information regarding the influence of in situ tectonic stress on long-time CO_{2} storage performance. In order to address this issue, an investigation of heterogeneous permeability fields including additional horizontal stress has been conducted.

#### 2. Theoretical Background

##### 2.1. Fracture Reactivation

The overburden stress, known as the total stress, induced by external loading, is sustained by the rock matrix and the formation fluid pressure. To evaluate the fluid pressure required to fracture reactivation in caprock, Terzaghi [15] suggested the law of effective stress as follows:where is the total stress (kPa), is the effective stress (kPa), and is the fluid pressure (kPa).

The tensile failure potential is affected by applying the assumption that tensile fracture propagation can occur when the fluid pressure exceeds the summation of least principal effective stress and tensile strength [16]:where is the critical fluid pressure (kPa), is the least principal effective stress (kPa), and is the tensile strength (kPa). Most of sedimentary rocks have a relatively low tensile strength, typically few MPa or less. Fjaer et al. [17] indicated that it is a standard approximation for numerous applications that the tensile strength is zero. In the case that the pressure build-up by injected fluid exceeds the least principal effective stress, the fracture propagates in the direction being perpendicular to least principal effective stress. This process is similar to the tensile reactivation of natural fractures.

##### 2.2. Modification of Fracture Permeability

A hydromechanical fracture permeability model was based on the modified Barton–Bandis model, as shown in Figure 1 [18]. Initially, effective stress () is much larger than fluid pressure (point A), and there are no activated fractures in the rock matrix. The effective stress is decreased along the CO_{2} injection (path AB). On this path, the range of fracture permeability is small, and the reactivation is reversible. When the effective stress is below point B, the fracture opens abruptly, and the fracture permeability increases drastically from the initial value to maximum value () on path BC. The threshold value (point B) is named the fracture opening stress (). Although the effective stress decreases consistently, the fracture permeability keeps its maximum value (path DE). Effective stress increment and the fracture aperture reduction due to decline of fluid pressure lead to the decrease of fracture permeability. However, fracture permeability change does not occur on the previous path BC. Instead, it follows a new hyperbolic curve (path FG) because the fracture cannot be closed completely. Fracture closure permeability is always bigger than the residual value of fracture closure permeability (). If the increment of fluid pressure induces another reduction of the effective stress, fracture permeability will increase along the path GFED to the maximum fracture permeability. Fracture permeability () can be calculated as follows [18, 19]:where is the fracture closure permeability (md), is the initial fracture aperture (m), is the joint closure under a normal fracture effective stress, is the initial normal fracture stiffness (kPa/m), and is the maximum fracture closure (m).

##### 2.3. Permeability-Porosity Relationship

Walsh and Brace [20] indicated that permeability is regarded as proportional to simple integer powers of the related pore geometry parameters such as hydraulic radius, porosity, and tortuosity. These parameters are normally regarded to be related to each other through the power-law relationship (), leading to power-law dependence of permeability on porosity [21]. Experimental evidence indicates that power-law exponent does not always match porosity changes. However, for elastic or plastic compaction, chemical alteration, or microcracking, it is usually regarded that there is a power-law relationship between these parameters [22].

##### 2.4. Porosity-Geomechanical Properties Relationship

To analyse the petrophysical and geomechanical properties, sedimentary rocks were studied using samples cored from several wells from the Krishna-Godavari basin on the east coast of India [23]. Figure 2 shows best fit linear curves for uniaxial compressive strength, effective porosity, and Young’s modulus () of core samples. As could be seen that the linear relationship between Young’s modulus strength, and uniaxial compressive strength has reasonable correlation coefficients ( = 0.93). Similarly, the correlation coefficients of the best fit linear curves between uniaxial compressive and effective porosity are high ( = 0.89). Using linear equations, it is possible to predict the geomechanical properties of sedimentary rocks from petrophysical properties. The results of the regression analysis indicate that the petrophysical and geomechanical properties have a positive correlation to each other. Figure 3 shows the relationship between the geomechanical properties of and Poisson’s ratio () [24]. Rock physics templates show the mapping of lines of mineralogical mixtures and constant porosity on the crossplots of geomechanical properties [25]. Solid lines represent the clay content, and dashed lines indicate constant values of porosity. The intersection of these lines is overlaid with points, representing chosen values of porosity, from 0 to 0.20. Each color of dots represents constant clay content, explained by the legend. From each constant porosity line, it can be seen that the increase in porosity will decrease the values of and . This may imply that porosity has negative correlation with and .

**(a)**

**(b)**

**(a)**

**(b)**

#### 3. Numerical Simulation

GEM was used to assess the fracture propagation during the CO_{2} injection, which is a compositional simulator developed by the Computer Modeling Group [18]. The geomechanics module in GEM was applied to compute the deformation in the rock and the stress change. All the rocks in the matrix were assumed to be elastic during the injection process. To describe the fracture permeability change, dual permeability system and Barton–Bandis model were used. For the in situ stress, the weight of the overburden rock causes the vertical stress. The initial fluid pressure was induced by the hydrostatic gradient of 9.81 MPa/km with an atmospheric pressure of 0.1 MPa from the surface. Constant pressure condition is supposed at the horizontal outer boundary to minimize the boundary effects. All of the simulations were organized in a constant temperature. The static geomechanical and hydraulic properties of the base case are shown in Table 1, with the typical rock types of the permeable caprock and the formations representing sandstone and shale [26].

##### 3.1. Model Description

In order to assess the effects of the geomechanical rock properties of the fracture reactivation, several simulations were conducted on two-dimensional cross-sectional model. The model considers permeability and porosity of rock matrix, Poisson’s ratio, and Young’s modulus, which are regarded to be the most important properties for the investigation of the geomechanical state. The formation model was constructed, as shown in Figure 4, including an aquifer, caprock, and overburden rock. The formation model locates vertically from 1,000 to 1,330 m depth. Model of 2,020 20 330 m^{3} was divided into 101 1 33 grid blocks. Pure CO_{2} was injected through an injection well at 1,300 m depth and a constant rate of 120,000 m^{3}/day for 10 years. Table 2 shows component properties of CO_{2} used in the simulation [27]. CO_{2} existed as a supercritical state at the given temperature and pressure ( = 50°C and ). At first, the injected CO_{2} was sequestered under the caprock and then leaked through the fracture in the caprock after the fracture reactivation occurred. During the simulation, the vertical deformation and effective stress were measured.

Porosity-permeability relationship is calculated using a power law. To generate geostatistical hydrogeomechanically heterogeneous fields, the negative correlation between porosity and Young’s modulus is applied. Poisson’s ratio is negatively related to porosity, too. Table 3 shows physical data of geological caprocks from saline CO_{2} storage in the Paradox Basin [28–30]. Geomechanical properties of simulation values are modeled based on this formation data.

##### 3.2. Heterogeneity Description

In this study, a stochastic simulation method is used for heterogeneous reservoir model. To describe formation heterogeneity, Dykstra–Parsons coefficient (), which is the most common measure of heterogeneity in petroleum industry, is used [31, 32]. It is defined aswhere is the median permeability value and is the permeability at 84.1% probability (one standard deviation). The model with same also has the same correlation length (Table 4). Averages of all variables are same as homogeneous model (Table 1). Each model has same mean and variance of permeability, but all of the realizations have their own permeability distribution. The result graphs are drawn from the average of 25 realizations. The number of realizations was chosen after testing the sensitivity of results. Though variations were sometimes observed, 25 were chosen for the final number of realizations for analysis due to computational constraints [11].

To evaluate how heterogeneity considerations for petrophysical and geomechanical properties affect output results, four cases are designed as below: Case A: heterogeneous petrophysical and homogeneous geomechanical properties Case B: homogeneous petrophysical and heterogeneous geomechanical properties Case C: heterogeneous petrophysical and geomechanical properties Case D: homogeneous petrophysical and geomechanical properties

The effects of heterogeneity are revealed in comparison to homogeneous model. The results are compared in terms of vertical deformation and the amount of leaked CO_{2}.

#### 4. Results and Discussion

##### 4.1. CO_{2} Storage

Figure 5 shows gas saturation of single realization and CO_{2} leakage regime. The figure is about Case C which has petrophysically heterogeneous properties and geomechanically heterogeneous properties and Dykstra–Parsons coefficient of 0.7. Increased formation pressure causes fractures in caprock, and CO_{2} escapes from the aquifer to overburden rock. Figure 6 depicts the cumulative amount of remaining CO_{2} in aquifer during 10-year injection. These graphs show the total amount of gas sequestered in the aquifer. Each curve is generated by averaging results from 25 heterogeneous models. The rest of gas has escaped to overburden formation. There could be further escape after injection, until the system reaches equilibrium. Figure 6(a) indicates simulation model with petrophysically heterogeneous and geomechanically homogeneous properties. By investigating the results of Case A, it can be seen that heterogeneous model shows less CO_{2} moles in the aquifer than homogeneous model. Also, the amount of leaked CO_{2} depends on heterogeneity. In other words, more petrophysically heterogeneous model presents more CO_{2} leakage. In heterogeneous caprock, there could be relatively weak shale which has high porosity and permeability. Injected fluids tend to flow through high permeable zone and partially increase formation pressure. Concentration of pressure incurs fracture reactivation and earlier leakage of CO_{2}. Figure 6(b) shows simulation model with petrophysically homogeneous and geomechanically heterogeneous properties. In early stage of injection, heterogeneous models show more leakage than homogeneous model. However, as time goes on, the leakage is mitigated compared to homogeneous model. This difference is caused by imbalance of geomechanical properties. Being different from the heterogeneous caprock, the failure of homogeneous caprock occurs abruptly at the center. In heterogeneous caprock, the leakage does not occur above the injector. Initially, CO_{2} gas leaks through the location region with high Young’s modulus and Poisson’s ratio because this spot could not elastically disperse the stress. As injection continues, further escape of CO_{2} is mitigated by the heterogeneity of caprock. Compared to homogeneous model, heterogeneous model shows unstable increase of cumulative CO_{2} moles. Also, more heterogeneous model has larger difference compared with homogeneous model. Figure 6(c) shows the cumulative CO_{2} moles in petrophysically heterogeneous and geomechanically heterogeneous caprock. This model also shows more leakage in more heterogeneous model at early stage of injection, but after 4 years of injection, more heterogeneous model shows slightly more CO_{2} storage under the influence of geomechanical heterogeneity. After that, each model reveals similar leakage along continuous injection due to integrated effect of heterogeneous properties. Heterogeneous models show less difference compared to above two cases due to combined effect of petrophysical and geomechanical properties.

**(a)**

**(b)**

**(c)**

In order to assess the range on uncertainty, a box-whisker plot format is used to illustrate minimum, 25th, 50th, 75th percentiles, and maximum values for sequestered CO_{2} after 10-year injection. As could be seen from Figure 7(a), the range of uncertainty is getting larger along the Dykstra–Parsons coefficient in Case A, and the amount of sequestered CO_{2} is getting smaller compared with homogeneous model. On the contrary to Figures 7(a) and 7(b) shows that geomechanical heterogeneity has favorable effect in CO_{2} trapping. Similar to Case A, the range of uncertainty is also increasing with Dykstra–Parsons coefficient getting bigger. Likewise, Case C shows the increase of uncertainty range as shown in Figure 7(c).

**(a)**

**(b)**

**(c)**

Figure 8, which shows the box-whisker plot comparing four cases, presents the effect of heterogeneity of each property after 10-year injection for most heterogeneous models ( = 0.7). Case A indicates 3.9% lower sequestered CO_{2}. Heterogeneities of permeability and porosity have a negative influence on CO_{2} trapping because CO_{2} leakage occurs through the highly permeable shale first in petrophysically heterogeneous formation. However, in geomechanically heterogeneous model, leakage also occurs through the weak point in first 4 years. Then, geomechanical heterogeneity mitigates further escape of CO_{2} because geomechanically heterogeneous caprock includes elastic formation which has low Young’s modulus and Poisson’s ratio, and this section helps caprock endure the stress elastically. After 10 years of injection, Case B shows 4.6% higher amount of stored CO_{2}. Due to combined effects of these two properties, Case C reveals only 0.3% difference compared with homogeneous case because composite effects of geomechanical and petrophysical properties are offset.

##### 4.2. Vertical Displacement

In Figure 9, profiles of vertical displacements of caprock for three cases are compared. Continuous injection of CO_{2} increases fluid pressure in formation. Once the fluid pressure exceeds the least principal effective stress, the mechanical failure of caprock occurs, and CO_{2} leaks out of the storage formation. The vertical deformation of the caprock was calculated at the top of the caprock after 2 years of injection, just before the caprock cracks. Again, each curve indicates the average of 25 models. During CO_{2} injection into the formation, increased stress makes the surrounding formation rocks expand and push the caprock upwards. The caprock is bended more in the region just above the injector. As seen in Figure 9(a), for the case with heterogeneous petrophysical properties, all of the curves show higher vertical displacements than homogeneous case. Especially, values of vertical displacements are in order of . This result suggests that top of the caprock of homogeneous model collapses earlier than other models because of pressure concentration. Similarly, Figure 9(b) shows higher vertical displacement of geomechanically heterogeneous model than homogeneous model regardless of . Using petrophysically heterogeneous and geomechanically heterogeneous properties, Figure 9(c) depicts the combined effect of two properties.

**(a)**

**(b)**

**(c)**

Figure 10 reveals that the distinction of vertical displacement depends on at the center of the caprock. Heterogeneous models at Figure 10(a) show 4.4%, 5.1%, and 6.3% higher vertical displacement compared with homogeneous case, respectively. Also, Case B shows 4.3%, 3.9%, and 4.1% higher displacement, and Case C shows 4.5%, 5.5%, and 6.2% higher displacement as could be seen in Figures 10(b) and 10(c). In comparison among the effect of heterogeneity for most heterogeneous models ( = 0.7), Cases A and C show similar range of uncertainty, and Case B shows less deformation but it is still larger than the homogeneous model. As could be seen in Figure 10(a), the vertical displacement of caprock increases along , and the range of uncertainty is also increasing. Figure 10(b) also indicates that geomechanically heterogeneous caprock shows higher value of vertical displacement than homogeneous model, but uncertainty range and vertical displacement of heterogeneous caprock are relatively constant compared with Case A. It is well known that geomechanical properties have influence on rock deformation. Khajeh et al. [11] indicated that heterogeneous geomechanical properties result in significant differences in range of uncertainty in geomechanical response. In this study, geomechanical properties are correlated with petrophysical properties as a function. In this limited situation, heterogeneity in rock mechanical properties has less effect on the geomechanical change of formation in comparison to heterogeneity consideration of petrophysical properties.

**(a)**

**(b)**

**(c)**

The amounts of sequestered CO_{2} and range of uncertainty in vertical displacement are larger for models with higher Dykstra–Parsons coefficient. It indicates that heterogeneity for petrophysical properties result in unfavorable effects in flow output variables. In contrast, the heterogeneity of rock mechanical properties can reduce additional escape of stored CO_{2}. In terms of geomechanical change, heterogeneities of both petrophysical and geomechanical properties can affect the uplift of caprock. Highly heterogeneous caprock will be transformed more. Despite the heterogeneity of geomechanics affects the displacement of caprock, Case B shows relatively small range of uncertainty. It is obvious that geomechanical properties affect mechanical change of formation. However, there is limitation of displacement when reservoir properties are correlated. Since both geomechanical and flow responses are of importance to predict and optimize CO_{2} injection performance, ignoring heterogeneity effects for storage site properties may result in inaccurate analyses.

##### 4.3. Effects of Horizontal Tectonic Stress

In order to assess the effect of tectonic stresses, horizontal stress is added to highly heterogeneous model (c) ( = 0.7) and homogeneous model along the *x*-axis. The simulation model has initial stress of 15,503.6 kPa with hydrostatic gradient. To add a 5,000 kPa at *x*-axis, horizontal compression is loaded to whole grids. Figure 11 indicates the sequestered CO_{2}. Compared to cases considering normal stress only, the model with tectonic stress shows 24% more CO_{2} storage capacity just before the caprock failure during 4-5 years of injection. Due to the effect of additional horizontal compressive stress, the propagation of fractures is mitigated and a larger amount of CO_{2} is trapped. After the leakage occurs, the storage efficiency rapidly decreases for a homogeneous model with tectonic stress. Comparing two models including tectonic stress, the difference is up to 8% after 10 years of injection, and it is possible to get bigger. Vertical deformation at the center of the caprock is also increased by 0.6% at heterogeneous model and 1.1% at homogeneous model, respectively (Figure 12). Thus, sealing capacity could be better in more heterogeneous rocks under the high horizontal tectonic stress condition, but unpredicted high stress could make early fracture reactivation of caprock and bring unfavorable effect on long-term CO_{2} storage relatively homogeneous formation.

#### 5. Conclusion

Comparisons among storage sites with heterogeneous caprock showed the significance of coupled simulation of fluid flow and geomechanical deformation with respect to CO_{2} storage integrity. In particular, they have shown the importance of correlation between hydraulic and geomechanical properties in controlling CO_{2} leakage through caprock. Four major factors affecting mechanical caprock failure considered in this study are permeability, porosity, Young’s modulus, and Poisson’s ratio. These four properties are correlated as linear relationship based on field data. Heterogeneous petrophysical properties have unfavorable effects on CO_{2} storage, and geomechanical heterogeneity could mitigate further escape of CO_{2} after leakage occurs. Heterogeneity effects of these two properties can be offset for CO_{2} trapping. Heterogeneity of these properties also affects on vertical displacement of caprock. Results of analysis show that heterogeneity of four properties induces further displacement of caprock, but geomechanical properties have limited effects when it is correlated with other properties as a function. Subsurface storage needs to be evaluated in terms of heterogeneity to avoid any risk when operating CO_{2} sequestration. The hydrogeomechanical correlations are applied for more accurate investigation of CO_{2} leakage. We considered only four properties because correlations between hydraulic properties and other geomechanical properties have not been thoroughly understood. However, investigation including other formation properties such as rock compressibility and yield stress could provide better understanding of geomechanical study. Another improvement on the analysis of caprock stability can be achieved by the evaluation of 3D model. It could be an effective way to evaluate more realistic CO_{2} flow regime.

#### Conflicts of Interest

The authors declare that they have no conflicts of interest.

#### Acknowledgments

This work was supported by the Energy Efficiency & Resources of the Korea Institute of Energy Technology Evaluation and Planning (KETEP) grant funded by the Korea Government Ministry of Knowledge Economy (no. 20152520100760).