Abstract

Shale is a typical medium of transverse isotropy with a vertical axis of symmetry (VTI), and its strong anisotropy is mainly due to the combined effect of intrinsic anisotropy and that induced by horizontal fractures. To calculate the anisotropy parameters of shale, a physical rock model is built based on Hudson's thin-coin fracture model and Schoenberg's linear-sliding model, and an approximate theoretical calculation method for Thomsen's anisotropy parameters of VTI media with horizontal fractures is proposed. These calculation results using the proposed method confirm that this anisotropy contributed to by horizontal fractures cannot be ignored to the overall anisotropy of shale. To simplify Rüger's formula that is an approximate theoretical formula for calculating the anisotropic reflection coefficients of VTI media, a new four-term approximate formula is derived in a standard reflectivity form based on Rüger's and Aki–Richards’ formulas. The simulation results of a VTI theoretical model and logging data of shale reservoirs show that there is only a small difference between the newly derived four-term formula and Rüger's formula for incidence angles <40°, and the new four-term formula can correctly reveal the seismic amplitude-versus-offset (AVO) characteristics of VTI media and fully retain the corresponding anisotropic seismic responses. Compared to Rüger's formula, the proposed new formula only has four terms of unknown parameters and can directly decouple Thomsen's anisotropy parameter ε from them, which helps to alleviate the ill-posed problems of simultaneous inversion of multiple parameters and enhance its application potential in seismic inversion of VTI media as shale.

1. Introduction

Anisotropy is one of the most significant geophysical features of shale reservoirs that can generally be approximated as a special medium of transverse isotropy with a vertical axis of symmetry (VTI), and the difference of seismic wave velocity measured by vertical bedding and parallel bedding can be as high as 60% (Wang 2002; Lonardelli et al.2007). The VTI anisotropy of shale is not only related to its inherent properties, and the directional arrangement of clay minerals is an important reason, as well as the organic matter and the pores with small aspect ratios in shale are also important factors affecting its inherent anisotropy strength (Kwon et al.2004; Sondergeld & Rai 2011; Deng et al.2014; Zhang et al.2017). The effects of porosity, pore structure and organic matter content are analyzed on elastic parameters including P- and S-wave velocities using a physical rock model (Deng et al.2015). The elastic anisotropy characteristics of shale and their influencing factors are studied in Wufeng Formation of Upper Ordovician (O3w) and Longmaxi Formation of Lower Silurian (S1l) using systematic physical experiments and theoretical analysis (Liu et al.2019). It is believed that the anisotropy of shale has a good correlation with its clay content, and the degree of directional arrangement of clay has an important impact on this relationship. The mineral composition, microstructure and micropore volume of Longmaxi Formation shale reservoirs at different depths in southeastern Sichuan are systematically studied using the thin slice method, X-ray diffraction, high-resolution argon ion polishing scanning electron microscopy (SEM), mercury intrusion, nitrogen adsorption, and other methods, combined with the analysis of total organic carbon content and porosity (Liu et al.2020). The influence of sedimentary and diagenetic processes is analyzed on the elastic parameters of the Wufeng-Longmaxi Formation shale (Zhao et al.2020). This research mainly highlights the inherent anisotropy of shale. However, the analysis of fracture-induced anisotropy is not in-depth, and the corresponding petrophysical model has not been completely established.

Shale is generally characterized by coupling of intrinsic anisotropy and fracture-induced anisotropy. Fracture development is very important for the storage and fracturing of shale gas. Especially, the horizontal bedding fracture is a structural weak plane during hydraulic fracturing, which is easy to form artificial network cracks (Engelder et al.2009). The different types of fracture cause also large variations in physical rock properties. For example, the horizontal bedding fractures affect the polarization anisotropy strength of shale (in VTI media) (Guo & Liu 2018). By analyzing this fracture-induced anisotropy, the equivalent static theory for media with isolated fractures is established (Hudson 1980). The linear-sliding model equates the fracture to a plane with zero thickness, and considers that the equivalent compliance tensor of fracture media is only related to two compliance parameters of fracture (Schoenberg 1980, 1983; Schoenberg & Sayers 1995). The relationship is showed between the anisotropy parameters in another special medium of transverse isotropy with a horizontal axis of symmetry (HTI) and the petrophysical parameters of the Schoenberg fracture model (Bakulin et al.2000). However, due to the lack of theoretical calculation formula of rock physics combining the intrinsic VTI anisotropy of shale and the induced anisotropy of horizontal fractures, it is difficult to further characterize the mineral composition of shale and the geological characteristics of horizontal fractures by using the anisotropy parameters from seismic inversion. Therefore, according to the relationship between the Hudson and Schoenberg fracture models, it is necessary to establish the petrophysical theoretical calculation formula of VTI anisotropy parameters of shale with horizontal fractures for the prediction of horizontal fracture density in the subsurface strata.

Anisotropic seismic inversion is an important approach to obtain anisotropy parameters and fracture density, and the azimuth AVO characteristics have been proved to be applicable to high-angle fracture prediction (Wang et al.2017; Shi et al.2018; Ma et al.2019; Huang et al.2020; Pan et al.2020; Li et al.2023). Rüger's approximate formula is the most widely known and used for reflection and transmission of VTI and HTI media (Rüger 1997, 1998; Rüger & Tsvankin 1997). However, it is a five-item approximation formula, resulting in many parameters to be inverted, large sensitivity differences between parameters, and highly ill-posed problems from seismic amplitude-versus-offset or angle (AVO/AVA) inversion for VTI media using PP-waves only (Zhang et al.2020; Li et al.2022). This can pose a great challenge for VTI inversion of anisotropy parameters and the corresponding prediction of horizontal fractures (low-angle fractures) in shale reservoirs. Therefore, it is necessary to further simplify Rüger's approximate formula and reduce the number of parameters to be inverted in the approximate formula for VTI media.

In this paper, the systematic tests of P- and S-wave velocity in two directions of shale reservoir cores in the O3w-S1lFormations in the Sichuan Basin are introduced to confirm that shale has the intrinsic anisotropy characteristics of VTI media. Then, in view of the widespread development of low-angle fractures in shale, a petrophysical approximate calculation method for VTI anisotropy parameters of shale with horizontal fractures is constructed, and a more simplified approximate formula for calculating the reflection coefficients of VTI media is derived to carry out the corresponding AVO theoretical analysis of shale models.

2. Shale anisotropy and its corresponding petrophysical calculation

2.1. Shale intrinsic anisotropy

The properties of an elastic rock media are described by the stiffness matrix C, which determines the relationship between stress and strain. Transverse isotropic media can be divided into VTI media and HTI media. Shales generally exhibit this anisotropy of VTI, and their stiffness matrix is as follows (Rüger 1997):

(1)

The five elastic stiffness coefficients are as follows:

(2)

where ρ is the rock bulk density, |${V}_{p0}$| and |${V}_{p90}$| are the vertical and horizontal P-wave velocities, respectively, and |${V}_{s0}$| and |${V}_{s90}$| are the vertical and horizontal S-wave velocities, respectively.

A simple expressions for the elastic constants is proposed using |$\rho $|⁠, |${V}_{p0}$|⁠, |${V}_{s0}$| and five other constants (Thomsen 1986):

(3)

where |$\varepsilon $| and |$\gamma $| represent Thomsen's anisotropy parameters of quasi P-waves and quasi S-waves, respectively, and Thomsen's anisotropy parameter |$\delta $| determines the complexity of the front shape of quasi P- and S-waves.

We carried out the systematic sample collection in the Wufeng-Longmaxi (O3w-S1l) shale of well Dys1 in the block DX of Sichuan Basin. The selected vertical and horizontal core samples were used for determining the P- and S-wave velocities experimentally. According to the actual drilling pressure data of shale reservoirs in the area, the maximum confining pressure was 60 MPa (shown as figure 1). Figure 1b–e show the experimental |${V}_{p0}$| and |${V}_{s0}$|⁠, and the logging velocities in the vertical direction (red and blue curves). The experimental velocities are is well consistent with the logging velocities. It can be seen that |${V}_{p90}$| is significantly larger than |${V}_{p0}$|⁠, and their difference tends to increase from the bottom to the top of this well (shown as figure 1b and c). Similarly, |${V}_{s90}$| is significantly higher than |${V}_{s0}$|⁠, and the difference tends to increase the bottom to the top of this well (shown as figure 1d and e). The anisotropy parameters |$\varepsilon $| and |$\gamma $| are calculated using equation (3) and shown as the scattered data points in figure 1f and g. The differences in the velocities in the two directions are as high as 15–60%, with a trend of increasing from the bottom to the top of this well.

(a) Experimental clay mineral content (scatter plot) and logging data (curve plot) interpretation clay mineral content, (b) experimental ${V}_{p0}$ and logging ${V}_{p0}$, (c) experimental ${V}_{p90}$, (d) experimental ${V}_{s0}$ and logging ${V}_{s0}$, (e) experimental ${V}_{s90}$, (f) Thomsen's anisotropy parameter $\varepsilon $, (g) Thomsen's anisotropy parameter $\gamma $, and (h) Thomsen's anisotropy parameter $\delta $.
Figure 1.

(a) Experimental clay mineral content (scatter plot) and logging data (curve plot) interpretation clay mineral content, (b) experimental |${V}_{p0}$| and logging |${V}_{p0}$|⁠, (c) experimental |${V}_{p90}$|⁠, (d) experimental |${V}_{s0}$| and logging |${V}_{s0}$|⁠, (e) experimental |${V}_{s90}$|⁠, (f) Thomsen's anisotropy parameter |$\varepsilon $|⁠, (g) Thomsen's anisotropy parameter |$\gamma $|⁠, and (h) Thomsen's anisotropy parameter |$\delta $|⁠.

Figure 2 reveals the clear relationships between the experimental |$\varepsilon $|⁠, |$\gamma $|and the experimental content of clay minerals |${V}_{clay}$|⁠, which can be fitted by the following polynomials:

(4)
(5)
(a) Thomsen's anisotropy parameter $\varepsilon $ vs. the content of clay minerals. (b) Thomsen's anisotropy parameter $\gamma $ vs. the content of clay minerals.
Figure 2.

(a) Thomsen's anisotropy parameter |$\varepsilon $| vs. the content of clay minerals. (b) Thomsen's anisotropy parameter |$\gamma $| vs. the content of clay minerals.

Thomsen's anisotropy parameters |$\varepsilon $| and |$\gamma $|(curves in figure 1f and g) can be calculated by using the logging data interpretation curve of clay mineral content (black curve in figure 1a). The calculated curve is very consistent with the experimental data.

For the same shale reservoir, there is an exponential relationship between the anisotropy parameters and the content of clay minerals (Deng et al.2015). The results presented in this paper are essentially consistent with those findings. The values of |$\delta $| can be obtained by using the positive correlation between |$\delta $| and |$\varepsilon $|⁠. On this basis, the stiffness matrix of the background media with intrinsic VTI characteristics can be obtained from equation (3). Owing to the clear positive correlations between the shale intrinsic anisotropy parameters and the content of clay minerals, the experimental anisotropy parameters |$\varepsilon $| and |$\gamma $| are mainly affected by mineral components. By analyzing their sensitivity to geophysical parameters, we found that they have a clear positive correlation with the experimental density ρ, as shown in figure 3. This demonstrates that we can predict the shale intrinsic anisotropy from the density from the inversion of seismic data.

(a) Thomsen's anisotropy parameter $\varepsilon $ vs. density. (b) Thomsen's anisotropy parameter $\gamma $ vs. density.
Figure 3.

(a) Thomsen's anisotropy parameter |$\varepsilon $| vs. density. (b) Thomsen's anisotropy parameter |$\gamma $| vs. density.

2.2. Anisotropy induced by horizontal fractures of shales

In addition to the strong VTI anisotropy caused by clay directivity, the actual shale strata often develop many low-angle or horizontal fractures with the mixed characteristics of intrinsic and fracture-induced anisotropy, which can often also be seen in core and imaging logging data (as shown in figure 4). Therefore, the influence of horizontal fractures on the equivalent anisotropy media needs to be considered.

Horizontal fractures in cores and image logging data.
Figure 4.

Horizontal fractures in cores and image logging data.

Based on the equivalent static theory for the media with isolated fractures and coin-shaped fractures, the following formulas for fracture density and porosity are proposed (Hudson 1980):

(6)

where e is the fracture density, |${\phi }_f$| is the fracture porosity, |${V}_f$| is the fracture volume, |${V}_a$| is the wellbore volume, which can be calculated from FMI image logging data, and |$\chi $| is the crack aspect ratio according to the Hudson coin fracture model.

The elastic matrix for fractured media is formulated by the linear slip model (Schoenberg 1980, 1983; Schoenberg & Sayers 1995). So, the compliance tensor of fracture reservoir is approximately equivalent to the sum of the compliance tensors of the isotropic background media and the disturbed fractures. The equivalence implies that it remains a strongly anisotropic media of VTI. Therefore, based on the Schoenberg's linear slip theory, the compliance matrix of the actual shale is approximately equivalent to the sum of the compliance matrices of the background media and the disturbed fractures (Schoenberg & Sayers 1995):

(7)

where |${S}_b$| is the compliance tensor of the isotropic background media, and |${S}_f$| is the compliance tensor of the disturbed fractures. Since horizontal fractures are present, the compliance tensor |${S}_f$| is as follows:

(8)

where

(8a)

and M is the modulus of P-waves, |$\lambda $| and |$\mu $| are the background Lame parameters and the background shear modulus of the shale inherent anisotropy media, respectively, |${Z}_N$| is the normal compliance of horizontal fractures, |${Z}_T$| is the tangential compliance of horizontal fractures, |${\Delta }_N$| is the normal weakness of horizontal fractures, and |${\Delta }_T$| is the tangential weakness of horizontal fractures.

We assume after the Schoenberg's linear slip theory that these fractures are filled with a weak media and their bulk and shear moduli are |$k^{\prime}$| and |$\mu ^{\prime}$|⁠, respectively. Then, the fracture normal weakness |${\Delta }_N$| and tangential weakness |${\Delta }_T$| can be calculated as follows (Hudson 1980; Bakulin et al.2000):

(9)

where |${g}_0 = {{{\mu }_0} {/ {\vphantom {{{\mu }_0} {( {{\lambda }_0 + 2{\mu }_0} )}}} } {( {{\lambda }_0 + 2{\mu }_0} )}}$|⁠, |${\mu }_0$| and |${\lambda }_0$| are the shear modulus and the Lame parameter of the background media without fractures, respectively.

Figure 5 shows the development of horizontal fractures or low-angle fractures during the actual drilling. The normal and tangential weakness of fractures are calculated by equation (9). It can be seen in figure 5 that in the fracture development section, the values of normal and tangential weakness can reach 0.1–0.24. Using the Schoenberg linear slip theory, the approximate compliance matrix of VTI media with horizontal fractures can be obtained as follows (Schoenberg & Sayers 1995):

(10)
Horizontal fractures observed in image logging data and calculated fracture normal weakness ${\Delta }_N$ and tangential weakness ${\Delta }_T$.
Figure 5.

Horizontal fractures observed in image logging data and calculated fracture normal weakness |${\Delta }_N$| and tangential weakness |${\Delta }_T$|⁠.

By inverting the equivalent compliance matrix S, the equivalent stiffness matrix |${C}_{VTI}$| of VTI shale media with horizontal fractures can be obtained. In view of the complexity of the analytical solution, only the numerical solutions are used in this study. Finally, the equivalent Thomsen's anisotropy parameters |$\tilde{\varepsilon }$| and |$\tilde{\delta }$| that induced by horizontal fractures can be calculated by using the stiffness matrix |${C}_{VTI}$| and equation (3).

In summary, the calculation process of Thomsen's anisotropy parameters for shale VTI media with horizontal fractures is as follows:

  • Calculating the intrinsic Thomsen's anisotropy parameters |$\varepsilon $|⁠, |$\gamma $|⁠, and |$\delta $|⁠. First, |$\varepsilon $| and |$\gamma $| the can be calculated based on the clay mineral content |${V}_{clay}$| of logging curves using equations (4) and (5), and then |$\delta $| can be estimated using the positive correlation relationship between |$\delta $| and |$\varepsilon $|(Deng et al.2015).

  • Calculating the normal weakness |${\Delta }_N$| and the tangential weakness |${\Delta }_T$| of low-angle fractures. First, the density e of low-angle fractures can be calculated by equation (6), and then |${\Delta }_N$| and |${\Delta }_T$| can be obtained by equation (9).

  • The equivalent compliance matrix S can be obtained by the sequential application of equations (7), (8), and (10).

  • Calculating the equivalent Thomsen's anisotropy parameters |$\tilde{\varepsilon }$|⁠, |$\tilde{\lambda }$|⁠, and |$\tilde{\delta }$|⁠. First, the inverse matrix |${C}_{VTI}$| of the compliance matrix S can be estimated by the numerical computation, and then |$\tilde{\varepsilon }$|⁠, |$\tilde{\lambda }$| and |$\tilde{\delta }$| can be calculated using the five elastic stiffness coefficients of |${C}_{VTI}$|⁠.

Figure 6 shows the calculation results. The pink curves in the figure are the calculated |$\varepsilon $| and |$\delta $| of VTI media without the effect of horizontal fractures, which are positively correlated with the content of clay minerals. The black curves in figure 6 are the calculated |$\tilde{\varepsilon }$| and |$\tilde{\delta }$| of VTI media with the effect of horizontal fractures. The maximum difference of the two calculated values is 0.1. This shows that the VTI anisotropy induced by horizontal fractures contributes significantly to the overall VTI anisotropy of shale and cannot be ignored, but should rather be considered jointly with the intrinsic anisotropy in the theoretical calculations of shale anisotropy parameters.

Comparison of calculated Thomsen's anisotropy parameters for VTI media without and with horizontal fractures.
Figure 6.

Comparison of calculated Thomsen's anisotropy parameters for VTI media without and with horizontal fractures.

3. Approximate formula of VTI reflection coefficient and its AVO analysis

3.1. Four-term approximate formula of VTI reflection coefficient

It is well known that the Rüger's approximate formula of reflection and transmission of VTI media derived by Rüger is widely used (Rüger 1997, 1998; Pšenčík & Martins 2001; Chen 2015). This formula takes the isotropic media as the background and introduces Thomsen's parameters to realize the linear approximation of the reflection and transmission coefficients of anisotropy media. Rüger's approximate formula of reflection coefficients can be divided into the isotropic background part and the anisotropic disturbance part (Rüger 1997, 1998):

(11)

where |$\theta $| represents the incident angle. The formulas for the isotropic background reflection term |${R}_{ISO}\! ( \theta )$| and the anisotropic disturbance reflection term |${R}_{ANI}( \theta )$| are given as equations (12) and (13), respectively:

(12)
(13)

where |$\Delta $| is the difference between the lower and upper medium parameters, |$\overline \bullet $| is the average of the parameters on both sides of the medium interface, and |$\varepsilon $| and |$\delta $| are Thomsen's anisotropy parameters of VTI media.

The isotropic background reflection coefficient |${R}_{ISO}( \theta )$| can be rewritten in different forms for some particular tasks. The petrophysical analysis has shown that the velocity and density terms are the key parameters that can characterize the anisotropy of shale in Section 2. Therefore, in this paper the |${R}_{ISO}( \theta )$| will be rewritten as the Aki–Richards approximate formula for |${V}_{p0}$|⁠, |${V}_{s0}$|⁠, and |$\rho $| (Aki & Richards 1980):

(14)

where |$g = {( {\frac{{\overline {{V}_{s0}} }}{{\overline {{V}_{p0}} }}} )}^2$| is the square of the ratio of |$\overline {{V}_{p0}} $| and |$\overline {{V}_{s0}} $| as a known quantity, |$\frac{{\Delta {V}_{p0}}}{{\overline {{V}_{p0}} }}$| is the reflectivity of |${V}_{p0}$|⁠, |$\frac{{\Delta {V}_{s0}}}{{\overline {{V}_{s0}} }}$| is the reflectivity of |${V}_{s0}$|⁠, and |$\frac{{\Delta \rho }}{{\overline \rho }}$| is the reflectivity of |$\rho $|⁠.

Using |${\sin }^2\theta {\tan }^2\theta = {\tan }^2\theta - {\sin }^2\theta $| and |${\sec }^2\theta = {\tan }^2\theta + 1$|⁠, after reorganizing the isotropy and anisotropy terms, a new approximate formula for reflection coefficients of VTI media is obtained as in equation (15):

(15)

Compared to equation (11), there are only four terms in equation (15). The reflection coefficient terms can be further approximated using equation (16):

(16)

Considering that |${V}_{p0}$| and |${e}^\varepsilon $| of adjacent strata have small difference and are differentiable, there are |$\left \{ { _{\scriptscriptstyle{e}^\varepsilon \approx \overline {{e}^\varepsilon }}^{\scriptscriptstyle\overline {{V}_{p0}} \approx {V}_{p0}}} \right.$|⁠, and then |$\overline {{V}_{p0}} {e}^\varepsilon \approx \overline {{V}_{p0}} \overline {{e}^\varepsilon } \approx \overline {{V}_{p0}{e}^\varepsilon } $|⁠. Equation (16) can then be rewritten as equation (17):

(17)

Assuming that g as a known coefficient has independence from elastic parameters, the following equation (18) can be derived by following the derivation method of equation (17):

(18)

Substituting these equations (17) and (18), equation (15) can be rewritten as follows:

(19)

It can be seen from equation (19) that it has the standard reflectivity form, and the VTI reflection coefficients are obtained by the weighted sum of the reflectivity of four parameters. Therefore, we can use the traditional workflow of anisotropy and isotropy AVO/AVA inversion to obtain |${V}_{p0}{e}^\varepsilon $|⁠, |${V}_{p0}$|⁠, and |${V}_{s0}{e}^{\frac{{\varepsilon - \delta }}{{8g}}}$|⁠, and then calculate the anisotropy parameters |$\varepsilon $| and |$\delta $| of VTI media. In theory, if only estimating |$\varepsilon $|⁠, it can be achieved by directly using the equation (19) to perform four-parameter AVA inversion.

3.2. AVO of VTI reflection coefficients of shale

To test the AVO/AVA accuracy of the new four-term approximation formula derived in this paper (equation 19), we compare it to the Rüger's approximate formula and the Aki–Richards’ approximate formula. We design a VTI model of shale reservoir sandwiched between roof strata of isotropic mudstone and floor strata of isotropic limestone. The parameters of this model are listed in Table 1. The calculation results of incident angles from 0° to 40° are shown in figure 7.

Comparisons between Rüger's approximate formula and the new four-term approximate formula on RVTI. (a) The top-interface reflection coefficients from the shale reservoir model. (b) The bottom-interface reflection coefficients from the shale reservoir model.
Figure 7.

Comparisons between Rüger's approximate formula and the new four-term approximate formula on RVTI. (a) The top-interface reflection coefficients from the shale reservoir model. (b) The bottom-interface reflection coefficients from the shale reservoir model.

Table 1.

Elasticity and anisotropy parameters of a VTI model for shale reservoir

VpVsρεδ
Medium parameters(m s−1)(m s−1)(g cm−3)εδ
Mudstone (top layer)440023002.65--
Shale reservoir400021002.550.300.20
Limestone (bottom layer)550035002.70--
VpVsρεδ
Medium parameters(m s−1)(m s−1)(g cm−3)εδ
Mudstone (top layer)440023002.65--
Shale reservoir400021002.550.300.20
Limestone (bottom layer)550035002.70--
Table 1.

Elasticity and anisotropy parameters of a VTI model for shale reservoir

VpVsρεδ
Medium parameters(m s−1)(m s−1)(g cm−3)εδ
Mudstone (top layer)440023002.65--
Shale reservoir400021002.550.300.20
Limestone (bottom layer)550035002.70--
VpVsρεδ
Medium parameters(m s−1)(m s−1)(g cm−3)εδ
Mudstone (top layer)440023002.65--
Shale reservoir400021002.550.300.20
Limestone (bottom layer)550035002.70--

The top and bottom reflection coefficients of the shale reservoir calculated by the new formula (red dotted curves) and those by Rüger's approximation (blue curves) vary significantly with the incident angles, as shown in figure 7. Furthermore, there is a great difference between the anisotropic reflection coefficients calculated by the new formula and the isotropic reflection coefficients calculated by the Aki–Richards’ approximation (black curve), especially when the incident angle exceeds 20°. However, the results calculated by the new four-term formula are well consistent with those calculated by Rüger's approximation. The average relative error of the bottom-interface reflection coefficients is ∼6.61% between the new four-term formula and Rüger's approximation, and the average relative error of the top-interface reflection coefficients is about 0.17% between the new four-term formula and Rüger's approximation. Both are much smaller than the difference between the anisotropic and the isotropic reflection coefficients. The calculation results show that the new four-term formula can be adopted as a concise version of Rüger's approximation and can be used to calculate the AVO/AVA curves of VTI media.

Similarly, based on the logging anisotropy and elasticity parameters of the reservoir section of a shale gas well, we calculated a sequence of reflection coefficients and used the convolution model to synthesize these prestack AVA gathers, as shown in figures. 8 and 9. Here, we choose the Ricker wavelet with a dominant frequency of 30 Hz to produce synthetic seismic records for further analysis. The characteristics of the prestack AVA traces (figure 8a) synthesized by the new four-term formula are consistent with the results of Rüger's approximation (figure 8b), and the residual waveform energy of two synthetics is very small (figure 8d). However, there are obvious differences between the AVA gathers synthesized by the Aki–Richards’ approximation (figure 8 (c)) and the ones synthesized by the new 4-term formula, and the residual waveform energy of two synthetics is large (figure 8 (e)), especially when the incidence angle exceeds 20°. It shows that the prestack AVA gathers synthesized by the new 4-term formula derived in this paper fully retain the contribution of the VTI anisotropy parameters of the shale reservoir section of this well to the overall reflection coefficients. As an important supplement, it can be seen from figure 8f and g that there is no significant difference between the reflection coefficient gathers versus the incidence angle θ calculated by Rüger's formula and that calculated by the new 4-term approximation formula, and their residual energy that shown in figure 8h is very small and almost negligible.

Prestack AVA synthetic gathers based on the logging elasticity and anisotropy parameters of the reservoir section of a shale gas well: (a) prestack AVA gathers synthesized by the new four-term approximation, (b) prestack AVA gathers synthesized by Rüger's approximation, (c) prestack AVA gathers synthesized by Aki–Richards’ approximation, (d) waveform residuals subtracted from (a) and (b), (e) waveform residuals subtracted from (a) and (c), (f) reflectivity gathers versus the incidence angle θ by Rüger's approximation, (g) reflectivity gathers versus the incidence angle θ by the new four-term approximation, and (h) reflectivity residuals subtracted from (f) and (g)
Figure 8.

Prestack AVA synthetic gathers based on the logging elasticity and anisotropy parameters of the reservoir section of a shale gas well: (a) prestack AVA gathers synthesized by the new four-term approximation, (b) prestack AVA gathers synthesized by Rüger's approximation, (c) prestack AVA gathers synthesized by Aki–Richards’ approximation, (d) waveform residuals subtracted from (a) and (b), (e) waveform residuals subtracted from (a) and (c), (f) reflectivity gathers versus the incidence angle θ by Rüger's approximation, (g) reflectivity gathers versus the incidence angle θ by the new four-term approximation, and (h) reflectivity residuals subtracted from (f) and (g)

The logging anisotropy and elasticity parameters of the reservoir section of a shale gas well including δ, ε, Vp0, Vs0, and ρ from the left to right column, respectively.
Figure 9.

The logging anisotropy and elasticity parameters of the reservoir section of a shale gas well including δ, ε, Vp0, Vs0, and ρ from the left to right column, respectively.

Therefore, the new 4-term formula derived in this paper can be used for the inversion of anisotropy parameters of VTI media in theory. In addition, compared to Rüger's approximation, the new formula has only four reflectivity terms, which alleviates the ill-posed problem of simultaneous solution of multiple parameters and reduces the number of the division groups of prestack angle gathers required for AVO/AVA inversion.

Furthermore, based on the derived equation (19) and its numerical simulation results, it can be seen that two conditions need to be met at least when applying equation (19) for seismic VTI inversion: (i) the incident angle of prestack seismic data is >20°, and the relatively accurate low-frequency initial model and seismic wavelet can be obtained; and (ii) the actual prestack gathers near the borehole can be well calibrated with the synthesized prestack gathers from logging data on VTI-AVO features.

4. Conclusions

In this paper, we introduce a theoretical calculation method of VTI anisotropy parameters of shale with horizontal fractures through the petrophysical modeling, and propose a new four-term approximate formula of VTI reflection coefficients. Using the proposed method, we compared the calculated results of Thomsen's anisotropy parameters for VTI media without and with horizontal frames based on actual logging data. Furthermore, we analyze the AVA response curves based on a VTI theoretical model for shale reservoirs and some prestack synthetic gathers based on logging data of shale reservoirs. Based on these studies, we have gained three insights in the following:

  • The O3W-S1l shales in Sichuan Basin have the obvious anisotropy of VTI, and the proposed theoretical calculation method of anisotropy parameters of VTI media is applicable to shale reservoirs with low-angle or horizontal fractures.

  • The testing of AVO confirms that the calculation accuracy of the new four-term approximation formula is close to that of Rüger's formula in the calculation of VTI reflection coefficients, which means it can be used as an alternative approximate method to characterize the AVO responses of VTI media.

  • The new four-term approximate formula for VTI reflection coefficients has a standard reflectivity form and only four unknown parameters, which is beneficial for obtaining Thomsen's anisotropy parameters of shale reservoirs using the traditional AVO inversion frameworks.

Acknowledgements

This research is supported by National Natural Science Foundation of China (grant nos. U19B6003, 42030103). Thanks to the Sinopec Exploration Company for providing these research data.

Conflict of interest statement. None declared.

References

Aki
K.
,
Richards
P.G.
,
1980
.
Quantitative Seismology: Theory and Methods
,
Freeman
,
New York
.

Bakulin
A.
,
Grechka
V.
,
Tsvankin
I.
,
2000
.
Estimation of fracture parameters from reflection seismic data-Part I: HTI model due to a single fracture set
,
Geophysics
,
65
,
1788
1802
.

Chen
H.
,
2015
.
Study on Methodology of Pre-stack Seismic Inversion for Fractured Reservoirs Based on Rock Physics
:
China University of Petroleum (East Chian)
,
Shandong
.

Deng
J.
,
Shen
H.
,
Xu
Z.
,
Ma
Z.
,
Zhao
Q.
,
Li
C.
,
2014
.
Dynamic elastic properties of the Wufeng-Longmaxi formation shale in the southeast margin of the sichuan basin
,
Journal of Geophysics and Engineering
,
11
,
1
11
.

Deng
J.
,
Wang
H.
,
Zhou
H.
,
Liu
Z.
,
Song
L.
,
Wang
X.
,
2015
.
Microtexture, seismic rock physical proporties and modeling of Longmaxi Formation shale
,
Chinese Journal of Geophysics (in Chinese)
,
58
,
2123
2136
.

Engelder
T.
,
Lash
G.
,
Uzcátegui
R.
,
2009
.
Joint sets that enhance production from Middle and Upper Devonian gas shales of the Appalachian Basin
,
AAPG Bulletin
,
93
,
857
889
.

Guo
Z.
,
Liu
X.
,
2018
.
Seismic rock physics characterization of anisotropic shale-a Longmaxi Shale case study
,
Journal of Geophysics and Engineering
,
15
,
512
526
.

Huang
X.
,
Eikrem
K.S.
,
Jakobsen
M.
,
Nævdal
G.
,
2020
.
Bayesian full-waveform inversion in anisotropic elastic media using the iterated extended Kalman filter
,
Geophysics
,
85
,
125
139
.

Hudson
J.A.
,
1980
.
Overall properties of a cracked solid
,
Mathematical Proceedings of the Cambridge Philosophical Society
,
88
,
371
384
.

Kwon
O.
,
Kronenberg
A.K.
,
Gangi
A.F.
,
Johnson
B.
,
Herbert
B.E.
,
2004
.
Permeability of illite bearing shale: 1. Anisotropy and effects of clay content and loading
,
Journal of Geophysical Research
,
109
,
1
19
.

Li
Q.
,
Wang
H.
,
Yang
X.
,
Ma
S.
,
Liu
X.
,
Li
J.
,
2022
.
Seismic inversion and fracture prediction in tilted transversely isotropic media
,
Journal of Geophysics and Engineering
,
19
,
1320
1339
.

Li
S.
,
Wen
X.
,
Pu
Y.
,
Miao
Z.
,
Yuan
M.
,
2023
.
Fracture prediction based on attenuative anisotropy theory and its application to a shale gas reservoir
,
Journal of Geophysics and Engineering
,
20
,
196
210
.

Liu
W.
et al.
2020
.
Microstructure characteristics of Wufeng-Longmaxi shale gas reservoirs with different depth, southeastern Sichuan Basin
,
Petroleum Geology & Experiment
,
42
,
378
386
.

Liu
Z.
,
Zhang
F.
,
Li
X.
,
2019
.
Elastic anisotropy and its influencing factors in organic-rich marine shale of southern China
,
Science China (Earth Sciences)
,
62
,
1805
1818
.

Lonardelli
I.
,
Wenk
H.
,
Ren
Y.
,
2007
.
Preferred orientation and elastic anisotropy in shales
,
Geophysics
,
72
,
D33
D40
.

Ma
Z.
,
Yin
X.
,
Zong
Z.
,
Song
M.
,
Wang
Y.
,
Chen
Z.
,
Li
J.
,
Liu
X.
,
2019
.
Azimuthally variation of elastic impedances for fracture estimation
,
Journal of Petroleum Science and Engineering
,
181
,
1
13
.

Pan
X.
,
Li
L.
,
Zhang
G.
,
2020
.
Multiscale frequency-domain seismic inversion for fracture weakness
,
Journal of Petroleum Science and Engineering
,
195
,
107845
.

Pšenčík
I.
,
Martins
J.L.
,
2001
.
Properties of weak contrast PP reflection/transmission coefficients for weakly anisotropic elastic media
,
Studia Geophysica et Geodaetica
,
45
,
176
199
.

Rüger
A.
,
1997
.
P-wave reflection coefficients for transversely isotropic models with vertical and horizontal axis of symmetry
,
Geophysics
,
62
,
713
722
.

Rüger
A.
,
1998
.
Variation of P-wave reflectivity with offset and azimuth in anisotropic media
,
Geophysics
,
63
,
935
947
.

Rüger
A.
,
Tsvankin
I.
,
1997
.
Using AVO for fracture detection: analytic basis and practical solutions
,
The Leading Edge
,
16
,
1429
1434
.

Schoenberg
M.
,
1980
.
Elastic wave behavior across linear slip interfaces
,
The Journal of the Acoustical Society of America
,
68
,
1516
1521
.

Schoenberg
M.
,
1983
.
Reflection of elastic waves from periodically stratified media with interfacial slip
,
Geophysical Prospecting
,
31
,
265
292
.

Schoenberg
M.
,
Sayers
M.C.
,
1995
.
Seismic anisotropy of fractured rock
,
Geophysics
,
60
,
204
211
.

Shi
P.
,
Yuan
S.
,
Wang
T.
,
Wang
Y.
,
Liu
T.
,
2018
.
Fracture identification in a tight sandstone reservoir: a seismic anisotropy and automatic multisensitive attribute fusion framework
,
IEEE Geoscience and Remote Sensing Letters
,
15
,
1525
1529
.

Sondergeld
C.H.
,
Rai
C.S.
,
2011
.
Elastic anisotropy of shales
,
The Leading Edge
,
30
,
324
331
.

Thomsen
L.
,
1986
.
Weak elastic anisotropy
,
Geophysics
,
51
,
1954
1966
.

Wang
Y.
,
Liu
Y.
,
Zhang
M.
,
2017
.
Seismic Equivalent Medium Theory for Fractured Anisotropy. 1st edn (in Chinese)
:
Science Press
,
Beijing
.

Wang
Z.
,
2002
.
Seismic anisotropy in sedimentary rocks, part 2: laboratory data
,
Geophysics
,
67
,
1423
1440
.

Zhang
F.
,
Li
X.
,
Qian
K.
,
2017
.
Estimation of anisotropy parameters for shales based on an improved rock physics model, part 1: theory
,
Journal of Geophysics and Engineering
,
14
,
143
158
.

Zhang
F.
,
Wang
L.
,
Li
X.
,
2020
.
Characterization of a shale-gas reservoir based on a seismic amplitude variation with offset inversion for transverse isotropy with vertical axis of symmetry media and quantitative seismic interpretation
,
Interpretation
,
8
,
SA11
SA23
.

Zhao
L.
,
Wang
Y.
,
Liu
X.
,
Zhang
J.
,
Liu
Y.
,
Qin
X.
,
Li
K.
,
Geng
J.
,
2020
.
Depositional impact on the elastic characteristics of the organic shale reservoir and its seismic application: a case study of the Longmaxi-Wufeng Shale in the Fuling gas field, Sichuan Basin
,
Geophysics
,
85
,
B23
B33
.

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