Subsurface Fluid Injection and Energy Storage
Open Access
Numéro
Oil Gas Sci. Technol. – Rev. IFP Energies nouvelles
Volume 74, 2019
Subsurface Fluid Injection and Energy Storage
Numéro d'article 81
Nombre de pages 14
DOI https://doi.org/10.2516/ogst/2019047
Publié en ligne 18 novembre 2019

© L. Li & J. Ma, published by IFP Energies nouvelles, 2019

Licence Creative Commons
This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1 Introduction

The 4D seismic method is a remarkable approach with which to monitor the distribution of CO2 plume in the reservoir during CO2-EOR or aquifer storage (Ivanova et al., 2012; White, 2009, 2013). The map of CO2 plume from 4D seismic monitoring data is used to evaluate the dead oil zone and the efficiency of CO2 flooding. The estimation of CO2 volume in the reservoir is our ultimate goal and one of the most powerful evidence to prove the safety of CO2 geological storage.

Conventional 4D seismic interpretation methods (Avseth et al., 2006) assume that porosity never be changed before and after fluids (CO2 and water) injection. The bulk modulus of the mineral matrix and bulk modulus of the porous rock frame will not be changed as well. Then fluids substitution methods based on Gassmann’s equation is used to make rock physics model (Mavko et al., 1998), and further to make 4D seismic forward models for interpretation.

However, more scientific evidence proved that reservoir pore system or pore structure will be changed after CO2 injection. Berrezueta et al. (2013) observed that intergranular clay matrix detachment and partial removal in CO2 flooded sandstone samples that lead to porosity increased 3.75%. Geochemical data observed in CO2-injected fields (Hovorka, 2009; Karamalidis et al., 2010; Kharaka et al., 2006; Thordsen et al., 2012) indicate that injection lowers the pH of the brine by dissolving a few tens of milligrams per litre of calcium (Ca), iron (Fe), and manganese (Mn) from the matrix as the CO2 front moves through the rock-water system; this phenomenon also indicates that dissolution occurs not only in the injection well but also in the reservoir. Tambach et al. (2015) simulated CO2 dissolution in the K12-B area in Holland and found that CO2 dissolution can increase the porosity and decrease the formation pressure. Vanorio (2015) performed time-lapse petrophysical experiments to observe the effects of fluid–rock interactions on petrophysical parameters and observed that the P- and S-wave velocities of rocks decreased before and after CO2 injection under the same pressure regardless of whether the injection was performed in a carbonate or a sandstone; this phenomenon is a result of the coupling between compaction and variations in the porosity. Fӧster et al. (2010) observed that approximately 1–5% of the porosity in the Ketzin CO2 sequestration area in Germany resulted from CO2 dissolution of clastic rock. Lu et al. (2016) experimentally examined the dissolution effect of CO2 on deposits and found that CO2 dissolution could increase the porosity and that the magnitude of the increase was associated with the pore throat radius. Several authors also found the porosity decreased after CO2 flooding (Skorpa et al., 2017).

Establishing a reasonable 4D seismic forward model is the key to correctly interpreting 4D seismic data, and 4D seismic responses result from variations in the fluid saturation or pressure in a deposit. Accordingly, the goals of using the 4D seismic method in this manner are to identify the effects of the CO2 saturation and pressure, to monitor the state of underground CO2, and to provide a theoretical understanding of the safety of geological CO2 sequestration.

The formation pressure and fluid saturation can vary after CO2 injection; in addition, the porosity also changes due to the fracturing either induced by high-pressure CO2 injection or due to CO2 dissolution of the rock. Thus, it is necessary to establish a 4D seismic forward model that considers variations in the fluid saturation, porosity, and formation pressure simultaneously, and accurate estimations of P- and S-wave velocities constitute the basis of establishing a 4D seismic forward model.

For present-day 4D seismic exploration, many P- and S-wave velocities estimation models are based on the Hertz-Mindlin model (Duffaut and Landrø, 2007; Mindlin, 1949; Saul and Lumley, 2015). This method considers the effects of both pressure and fluid. Duffaut and Landrø (2007) proposed the use of the Hertz-Mindlin expression to calculate P- and S-wave velocities and establish 4D seismic forward and inverse models in poorly consolidated sandstone. Saul and Lumley (2015) developed an inelastic model based on diagenesis and grain cement using the Hertz-Mindlin expression to describe the pressure sensitivity of rock.

However, none of these methods fully considers the effects of fluids. Instead, we propose that the effects of fluids should be considered in two aspects: (1) the effects of variations in the fluid composition, which can be achieved by the substitution of different fluids in the Gassmann equation, and (2) the effects of variations in the porosity on the seismic response after CO2 injection, which causes changes in the fluid. Ignoring variations in the porosity could result in misinterpretation; for example, the main factor influencing the seismic response could be the fluid composition, but neglecting porosity variations could result in the conclusion that pressure is the main factor.

The Digby equation assumes that cemented rock grains are present, whereas the Hertz-Mindlin expression is applicable to uncompacted sandstone reservoirs. The strata of interest in our study area, the Weyburn oil reservoir, can be divided into a marly layer and a vuggy layer; the former is composed of granular limestone, while the latter is finely crystalline to granular dolomite. Based on a petrological analysis and Scanning Electron Microscopy (SEM) inspection (Jensen, 2016; Njiekak et al., 2013), the Digby equation is suitable in the study area. Here, we improve upon the Digby equation to calculate the bulk modulus and shear modulus of dry rock and to estimate the P- and S-wave velocities of saturated rock under different pressures in combination with the Gassmann equation. The P- and S-wave velocities both pre- and post-CO2 injection are the basis of a 4D seismic forward model. The S-wave velocity is sensitive to the pressure rather than the fluid composition, whereas the P-wave velocity is affected by both the pressure and the fluid composition. Therefore, we focus mainly on variations in the P-wave velocity in this study. A 4D seismic forward model is established based on this relationship to investigate the effects of changes in the fluid composition and pressure on the Amplitude Versus Offset (AVO) gradient and intercept. In addition, the CO2 conditions studied here can be applied not only in the Weyburn deposit but also in other areas as well, and our motivation is to provide new 4D seismic interpretations.

2 Geological and geophysical background

The Weyburn deposit was discovered in 1954, and water was injected in 1965; more recently, CO2-Enhanced Oil Recovery (EOR) was initiated in 2000. The deposit is located in the Midale layer at depths of 1300–1500 m. The Midale layer can be divided into two layers: an upper marly dolomite (2–12 m) and a lower vuggy limestone (8–22 m). The marly layer has a high porosity (29%) and a low permeability (mean (M) = ~10 mD), whereas the vuggy layer has a lower porosity (10%) but a higher permeability (M = ~50 mD) (Brown, 2002). When injecting CO2, the injection pressure at the bottom of the well is 22 MPa, the original pore pressure is 15 MPa, and the pressure at the bottom of the production well is controlled to be 7–8 MPa (Ma and Morozov, 2010; White, 2013). According to research on the Weyburn deposit conducted by the British Geological Survey (Riding and Rochelle, 2005), the high-pressure injection of CO2 can result in fracturing and increase the total porosity by 0.5–1.9%.

Well 4-23 was drilled at the edge of Phase 1A in August 2000, approximately one month before CO2 injection (Fig. 1), to obtain well log curves (Fig. 2). The S-wave velocity data originate only from this well, which is therefore treated as the basis of the AVO model and deposit estimates. However, cores were collected from Wells 14-11 and 14-23 rather than from Well 4-23. Ultrasonic measurements were performed on four core samples from the marly and vuggy units (Brown, 2002). Their density, porosity, and permeability values are shown in Table 1. Our study is also based on these results.

thumbnail Fig. 1

The Phase 1A area of CO2 sequestration in the Weyburn field, Canada. Well 4-23 is the logging well, whereas Well 12-13 and Well 14-11 are the coring wells.

thumbnail Fig. 2

Well log data from Well 4-23. From left to right are the P-wave velocity (V p ), S-wave velocity (V s ), resistivity (Rt) log, Gamma Ray (GR) log, shale content (V sh), porosity (Por; the dashed line indicates the effective porosity, while the solid line indicates the total porosity), and density.

Table 1

Properties of the rock cores.

3 Theory

Digby (1981) proposed a method to calculate the bulk modulus and shear modulus of closely compacted dry rock based on changes in the pressure:

(1)

(2)where b can be expressed as:

(3)while d satisfies equation (4):

(4)where K dry and μ dry are the bulk and shear moduli, respectively, of dry rock; υ and μ ma are Poisson’s ratio and the shear modulus, respectively, of the solid grains; Φ is the porosity; C p is the coordination number; p is the differential pressure; a is the radius of the small, flat, circular regions of neighbouring bonded particles before deformation; b is the radius after deformation; and R is the grain radius. K dry, μ dry and μ ma have the same unit of measurement (GPa), while the unit of measurement for p is MPa, and υ, φ, C p , a/R and b/R are all dimensionless. The bulk and shear moduli of the rock matrix are calculated by using the Voigt-Reuss-Hill (V-R-H) method.

Murphy (1982) assumed C p to be a constant, summarized C p for different porosities, and concluded that C p is proportional to e 1−ϕ . The constant C p values for different porosities summarized by Murphy are shown in Table 2 (where the porosity and C p are dimensionless).

Table 2

Summary of information on various porosities and their corresponding C p values.

Table 2 shows the C p values corresponding to the porosities summarized by Murphy, where the porosity is listed discretely as 0.25, 0.3, 0.35,…, 0.7. In Figure 3, which displays a fitting diagram based on Murphy’s statement that C p is proportional to e 1−ϕ , we compare C p with e 1−ϕ and use the values obtained from linear fitting for subsequent calculations using the equation C p  = 11.759e 1−ϕ −12.748.

thumbnail Fig. 3

Linear fitting between C p and e 1−φ . Black dots are Murphy’s summarization. Black line is result of linear regression.

After fitting data with this formula, the corresponding C p values for a continuous change in the porosity can be calculated. However, C p varies with pressure (Wang et al., 2018). Based on this linear relationship, we propose a revised equation (5) for C p :

(5)where W is a weighting factor that represents several effects, such as the effects of the pressure and lithology. The P-wave velocity is usually measured in actual oil fields, and the best W value is obtained by subtracting the absolute value of the P-wave velocity in the actual logging data from the predicted P-wave velocity and establishing a minimum or zero difference (Eq. (6)):

(6)where is a measured P-wave velocity and (W) is a predicted value; the units of both and (W) are km/s. Substituting equation (5) into equations (1) and (2), we can obtain K dry(W) and μ dry(W), which can then be substituted into Gassmann’s equation, following which we can obtain the bulk and shear moduli of a saturated rock that contains W. Then, (W) can be represented as equation (7), where K sat(W) and μ sat(W) are the bulk and shear moduli, respectively, of a saturated rock that contains W; the units of both K sat(W) and μ sat(W) are GPa, and ρ is the density expressed in units of g/cm3. Equation (8) is an expansion of equations (6) and (7). Finally, we can calculate the shear modulus μ dry by substituting the predicted W and C p values into equation (2) to predict the S-wave velocity:

(7)

(8)

4 Calculation of elastic wave velocities

4.1 Method validation

As early as 2002, Brown performed the petrophysical tests on the rock of Weyburn field’s core samples in Colorado School of Mines. Four core samples were acquired from Marly and Vuggy units in two wells. For these four core samples, Brown measured the wave velocity once when pressure increased from low to high until the maximum pressure was reached. Then the velocity was measured again when pressure decreased from high to low. Two groups of wave velocity under two different pressures were obtained (Brown, 2002). For the compressional wave velocity, the average of the two measurements was selected. For the shear-wave velocity, because two values of V s1 and V s2 were measured for each core sample and each value of V s1 or V s2 was measured two times, there were four shear-wave velocities for each core sample. We took their average as the shear wave velocity in an isotropic medium and compared it with the predicted values (Brown, 2002).

Variations in the P-wave velocity during geological CO2 sequestration result from changes in both the fluid composition and the fluid pressure. To identify the influences of these two effects, the marly layer, which has the higher porosity between the two units, is selected as the subject of this study. Two marly rock cores, PC12-13 1410 and PC14-11 1420.2, are used to verify our method. Then, calculation models are developed to estimate the elastic wave velocities based on the rock testing results of rock core PC12-13 1410.

First, the S-wave velocity is estimated using the Digby equation in combination with the Gassmann expression based on the measured P-wave velocity, as shown in Tables 3 and 4 and Figures 4 and 5. The values of the weight coefficient W in Tables 3 and 4 are calculated using equation (8). The absolute errors of the estimated S-wave velocities are less than 5%, resulting in a prediction curve with a high accuracy, thereby indicating the efficiency of this method.

thumbnail Fig. 4

Comparison between the estimated and measured S-wave velocities in PC12-13 1410 from the marly layer. The error bar shows a 95% confidence interval for the mean of measured S-wave velocity.

thumbnail Fig. 5

Comparison between the estimated and measured S-wave velocities in PC14-11 1420.2 from the marly layer. The error bar shows a 95% confidence interval for the mean of measured S-wave velocity.

Table 3

Results for the marly layer in PC12-13 1410.

Table 4

Results for the marly layer in PC14-11 1420.2.

4.2 P- and S-wave velocity calculations

After verifying the accuracy of the proposed method, the physical parameters of rock core PC12-13 1410 are used for the following analysis because its average porosity (29%) is similar to that of the marly layer.

The P-wave velocity is calculated using the Digby equation with variable porosity and pressure. The Digby equation includes the pressure effect, and both the Digby equation and the coordination method proposed in our study include the porosity effect, which constitutes the theoretical basis of our study. Notably, based on the values of the weight coefficient shown in Tables 3 and 4, the coefficient remains constant under similar pressures with varied porosities (the porosities of these two rock cores are 0.29 and 0.33), implying that the weight coefficient is associated only with the pressure when the lithology is the same. On the other hand, the coordination number is associated with both the pressure and the porosity, enabling the S- and P-wave velocities of rock core PC12-13 1410 to be calculated under different porosities.

The elastic parameters employed herein are from Ma and Morozov (2010). To infer the effects of variations in the fluid on the P-wave velocity, P-wave velocities are calculated under different pressures, different porosities, and different CO2 saturations, as shown in Figure 6, the results of which reflect the variation in the P-wave velocity with a changing porosity. Different shapes represent different pore pressures, and different colours represent different CO2 saturations. The results show that P-wave velocity at various pressures can be distinguished under the different porosity and CO2 saturation. However, for a given pressure with variable porosity, the P-wave velocities under different CO2 saturations can overlap. For example, when the pore pressure is 8 MPa and the porosity is less than 0.24, the P-wave velocity under 5% CO2 saturation is greater than that under 20% CO2 saturation; however, when the porosity is greater than 0.24, the velocity pattern is reversed. This relationship implies that the effect of the porosity must be considered when examining variations in the fluid saturation.

thumbnail Fig. 6

Variations in the P-wave velocity with different porosities under different pressure and saturation conditions. The x-axis is the porosity, and the y-axis is the P-wave velocity. The triangle symbols represent a pore pressure of 8 MPa, and the star symbols represent a pore pressure of 15 MPa. The black, green, blue, and red colours represent the P-wave velocities associated with CO2 saturations of 5%, 10%, 15%, and 20%, respectively.

We focus on the effects of porosity variations on the S- and P-wave velocities and on the forward models before and after CO2 injection, following which P-wave velocities under different porosities and saturations are calculated for two fixed pressures. As shown in Figure 7, as the porosity varies, the P-wave velocities under different pressures cannot be identified because they overlap one another, that is, as the porosity increases by 1%, the P-wave velocity for a pore pressure of 8 MPa and a CO2 saturation of 10% (red triangle) overlaps the P-wave velocity for a pore pressure of 15 MPa and a CO2 saturation of 5% with the original porosity (black star). For a pore pressure of 8 MPa with increasing porosities and CO2 saturations, the P-wave velocities cannot be distinguished from those at a pore pressure of 15 MPa, highlighting the importance of considering porosity variations in the forward model. In Figure 6, the main factor influencing the P-wave velocity is still the fluid pressure; however, with varying porosity (Fig. 7), the P-wave velocities of two pressures cannot be distinguished from one another.

thumbnail Fig. 7

Variations in the P-wave velocity under different porosity, pressure and saturation conditions. The x-axis is the variation of porosity, and the y-axis is the P-wave velocity. The triangle symbols represent a pore pressure of 8 MPa, and the star symbols represent a pore pressure of 15 MPa. The black colour represents the P-wave velocity trend at a CO2 saturation of 5% and a porosity of 0.18–0.29, while the red, blue, green, grey and purple colours represent the P-wave velocity trends at CO2 saturations of 10%, 15%, 20%, 25%, and 30% and porosities that are 1–5%, respectively, higher than the data represented by the black colour.

Figure 8 shows the P-wave velocities at pore pressures of 8 MPa and 15 MPa. Similar to Figure 7, when the porosity is increased by 1% (red triangles), the P-wave velocity for a pore pressure of 8 MPa overlaps that for a pore pressure of 15 MPa.

thumbnail Fig. 8

Variations in the P-wave velocity under different pressures. The x-axis is the variation of porosity, and the y-axis is the P-wave velocity. The triangle symbols represent P-wave velocities at a pore pressure of 8 MPa and a CO2 saturation of 20%. In addition, the black colour represents the P-wave velocity under the original porosity (0.18–0.29), while the red, blue, green, grey and purple colours represent the P-wave velocities under porosities that are 1–5%, respectively, higher than the original porosity. The star symbols represent the P-wave velocities for a pore pressure of 15 MPa and a constant porosity, while the black, red, blue, and green colours represent the P-wave velocities under CO2 saturations of 5%, 10%, 15%, and 20%, respectively.

The P-wave velocity varies with the porosity. However, the effects of porosity variations on the S-wave velocity remain to be studied. Here, we study these effects with different pressures and CO2 saturations. In Figure 9, the x-axis represents the porosity, and the y-axis is the S-wave velocity; the figure shows data for different pressures and saturations and porosities. Clearly, the variation in the S-wave velocity is determined by the fluid pressure.

thumbnail Fig. 9

Variations in the S-wave velocity with the porosity. Different symbols represent different pore pressures, whereas different colours with the same symbol represent different CO2 saturations.

Here, we study the changes in the S-wave velocity with changes in the porosity. In Figure 10, the variation in the S-wave velocity is determined mainly by the pressure, although S-wave velocities at different pressures can still overlap in response to large variations in the porosity at high CO2 saturations (purple triangles). These results also demonstrate the significance of considering both the CO2 saturation and the porosity in the forward model.

thumbnail Fig. 10

Variations in the S-wave velocity with the porosity. The x-axis is the sample number, and the y-axis is the S-wave velocity. The triangle symbols represent the S-wave velocity at a pore pressure of 8 MPa and the original porosity (0.18–0.29). The black colour indicates a CO2 saturation of 5%, while the red, blue, green, grey, and purple colours represent CO2 saturations of 10%, 15%, 20%, 25%, and 30% and porosities that are 1–5%, respectively, higher than the original value. The star symbols represent a pore pressure of 15 MPa. The green, blue, red and black colours represent S-wave velocities with CO2 saturations of 5%, 10%, 15%, and 20%, respectively, for similar porosities.

5 Establishing a 4D seismic forward model

Based on the S- and P-wave velocities obtained under different pressures, CO2 saturations and porosities, the reflection coefficients are calculated using the Zoeppritz equation to compute the synthetic seismogram through the angle domain from convolution of the Ricker wavelet.

5.1 Two-layer model

A two-layer model is first developed for synthetic seismograms at 8 MPa with different saturations and porosities, and difference profiles are obtained by comparing these synthetic seismograms (Fig. 11) with the parameters shown in Table 5. As shown in Figure 11, the differences in the synthetic seismograms increase under similar pressures and saturations but with increasing porosity.

thumbnail Fig. 11

Synthetic seismograms of the two-layer model at a pore pressure of 8 MPa. The left-hand panels from left to right are the P-wave velocity, S-wave velocity, density, and synthetic seismogram. The right-hand panels from left to right are the P-wave velocity, S-wave velocity, density, synthetic seismogram and difference profile. In the plots of the P-wave velocity, S-wave velocity, and density, the curves are measured well log data, and the straight lines are block parameters. In each row, the synthetic seismogram on the far right is the difference profile between the synthetic seismograms on the left- and right-hand sides. a) The porosity is 0.18, while the CO2 saturation is 5% on the left and 20% on the right; b) The porosity is 0.24, while the CO2 saturation is 5% on the left and 20% on the right; c) The porosity is 0.34, while the CO2 saturation is 5% on the left and 20% on the right.

Table 5

Parameters used in the two-layer model shown in Figure 11.

The AVO gradient intercept figure (Fig. 12) corresponds to the two-layer model obtained above. As shown in Figure 12, variation in the porosity leads to variation in the AVO intercept, and variation in the CO2 saturation results in a change in the AVO gradient, which is consistent with the responses in the two-layer seismic model.

thumbnail Fig. 12

AVO gradient intercept plot. The stars, triangles, and circles represent porosities of 0.18, 0.24, and 0.34, respectively. The red and blue colours represent CO2 saturations of 5% and 20%, respectively.

To investigate the effects of porosity variations on the seismic response and AVO, synthesize seismograms for two-layer models with different porosities and pressures are established (Fig. 13) with the parameters shown in Table 6. In Figure 13a, the AVO trend of the difference between the synthetic seismograms increases in amplitude with an increase in the angle, whereas in Figure 13d, such a change is not obvious, implying the importance of considering variations in the porosity.

thumbnail Fig. 13

Two-layer models. The figures on the left show the P-wave velocity, S-wave velocity, density, and synthetic seismogram from left to right, whereas the figures on the right show the P-wave velocity, S-wave velocity, density, synthetic seismogram and difference between the synthetic seismograms. In the plots of the P-wave velocity, S-wave velocity and density, the curves are measured well log data, and the straight lines are block parameters. In each row, the difference between the synthetic seismograms on the left and right is shown on the far right. The conditions of the plots on the left and right of each row are as follows: a) Left: 15 MPa pore pressure, 0.2 porosity, and 5% CO2 saturation; right: 8 MPa pore pressure, 0.2 porosity, and 20% CO2 saturation; b) Left: 15 MPa pore pressure, 0.2 porosity, and 5% CO2 saturation; right: 8 MPa pore pressure, 0.21 porosity, and 20% CO2 saturation; c) Left: 15 MPa pore pressure, 0.2 porosity, and 5% CO2 saturation; right: 8 MPa pore pressure, 0.22 porosity, and 20% CO2 saturation; d) Left: 15 MPa pore pressure, 0.2 porosity, and 5% CO2 saturation; right: 8 MPa pore pressure, 0.23 porosity, and 20% CO2 saturation.

Table 6

Parameters used in the two-layer models shown in Figure 13.

To directly observe the effects of porosity and pressure variations on the AVO, the AVO gradient intercept is calculated based on the two-layer model data in Table 6. Figure 14 shows that variation in the fluid saturation leads to variation in the AVO gradient, where the variation in the AVO intercept is due to the variations in the pressure and porosity.

thumbnail Fig. 14

AVO gradient intercept plot. Circles represent the conditions under a 15 MPa pore pressure, 5% CO2 saturation, and 0.2 porosity. Stars represent the conditions under a 15 MPa pore pressure and 20% CO2 saturation. Triangles represent the conditions under an 8 MPa pore pressure and 20% CO2 saturation. The red, blue, green, and orange colours represent porosities of 0.2, 0.21, 0.22, and 0.23, respectively.

5.2 Well model

Based on the abovementioned two-layer models, a well model can be developed using data from Well 4-23 to compute synthetic seismograms and the differences between those seismograms. The well log curve from Well 4-23 is shown in Figure 2. First, a well model is derived using the original well log data to create a synthetic seismogram; then, the fluid is substituted based on the proposed method to predict the P- and S-wave velocities after CO2 injection to fully develop the well model. Before fluid substitution, the pore pressure is 15 MPa, and the CO2 saturation is 0; after fluid substitution, the pore pressure is 8 MPa, and the CO2 saturation is 20%. To consider porosity variations, two synthetic seismograms are generated after fluid substitution: one without changes in the porosity and one in which the porosity is increased by 1.5%. Figures 15a–15c show the synthetic seismograms of the well model before and after fluid substitution. The differences between the synthetic seismograms with varying porosity are greater than those between the synthetic seismograms with a consistent porosity, further demonstrating the effect of a changing porosity on the seismic response during the establishment of the forward model. However, further discussion of the unclear difference in Figure 15c is required.

thumbnail Fig. 15

Synthetic seismogram of the well model. The figures on the left show the P-wave velocity, S-wave velocity, density, and synthetic seismogram from left to right, whereas the figures on the right show the P-wave velocity, S-wave velocity, density, synthetic seismogram and difference between the seismograms (the synthetic seismogram on the left is subtracted from that on the right). In each row, the left panel shows the original synthetic seismogram, and the right panel shows the synthetic seismogram with a consistent porosity after fluid substitution.

6 Discussion

6.1 The effects of fluid variations on the P-wave velocity

Figure 7 shows that, at a constant pressure, the P-wave velocity decreases and then increases with increasing CO2 saturation and porosity. We propose two explanations for this observed pattern. First, the porosity increases to a critical value, at which point the carrier of the seismic wave converts from minerals to fluid (Mavko et al., 1998). When the porosity is lower than the critical value, the minerals impose a dominant influence on the P-wave velocity; therefore, with increasing CO2 saturation, the P-wave velocity decreases. After reaching the critical porosity, the main influence is derived from the fluid, and the variation in the P-wave velocity is associated with the properties of this fluid. Second, the properties of the fluid mixture vary with pressure. The bulk moduli and densities of the fluid mixture at different pressures are shown in Figures 16a and 16b, respectively. After fluid substitution, the pore pressure of the deposit is 8 MPa, which represents the bubble point. In Figure 16, the bulk moduli and densities of the fluid mixture exhibit their greatest changes below 8 MPa in each state. Due to great variations in the fluid properties, the P-wave velocity can vary substantially at high porosity values.

thumbnail Fig. 16

Properties of fluid mixtures under different pressures. The x-axis is the pore pressure. a) The y-axis is the bulk strain modulus, which varies with pressure, of the fluid mixture; b) The y-axis is the density, which varies with pressure, of the fluid mixture.

6.2 Analysis of the synthetic seismogram of the well model

According to the synthetic seismogram of the well model presented in Section 5.2, a distinct difference is observed between the original seismic traces and the seismic traces after fluid substitution, implying that CO2 injection causes variations in the fluid saturation and pore pressure, which can be identified based on the seismic response. However, no significant difference exists between the synthetic seismograms with and without considering porosity variations (Fig. 15c), and the porosity effect shown in the synthetic seismogram difference profile from the well model is not as clear as that from the two-layer model; this discrepancy is likely due to two reasons. First, according to the well log curve in Figure 2, the marly layer has a higher porosity than the vuggy layer. Since the seismic responses associated with porosity variations are due to the fluid effect, the porosity below 0.1 cannot significantly influence the P- and S-wave velocities or the seismic response, even if the porosity varies by 1.5%. Second, the deposit layer is thin, the marly layer is even thinner. Therefore, any variation in the marly layer is hardly observable due to the sampling restrictions of the preceding depth-time conversion procedure when synthesizing the seismogram.

7 Conclusion

Based on this study, we can conclude the following:

  1. Through calculations of the P- and S-wave velocities of the marly layer in the Weyburn deposit under different conditions, variations in the porosity can result in variations in the P- and S-wave velocities and even cause the velocities to overlap at different pressures. Therefore, in addition to variations in the pressure and fluid saturation, changes in the porosity must be considered when estimating 4D seismic P- and S-wave velocities, especially when the porosity of the deposit is relatively high.

  2. According to two layer synthetic seismogram, when comparing cases with and without porosity variations before and after fluid substitution, porosity variations may lead to greater seismic responses and can even change the AVO trend of the synthetic seismogram.

  3. Variations in porosity can provide a new possible explanation for observed 4D seismic patterns and can reduce the diversity of interpretations.

A 4D seismic forward model has been established. Inversion of 4D seismic data and interpretation 4D seismic data will be done on the basis of this 4D seismic forward model, which provide basis for the safety of CO2 geological storage.

References

  • Avseth P., Mukerji T., Mavko G. (2006) Quantitative seismic interpretation: Applying rock physics tools to reduce interpretation risk, Cambridge University Press, Cambridge, UK. [Google Scholar]
  • Berrezueta E., González-Menéndeza L., Breitner D., Luquot L. (2013) Pore system changes during experimental CO2 injection into detritic rocks: Studies of potential storage rocks from some sedimentary basins of Spain, Int. J. Greenhouse Gas Control 17, 411–422. [CrossRef] [Google Scholar]
  • Brown L.T. (2002) Integration of rock physics and reservoir simulation for the interpretation of time-lapse seismic data at Weyburn Field, Saskatchewan, Master’s Thesis, Colorado School of Mines, Golden, CO. [Google Scholar]
  • Digby P.J. (1981) The Effective elastic moduli of porous granular rocks, J. Appl. Mech. 48, 803–808. [Google Scholar]
  • Duffaut K., Landrø M. (2007) Vp/Vs ratio versus differential stress and rock consolidation – A comparison between rock models and time-lapse AVO data, Geophysics 72, 5, 81–94. [CrossRef] [Google Scholar]
  • Fӧster A., Schӧner R., Fӧster H.J., Norden B., Blaschke A.-W., Luckert J., Beutler G., Gaupp R., Rhede D. (2010) Reservoir characterization of a CO2 storage aquifer: The upper Triassic Stuttgart formation in the northeast German basin, Mar. Petrol. Geol. 27, 10, 2156–2172. [CrossRef] [Google Scholar]
  • Hovorka S. (2009) Frio brine pilot: The first U.S. sequestration test, South-West Hydrol. 8, 26–31. [Google Scholar]
  • Ivanova A., Kashubin A., Juhojuntti N., Kummerow J., Henninges J., Juhlin Ch, Lüth S., Ivandic M. (2012) Monitoring and volumetric estimation of injected CO2 using 4D seismic, petrophysical data, core measurements and well logging: A case study at Ketzin, Germany, Geophys. Prospect. 60, 957–973, doi: 10.1111/j.1365-2478.2012.01045.x. [Google Scholar]
  • Jensen G.K.S. (2016) Weyburn oilfield core assessment investigating cores from pre and post CO2 injection: Determining the impact of CO2 on the reservoir, Int. J. Greenhouse Gas Control 54, Part 2, 490–498. [CrossRef] [Google Scholar]
  • Karamalidis A., Hakala J.A., Griffith C., Hedges S., Lu J. (2010) Laboratory investigation of CO2 -rock-brine interactions using natural sandstones and brine samples from the SECARB Tuscaloosa injection zone. Presented at Geological Society of America Denver Annual Meeting, 31 October – 3 November, Denver, CO, USA. [Google Scholar]
  • Kharaka Y.K., Cole D.R., Hovorka S.D., Gunter W.D., Knauss K.G., Freifeld B.M. (2006) Gas-water-rock interactions in Frio Formation following CO2 injection: Implications for the storage of greenhouse gases in sedimentary basins, Geology 34, 577–580, doi: 10.1130/G22357.1. [Google Scholar]
  • Lu J., Nicot J.P., Mickler P.J., Ribeiro L.H., Darvari R. (2016) Alteration of Bakken reservoir rock duing CO2-based fracturing – An autoclave reaction experiment, J. Unconv. Oil Gas Res. 14, 2016, 72–85. [CrossRef] [Google Scholar]
  • Ma J., Morozov I.B. (2010) AVO modeling of pressure-saturation effects in Weyburn CO2 sequestration, Leading Edge 29, 2, 178–183. [CrossRef] [Google Scholar]
  • Mavko G., Mukerji T., Dvorkin J.P. (1998) The rock physics handbook: 147–159, Cambridge University Press, Cambridge, UK. [Google Scholar]
  • Mindlin R.D. (1949) Compliance of elastic bodies in contact, J. Appl. Mech. 16, 259–268. [Google Scholar]
  • Murphy W. (1982) Effects of microstructure and pore fluids on the acoustic properties of granular sedimentary materials, PhD Thesis, Stanford University, Stanford, CA. [Google Scholar]
  • Njiekak G., Schmitt D.R., Yam H., Kofman R.F. (2013) CO2 rock physics as part of Weyburn-Midale geological storage project, Int. J. Greenhouse Gas Control 16, Suppl. 1, S118–S133. [CrossRef] [Google Scholar]
  • Riding J.B., Rochelle C. (2005) The IEA Weyburn CO2 monitoring and storage project-final report of the European research team, British Geological Survey, Nottingam. [Google Scholar]
  • Saul M., Lumley D. (2015) The combined effects of pressure and cementation on 4D seismic data, Geophysics 80, 2, 135–148. [CrossRef] [Google Scholar]
  • Skorpa R., Todorovic J., Torsæter M. (2017) Porosity changes in mud-affected rock and cement upon reaction with CO2, Energy Procedia 114, 5266–5274. [Google Scholar]
  • Tambach T.J., Koenen M., Wasch L.J., Bergen F. (2015) Geochemical evaluation of CO2 injection and containment in a depleted gas field, Int. J. Greenhouse Gas Control 32, 2015, 61–80. [CrossRef] [Google Scholar]
  • Thordsen J.J., Kharaka Y.K., Thomas R.B., Ambats G., Abedini A., Manning M.A., Lu J. (2012) Natural heterogeneity and evolving geochemistry of Lower Tuscaloosa Formation brine in response to continuing CO2 injection at Cranfield EOR site, Mississippi, USA, American Geophysical Union 2012 Fall Meeting, December 3–7, 2012, San Francisco, CA, United States. [Google Scholar]
  • Vanorio T. (2015) Recent advances in time-lapse, laboratory rock physics for the characterization and monitoring of fluid-rock interactions, Geophysics 80, 2, WA49–WA59. [CrossRef] [Google Scholar]
  • Wang Y., Han D.H., Ren J.L., Zhang Y., Zhao L. (2018) Microstructure effects on static and dynamic moduli for two sandstones, Meeting, SEG 2018 DPRP, 20–22 May, Beijing, China, pp. 21–24. [Google Scholar]
  • White D.J. (2009) Monitoring CO2 storage during EOR at the Weyburn-Midale Field, Leading Edge 28, 7, 838–842. [CrossRef] [Google Scholar]
  • White D.J. (2013) Seismic characterization and time-lapse imaging during seven years of CO2 flood in the Weyburn field, Saskatchewan, Canada, Int. J. Greenhouse Gas Control 16S, 78–94. [CrossRef] [Google Scholar]

All Tables

Table 1

Properties of the rock cores.

Table 2

Summary of information on various porosities and their corresponding C p values.

Table 3

Results for the marly layer in PC12-13 1410.

Table 4

Results for the marly layer in PC14-11 1420.2.

Table 5

Parameters used in the two-layer model shown in Figure 11.

Table 6

Parameters used in the two-layer models shown in Figure 13.

All Figures

thumbnail Fig. 1

The Phase 1A area of CO2 sequestration in the Weyburn field, Canada. Well 4-23 is the logging well, whereas Well 12-13 and Well 14-11 are the coring wells.

In the text
thumbnail Fig. 2

Well log data from Well 4-23. From left to right are the P-wave velocity (V p ), S-wave velocity (V s ), resistivity (Rt) log, Gamma Ray (GR) log, shale content (V sh), porosity (Por; the dashed line indicates the effective porosity, while the solid line indicates the total porosity), and density.

In the text
thumbnail Fig. 3

Linear fitting between C p and e 1−φ . Black dots are Murphy’s summarization. Black line is result of linear regression.

In the text
thumbnail Fig. 4

Comparison between the estimated and measured S-wave velocities in PC12-13 1410 from the marly layer. The error bar shows a 95% confidence interval for the mean of measured S-wave velocity.

In the text
thumbnail Fig. 5

Comparison between the estimated and measured S-wave velocities in PC14-11 1420.2 from the marly layer. The error bar shows a 95% confidence interval for the mean of measured S-wave velocity.

In the text
thumbnail Fig. 6

Variations in the P-wave velocity with different porosities under different pressure and saturation conditions. The x-axis is the porosity, and the y-axis is the P-wave velocity. The triangle symbols represent a pore pressure of 8 MPa, and the star symbols represent a pore pressure of 15 MPa. The black, green, blue, and red colours represent the P-wave velocities associated with CO2 saturations of 5%, 10%, 15%, and 20%, respectively.

In the text
thumbnail Fig. 7

Variations in the P-wave velocity under different porosity, pressure and saturation conditions. The x-axis is the variation of porosity, and the y-axis is the P-wave velocity. The triangle symbols represent a pore pressure of 8 MPa, and the star symbols represent a pore pressure of 15 MPa. The black colour represents the P-wave velocity trend at a CO2 saturation of 5% and a porosity of 0.18–0.29, while the red, blue, green, grey and purple colours represent the P-wave velocity trends at CO2 saturations of 10%, 15%, 20%, 25%, and 30% and porosities that are 1–5%, respectively, higher than the data represented by the black colour.

In the text
thumbnail Fig. 8

Variations in the P-wave velocity under different pressures. The x-axis is the variation of porosity, and the y-axis is the P-wave velocity. The triangle symbols represent P-wave velocities at a pore pressure of 8 MPa and a CO2 saturation of 20%. In addition, the black colour represents the P-wave velocity under the original porosity (0.18–0.29), while the red, blue, green, grey and purple colours represent the P-wave velocities under porosities that are 1–5%, respectively, higher than the original porosity. The star symbols represent the P-wave velocities for a pore pressure of 15 MPa and a constant porosity, while the black, red, blue, and green colours represent the P-wave velocities under CO2 saturations of 5%, 10%, 15%, and 20%, respectively.

In the text
thumbnail Fig. 9

Variations in the S-wave velocity with the porosity. Different symbols represent different pore pressures, whereas different colours with the same symbol represent different CO2 saturations.

In the text
thumbnail Fig. 10

Variations in the S-wave velocity with the porosity. The x-axis is the sample number, and the y-axis is the S-wave velocity. The triangle symbols represent the S-wave velocity at a pore pressure of 8 MPa and the original porosity (0.18–0.29). The black colour indicates a CO2 saturation of 5%, while the red, blue, green, grey, and purple colours represent CO2 saturations of 10%, 15%, 20%, 25%, and 30% and porosities that are 1–5%, respectively, higher than the original value. The star symbols represent a pore pressure of 15 MPa. The green, blue, red and black colours represent S-wave velocities with CO2 saturations of 5%, 10%, 15%, and 20%, respectively, for similar porosities.

In the text
thumbnail Fig. 11

Synthetic seismograms of the two-layer model at a pore pressure of 8 MPa. The left-hand panels from left to right are the P-wave velocity, S-wave velocity, density, and synthetic seismogram. The right-hand panels from left to right are the P-wave velocity, S-wave velocity, density, synthetic seismogram and difference profile. In the plots of the P-wave velocity, S-wave velocity, and density, the curves are measured well log data, and the straight lines are block parameters. In each row, the synthetic seismogram on the far right is the difference profile between the synthetic seismograms on the left- and right-hand sides. a) The porosity is 0.18, while the CO2 saturation is 5% on the left and 20% on the right; b) The porosity is 0.24, while the CO2 saturation is 5% on the left and 20% on the right; c) The porosity is 0.34, while the CO2 saturation is 5% on the left and 20% on the right.

In the text
thumbnail Fig. 12

AVO gradient intercept plot. The stars, triangles, and circles represent porosities of 0.18, 0.24, and 0.34, respectively. The red and blue colours represent CO2 saturations of 5% and 20%, respectively.

In the text
thumbnail Fig. 13

Two-layer models. The figures on the left show the P-wave velocity, S-wave velocity, density, and synthetic seismogram from left to right, whereas the figures on the right show the P-wave velocity, S-wave velocity, density, synthetic seismogram and difference between the synthetic seismograms. In the plots of the P-wave velocity, S-wave velocity and density, the curves are measured well log data, and the straight lines are block parameters. In each row, the difference between the synthetic seismograms on the left and right is shown on the far right. The conditions of the plots on the left and right of each row are as follows: a) Left: 15 MPa pore pressure, 0.2 porosity, and 5% CO2 saturation; right: 8 MPa pore pressure, 0.2 porosity, and 20% CO2 saturation; b) Left: 15 MPa pore pressure, 0.2 porosity, and 5% CO2 saturation; right: 8 MPa pore pressure, 0.21 porosity, and 20% CO2 saturation; c) Left: 15 MPa pore pressure, 0.2 porosity, and 5% CO2 saturation; right: 8 MPa pore pressure, 0.22 porosity, and 20% CO2 saturation; d) Left: 15 MPa pore pressure, 0.2 porosity, and 5% CO2 saturation; right: 8 MPa pore pressure, 0.23 porosity, and 20% CO2 saturation.

In the text
thumbnail Fig. 14

AVO gradient intercept plot. Circles represent the conditions under a 15 MPa pore pressure, 5% CO2 saturation, and 0.2 porosity. Stars represent the conditions under a 15 MPa pore pressure and 20% CO2 saturation. Triangles represent the conditions under an 8 MPa pore pressure and 20% CO2 saturation. The red, blue, green, and orange colours represent porosities of 0.2, 0.21, 0.22, and 0.23, respectively.

In the text
thumbnail Fig. 15

Synthetic seismogram of the well model. The figures on the left show the P-wave velocity, S-wave velocity, density, and synthetic seismogram from left to right, whereas the figures on the right show the P-wave velocity, S-wave velocity, density, synthetic seismogram and difference between the seismograms (the synthetic seismogram on the left is subtracted from that on the right). In each row, the left panel shows the original synthetic seismogram, and the right panel shows the synthetic seismogram with a consistent porosity after fluid substitution.

In the text
thumbnail Fig. 16

Properties of fluid mixtures under different pressures. The x-axis is the pore pressure. a) The y-axis is the bulk strain modulus, which varies with pressure, of the fluid mixture; b) The y-axis is the density, which varies with pressure, of the fluid mixture.

In the text

Les statistiques affichées correspondent au cumul d'une part des vues des résumés de l'article et d'autre part des vues et téléchargements de l'article plein-texte (PDF, Full-HTML, ePub... selon les formats disponibles) sur la platefome Vision4Press.

Les statistiques sont disponibles avec un délai de 48 à 96 heures et sont mises à jour quotidiennement en semaine.

Le chargement des statistiques peut être long.