Experimental Investigation of Stress Sensitivity of Elastic Wave Velocities for Anisotropic Shale in Wufeng–Longmaxi Formation

Abstract

The shale of the Wufeng–Longmaxi formation in the Sichuan Basin is the preferred layer for shale gas exploration in China, and its petrophysical characteristics are the key to geological and engineering sweet spot prediction. However, the characteristics and impact mechanisms of its acoustic wave velocity and elastic anisotropy are currently unclear. In this paper, the Wufeng–Longmaxi shale is taken as the research object, and the P-wave and S-wave velocities of the samples are tested under the loading and unloading processes of confining pressure. The stress sensitivity variations in parameters such as wave velocity, wave velocity ratio, and anisotropy are discussed. P-wave and S-wave anisotropy parameters are correlated under different pressure conditions. X-ray diffraction, casting thin sections, scanning electron microscopy, micron CT scanning, and other analytical techniques are used to explore the mechanisms of stress sensitivity of elastic parameters. The research results indicate that: (1) the acoustic velocities of samples from different angles are V90° > V45° > V0°, and there is a positive correlation between the wave velocity and the confining pressure. After unloading the confining pressure, irreversible plastic deformation occurs due to the closure of some microfractures in the rock core, causing the wave velocity to be higher than the initial value. (2) The stress sensitivity coefficient of the P-wave (The mean is 3.00 m·s⁻¹·MPa⁻¹) is higher than that of the S-wave (the mean is 1.23 m·s⁻¹·MPa⁻¹), and the stress sensitivity coefficient of the compacted stage (the mean is 3.02 m·s⁻¹·MPa⁻¹) is higher than that of the elastic stage (the mean is 1.21 m·s⁻¹·MPa⁻¹). (3) The anisotropy of the P-wave and S-wave is negatively correlated with the confining pressure. When the confining pressure is loaded to 65 MPa, the change rate of the P-wave anisotropy coefficient is 23%, and its stress sensitivity is higher than that of S-wave anisotropy coefficient (the change rate is 13.7%). After unloading the confining pressure, the degree of anisotropy is reduced due to the closure of some microfractures. The empirical formula of P-wave and S-wave anisotropy parameters under different pressures is established through linear regression, which can provide a reference for mutual predictions. (4) The variation in wave velocity anisotropy with stress can be divided into stress and material anisotropy, which are related to the directional arrangement of microfractures and clay minerals, respectively. The quantitative characterization of shale anisotropy can be realized by evaluating the development degree of reservoir fractures and mineral components, providing a reference for logging interpretations, sweet spot prediction, and fracturing construction of shale gas reservoirs.

1. Introduction

Ultrasonic testing is usually one of the most effective methods for non-destructive characterization of rock physical characteristics, wellbore stability evaluation, reservoir geological characteristics inversion, and fracture identification [1]. The elastic and mechanical properties of organic shale usually exhibit transverse isotropic symmetry [2]. Due to the complex development of the bedding and stress environment, it usually exhibits high anisotropy, significantly affecting the propagation of ultrasonic waves and the direction of fracture propagation during fracturing [3]. However, due to limited experimental equipment and conditions, research on rock physics, mechanical response, and the anisotropic characteristics of unconventional reservoirs under complex stress is relatively lagging. This dramatically limits the deployment of horizontal wells, wellbore design, and the development of fracturing modification plans in production practice [4].

The elastic properties and mechanical behavior of shale have attracted great interest in geological engineering [5,6,7,8,9,10,11]. Vernik et al. [12,13] pointed out that there is a high correlation between in situ stress, porosity, permeability parameters, organic matter content, and shale acoustic wave velocity. Sondergeld [14] found that the acoustic wave velocity of shale increases with increasing effective stress. Through the experimental analysis of 17 saturated shales, Wang [15] found that the lower the porosity, the greater the anisotropy of the samples. In general, preferential clay minerals and microfractures can highly affect the anisotropy of shale [16,17,18]. Based on the experimental results of 51 samples, Ma et al. [19] pointed out that the P-wave velocity and S-wave velocity were exponentially correlated with the effective pressure. Ma et al. [20] found that the anisotropy parameter of shale decreases with increasing effective stress, and the anisotropy of favorable reservoirs is much larger than that of unfavorable reservoirs. Deng et al. [21] found that the P- and S-wave velocities increase with an increasing bedding angle or axial stress. Zhai et al. [22] conducted differential rock physics experiments on the Longmaxi Formation shale through pseudo triaxial stress loading experiments and concluded that applied isotropic stress to the shale can increase its P- and S-wave velocity, elastic modulus, and Poisson’s ratio; however, the anisotropy degree also weakens.

In recent years, significant attention has been paid to the rock mechanical properties and elastic characteristics of the Wufeng–Longmaxi formation shale in the Sichuan Basin of China [23,24,25,26,27,28]. However, the characteristics and mechanisms of acoustic anisotropy are not clear, and the relative relationship between P- and S-wave anisotropy is lacking. The influence of confining pressure changes on acoustic wave velocity, anisotropy, the stiffness coefficient, and other parameters of deep Wufeng–Longmaxi Formation shale require further study [29]. In order to clarify the influence of the law of elastic characteristics and stress sensitivity characteristics of deep shale reservoirs, this article takes the Wufeng–Longmaxi formation shale in the Sichuan Basin as the research object; tests the P- and S-wave velocities of samples during the loading and unloading processes under confining pressure; explores the stress-sensitive changes in the acoustic wave velocity, wave velocity ratio, elastic stiffness coefficient, and anisotropic parameters; establishes correlation equations for the anisotropy parameters of P- and S-waves under different pressure conditions; reveals the impact mechanism; and provides a reference for shale gas reservoir evaluation, dessert prediction, and fracturing construction.

2. Experimental Specimen and Test Procedure

2.1. Sample Description and Preparation

The experimental samples were from the deep Wufeng–Longmaxi formation shale in the Sichuan Basin. The TOC of the samples (sample numbers: L1, L2, L3, L4) was 0.6~2.7%, the porosity was 1.72~6.14%, the density was 2.53~2.74 g·cm⁻³, the brittle mineral content (quartz, feldspar, carbonate, pyrite) was 46~82%, and the content of clay minerals was 16~52%, mainly including illite and chlorite (

Table 1. Experimental samples.

Sample Number	Formation	Depth, m	TOC, wt%	Porosity, %	Density, g·cm−3	Type and Content of Minerals,%
						Quartz	Feldspar	Calcite	Dolomite	Pyrite	Clay
L1	Longmaxi	3936.3	0.6	1.72	2.74	38	8	0	0	2	52
L2	Longmaxi	3972.4	1.9	3.64	2.63	19	7	20	16	2	36
L3	Longmaxi	3522.5	1.8	6.14	2.53	35	16	17	0	2	30
L4	Wufeng	3553.2	2.7	4.35	2.59	18	3	50	11	2	16

).

The elastic and mechanical properties of shale usually show transversely isotropic symmetry (TI). In the laboratory, scientists often use the multi-core column method [12,30] or single plug method [31,32,33,34] to study the degree of anisotropy of rock [2]. In this paper, the multi-core column method was used (see Figure 1). The studied samples were drilled in three directions: parallel to the bedding direction (perpendicular to the symmetry axis), vertical to the bedding direction (parallel to the symmetry axis), and at an angle of 45° with the symmetry axis. The end face of the core was polished, and a sensor was installed at the top and bottom, which recorded the transmission time of the P-wave or S-wave.

2.2. Test Equipment and Procedures

This paper mainly discusses the stress sensitivity of the elastic wave velocity of dry shale. Before the experiment, the samples were evenly dried in an oven at 80 °C for 48 h and placed in a humid open-air environment for more than 24 h to obtain a core containing a small amount of water so as to avoid the damage of clay dehydration to the skeleton [35]. The BenchLab 7000 (BenchLab. Co., Begawan, Brunei) automatic porosity–permeability and velocity measurement system introduced by the NER Company in the United States was used to measure the rock density, physical properties, complex resistivity, and acoustic wave velocity under high temperature and pressure conditions (Figure 2). The maximum design confining pressure was 200 MPa, the pore pressure was 70 MPa, and the temperature was 150 °C. The main frequency of the longitudinal wave transducer of the velocity measurement module was 1 MHz, the main frequency of the transverse wave transducer was 0.5 MHz, and the experimental error was ±0.5%. The measurement principle was the ultrasonic pulse penetration method [36]. Nine acoustic wave velocities were measured for each group of samples under the same pressure conditions, including V_p0, V_s10, and V_s20 of the ultrasonic wave propagating perpendicular to the bedding direction; V_p90, V_s190, and V_s290 of the ultrasonic wave propagating parallel to the bedding direction; and the ultrasonic wave propagating at an angle of 45° with the symmetry axis of V_p45, V_s145, and V_s245, where V_s1 and V_s2 are mutually perpendicular polarization directions.

The process of shale drilling, fracturing, and gas injection is accompanied by stress fluctuations in the reservoir environment [37,38,39]. The stress sensitivity experiment simulated two stress change stages related to the depletion and recovery of the reservoir environment [40]. That is to say, the difference in stress loading and unloading values was used to study the hysteresis behavior of its elastic and mechanical properties. The buried depth of deep shale in the study area was more than 3500 m, and the measured formation pressure of the Wufeng–Longmaxi formation was 60~70 MPa. In order to simulate the high in situ stress field of the reservoir, the maximum confining pressure was set to 65 MPa. During the experiment, the effective pressure was increased from 5 MPa to 65 MPa, then the recovery pressure was increased to 5 MPa. The longitudinal wave and fast and slow shear wave velocities of the core were measured every 5 MPa. After each test point was pressurized, the loading state was maintained for more than 30 min to ensure the rock sample was compressed and stable. Then, the data were collected. Each group of experiments was measured three times to take the average value to reduce the experimental error.

2.3. Related Concepts and Formula Calculation

The following formula was used to calculate the P-wave and S-wave velocities of the sample:

$$V=\frac{L}{t-t_0}\times 1000=\frac{L}{\Delta t}\times 1000$$

(1)

where V is the P-wave or S-wave velocity of the sample in m/s, L is the sample length in mm, t is the propagation termination time of the P-wave or S-wave in the sample in m/s, t₀ is the starting time of the propagation of the P-wave or S-wave in the sample in m/s and ∆t is the propagation time of the P-wave or S-wave in the sample in m/s.

(1)

The elastic stiffness coefficient

The elastic stiffness coefficients of geological materials are usually determined by the wave velocity of ultrasonic testing [41]. Among them, the elastic stiffness coefficients are the fundamental physical quantity used to describe the elastic deformation form of the material under the action of an external force. The shale of transversely isotropic medium (VTI medium) can be described by six independent elastic stiffness coefficients representing the relationship between stress and strain of the material in different directions after being subjected to the applied stress. The relationship between the velocity and elastic stiffness coefficients is as follows [42]:

$$C_{11}=\rho V_{S{90}^{\circ }}^2$$

(2)

$$C_{33}=\rho V_{P{0}^{\circ }}^2$$

(3)

$$C_{44}=\rho V_{S{0}^{\circ }}^2$$

(4)

$$C_{66}=\rho V_{S{90}^{\circ }}^2$$

(5)

$$C_{12}=C_{11}-2C_{66}$$

(6)

$$C_{13}=\sqrt{\left(C_{11}+C_{44}-2\rho V_{P{45}^{\circ }}^2\right)\left(C_{33}+C_{44}-2\rho V_{P{45}^{\circ }}^2\right)}-C_{44}$$

(7)

where ρ is the rock density.

(2)

The anisotropy parameters

The Thomsen anisotropy parameter [43] is widely used to describe the macroscopic elastic anisotropy characteristics. The calculation formulas are written as follows:

$$\epsilon =\frac{C_{11}-C_{33}}{2C_{33}}$$

(8)

$$\gamma =\frac{C_{66}-C_{44}}{2C_{44}}$$

(9)

$$\delta =\frac{{\left(C_{13}+C_{44}\right)}^2-{\left(C_{33}-C_{44}\right)}^2}{2C_{33}\left(C_{33}-C_{44}\right)}$$

(10)

In the above, ε and γ describe the anisotropic strengths of the P-wave and S-wave, respectively. An ellipticity coefficient δ indicates the anisotropic parameters of P-wave and S-wave coupling and is closely related to seismic data processing [44]. This paper discusses only the anisotropic parameters of the P-wave and S-wave.

3. Experimental Results

3.1. The Effect of Confining Pressure on Acoustic Wave Velocity

The acoustic wave velocity measurement of the sample is shown in

Table 2. Acoustic wave velocity measurement results of shale samples (5 MPa).

Sample Number	Vp, m/s			Vs1, m/s			Vs2, m/s			Vp/Vs1			Vp/Vs2
	0°	45°	90°	0°	45°	90°	0°	45°	90°	0°	45°	90°	0°	45°	90°
L1	4298	4897	5454	2881	3211	3401	2771	3052	3413	1.49	1.53	1.60	1.55	1.60	1.60
L2	4352	4610	4737	2782	2840	2910	2824	2878	2867	1.56	1.62	1.63	1.54	1.60	1.65
L3	4411	4586	4850	2772	2884	2942	2720	2800	2827	1.59	1.59	1.65	1.62	1.64	1.72
L4	4450	4588	4821	2656	2730	2793	2701	2786	2849	1.68	1.68	1.73	1.65	1.65	1.69

(only the test value of the initial pressure point was selected). The change in the longitudinal and transverse wave velocity under the confining pressure is shown in Figure 3 and Figure 4, which show the velocity variation of the shale at 0°, 45°, and 90° angles under the pressure and after pressure relief. The experimental results show that V_P was significantly higher than V_S1 and V_S2. The V_S1 and V_S2 of the different samples overall follows V_s10 ≈ V_s20, V_s145 ≈ V_s245, and V_s190 ≈ V_s290 and confirms the transverse isotropic characteristics of the shale. The study of fast and slow shear waves is of great significance. S-wave S₂ usually demonstrated the effect of a fracture dip angle. The relative relationship between Vp/V_S2, V_S1, and V_S2 can be used to preliminarily judge and invert the dominant fracture dip angle [4]. Because the experimental data is too redundant, this paper only uses V_S1 for the analysis.

Under the same confining pressure, the P-wave and S-wave velocities showed the same variation, namely V_p90 > V_p45 > V_p0, V_s90 > V_s45 > V_s0. The P and S-waves increased exponentially with increasing confining pressure, which was caused by decreased pore structure space and increased particle contact stiffness after stress compression [45]. At the same time, both P- and S-waves had noticeable stress sensitivity. When the confining pressure was less than 40 MPa, the pressure increase had a relatively significant influence on the velocity, and the wave velocity gradually tended to become stable after exceeding 40 MPa. When the confining pressure increased from 5 MPa to 65 MPa, the change rate of the P-wave velocity relative to the initial wave velocity was 4.7~9.5% (0° sample), 3.4~6.0% (45° sample), and 2.2~4.7% (90° sample) (see

Table 3. The change rate of the wave velocity and wave velocity ratio after loading and unloading the samples.

Sample Number		Vp, %		Vs1, %		Vp/Vs1, %
		5~65 MPa	65~5 MPa	5~65 MPa	65~5 MPa	5~65 MPa	65~5 MPa
L1	0°	9.5	−1.7	4.7	0.7	4.6	−2.3
	45°	6.0	0.1	4.4	1.5	1.5	−1.4
	90°	4.7	1.8	4.0	−0.2	0.6	2.0
L2	0°	5.4	1.1	2.6	−0.6	2.7	1.6
	45°	3.4	−1.4	2.3	0.3	1.1	−1.6
	90°	3.2	−0.6	2.2	−1.3	1.0	0.7
L3	0°	6.1	0.9	3.0	0.1	3.0	0.9
	45°	3.5	−1.5	1.7	−0.8	1.7	−0.8
	90°	3.1	−0.4	1.8	−1.4	1.3	1.0
L4	0°	4.7	−0.8	3.8	−0.1	0.9	−0.7
	45°	4.0	−1.5	3.0	1.4	0.9	−2.9
	90°	2.2	−1.4	2.5	0.4	−0.3	−1.9

and Figure 5). Correspondingly, the change rate of the S-wave velocity relative to the initial wave velocity was 2.6~4.7% (0° sample), 1.7~4.4% (45° sample), and 1.8~4.0% (90° sample). The above-mentioned variation shows that the change rate of the P-wave velocity with confining pressure was significantly higher than that of the S-wave, presenting 0° > 45° > 90°.

When the confining pressure was reduced from 65 MPa to 5 MPa, the P- and S-wave velocity variation relative to the initial value was relatively complex (see Table 3 and Figure 5). The difference in wave velocity values between the unloading and loading process was slight, and some values were lower than the initial value. The change rate was generally less than 2%, indicating that the sample was still in the elastic deformation stage. Some wave velocity values increased relative to the initial value, and the change rate was still less than 2%, which may have been caused by the permanent closure of inelastic microfractures or micropores in organic matter during loading. When the sample returned to the initial stress state after pressure relief, the physical properties of the medium were difficult to recover to the initial level. Permanent damage to the reservoir is usually referred to as the stress-sensitive hysteresis effect [46,47]. The stress and strain of ideal elastic porous media are entirely reversible. The stress-sensitive hysteresis effect of the irrational elastic body of the shale reservoir produces a certain irreversible plastic deformation [40].

3.2. The Effect of Confining Pressure on the Wave Velocity Ratio

The P-wave and S-wave velocity ratio (V_p/V_s) is an important index for characterizing reservoir physical mechanics and sensitive parameters. It is usually used to index mineralogical content, organic matter richness, and fluid content [48,49,50], which is convenient for resevoir identification, formation evaluation, and geological model establishment. Figure 6 shows the variation in the P-wave and S-wave velocity ratio under confining pressure. Under the initial confining pressure, V_p/V_s was 1.49~1.73, and V_p90°/V_s90° > V_p45°/V_s45° > V_p0°/V_s0° or V_p90°/V_s90° > V_p45°/V_s45° = V_p0°/V_s0°(see Table 2). The wave velocity ratio increased with increasing confining pressure, indicating that the stress response of the P-wave was more significant than that of the S-wave. This further verifies that the stress sensitivity of the P-wave was stronger. With the increase in confining pressure, the relative relationship between the aspect ratio of different bedding angles was inconsistent. Two trends of V_p90°/V_s90° > V_p45°/V_s45° > V_p0°/V_s0° or V_p90°/V_s90° > V_p0°/V_s0° > V_p45°/V_s45° reflect that the wave velocity ratio also had significant anisotropy under the same pressure conditions. When the confining pressure increased from 5 MPa to 65 MPa, the change rate of the wave velocity ratio of the 0° sample was 0.9~4.6%, the 45°sample was 0.9~1.7%, and the 90°sample was −0.3~1.3% (see Table 3). All of these are expressed as V_p0°/V_s0° > V_p45°/V_s45° > V_p90°/V_s90°, that is, the increase in the wave velocity ratio V_p/V_s of the shale perpendicular to the bedding was the largest. When the confining pressure decreased from 65 MPa to 5 MPa, part of the wave velocity ratio V_p/V_s increased relative to the initial value (see Table 3), which was caused by the increase in wave velocity due to irreversible plastic deformation.

3.3. The Effect of Confining Pressure on the Elastic Stiffness Coefficients

The calculation results of the elastic stiffness coefficient of the experimental sample are shown in

Table 4. Calculation results of the elastic stiffness coefficients of shale samples (5 MPa).

Sample Number	C11, GPa	C33, GPa	C44, GPa	C66, GPa	C12, GPa	C13, GPa
L1	81.50	50.62	22.74	31.70	18.10	16.95
L2	59.02	49.81	20.35	22.27	14.47	16.38
L3	59.51	49.23	19.45	21.90	15.72	12.73
L4	60.19	51.29	18.27	20.20	19.78	16.47

(only the test values of the initial pressure point were selected). The change trend and change rate of the elastic stiffness coefficient C_ij with the confining pressure are shown in

Table 5. The change rate of elastic stiffness coefficients after loading and unloading tests of samples.

Sample Number	C11, %		C33, %		C44, %		C66, %		C12, %		C13, %
	5~65 MPa	65~5 MPa	5~65 MPa	65~5 MPa	5~65 MPa	65~5 MPa	5~65 MPa	65~5 MPa	5~65 MPa	65~5 MPa	5~65 MPa	65~5 MPa
L1	9.5	3.7	20.0	−3.3	9.7	1.3	8.2	−0.3	14.0	18.0	21.6	−12.9
L2	6.6	−1.2	11.1	2.1	5.2	−1.1	4.4	−2.7	13.1	3.2	6.4	−16.1
L3	6.3	−0.9	12.6	1.9	6.1	0.1	3.6	−2.8	13.7	4.6	3.2	−26.9
L4	4.4	−2.9	9.7	−1.6	7.8	−0.2	5.1	0.9	2.9	−10.5	14.5	−11.7

and Figure 7. The relationship between elastic stiffness coefficients was C₁₁ > C₃₃ > C₆₆ > C₄₄ > C₁₂ > C₁₃, and the L2 sample appeared as C₁₃ > C₁₂, reflecting that the elastic stiffness coefficient of the shale samples parallel to the bedding was higher than that of the vertical bedding samples, and the elastic stiffness coefficient showed a specific anisotropy. With an increase in the confining pressure, the elastic stiffness coefficient (C₁₁, C₃₃, C₄₄, C₆₆, C₁₂) increased monotonously, and the growth was more evident under a low confining pressure, while C₁₃ showed a fluctuating upward trend. When the confining pressure increased from 5 MPa to 65 MPa, the change rate of the elastic stiffness coefficient was the highest in C₃₃ (average 13.4%), followed by C₁₃ (average 11.4%), C₁₂ (average 10.9%), and C₆₆. The change rate was the lowest (average 5.3%), that is, C₃₃ was the most sensitive to the confining pressure, and C₆₆ was the lowest. When the confining pressure decreased from 65 MPa to 5 MPa, some elastic stiffness coefficients decreased to varying degrees, among which C₁₃ had the most significant decrease of 11.7~26.9%, further confirming that the sample underwent irreversible deformation under the action of the external force.

4. Discussion

4.1. The Stress Sensitivity of the Acoustic Wave Velocity

The relationship between the rock acoustic wave velocity and the loading stress is called wave velocity stress sensitivity [51]. The stress sensitivity of the wave velocity can be used to characterize fracture closures and propagation in unconventional reservoirs. In this paper, the stress sensitivity coefficient of the wave velocity proposed by Wang et al., 2016 [52] was used to characterize the sensitivity of the wave velocity to a change in confining pressure, and the correlation curve of the wave velocity and stress was drawn. The least square method was used to obtain a fitting line. The slope of the fitting line was defined as χ, that is, the stress sensitivity coefficient of the acoustic wave velocity. A larger χ value indicates a stronger stress sensitivity [52].

Taking the experimental core L2 as an example, the curve of the wave velocity changing with stress was plotted (Figure 8). Although the P- and S-wave velocities increased exponentially with increasing confining pressure, the curve can be divided into rapid growth and stable growth sections, with 40 MPa as the node, defined as the compacted and elastic stages [51]. The stress sensitivity coefficients of the P and S-waves were obtained through linear fitting and calculation, respectively. As shown in

Table 6. Stress sensitivity coefficients of the acoustic wave velocity.

Sample Number		Stress Sensitivity Coefficient for P-Wave (m·s−1·MPa−1)		Stress Sensitivity Coefficient for S-Wave (m·s−1·MPa−1)
		Compacted Stage	Elastic Stage	Compacted Stage	Elastic Stage
L1	0°	7.02	3.45	2.61	1.19
	45°	6.99	1.75	2.54	0.10
	90°	4.59	1.23	1.79	1.18
L2	0°	5.33	1.09	1.44	0.59
	45°	3.11	0.70	1.32	0.55
	90°	2.52	1.15	1.25	0.52
L3	0°	6.84	2.30	2.05	0.21
	45°	3.47	2.11	1.38	0.55
	90°	2.08	1.71	1.11	0.85
L4	0°	4.50	2.25	2.55	0.82
	45°	2.59	1.27	1.92	0.99
	90°	2.09	1.90	1.48	0.63

, the stress sensitivity coefficient of the P-wave was higher than that of the S-wave. The stress sensitivity coefficient of the compacted stage was higher than that of the elastic stage. The relationship between the stress sensitivity coefficient of the compacted stage with different bedding angles was 0° > 45° > 90°, that is, perpendicular to the bedding direction > 45° bedding direction > parallel to the bedding direction. With increasing confining pressure, the rock samples entered the elastic stage, the stress sensitivity coefficient gradually decreased and tended to be stable, and the difference among the rock samples with different bedding angles was small.

The bedding angle of shale is an essential factor affecting the P-wave and S-wave velocities. Shale has the characteristics of transverse isotropy and is regarded as a homogeneous elastic medium in the same lamina. The shale matrix in the direction of parallel bedding is usually highly cemented, indicating that it is smooth and the mineral filling is complete, which makes the wave velocity propagate faster. When the wave velocity propagates in the direction perpendicular to the bedding, the acoustic wave passes through the bedding interface, equivalent to passing through the lithologic mutation surface or the sedimentary discontinuity surface [53]. At the same time, with an increase in the bedding angle, the cementation degree of the shale matrix decreases, and the surface fluctuation degree also increases. Therefore, the instantaneous velocity attenuation of acoustic energy occurs [29,54], which is the fundamental reason why the P-wave and S-wave velocities always maintained V_p/s90° > V_p/s45° > V_p/s0° under the same confining pressure. Secondly, the P-wave velocity was greater than the S-wave velocity because the ability of the solid material to resist the shear strain (shear modulus) was less than its ability to resist the tensile strain (Young’s modulus).

Under a low confining pressure, the degree of fracture development determines the trend in the velocity increase. When the stress is applied to the rock sample, the fracture preferentially closes, and the softer organic matter is compacted. At the same time, the maximum stress sensitivity coefficient of the rock sample (0°) perpendicular to the bedding direction occurs because the microfractures are usually developed along the bedding direction (Figure 9), and the confining pressure is loaded perpendicular to the bedding direction. In this study, the normal force on the fracture surface was the largest. Compared with the horizontal bedding direction (90°), the fracture was more likely to close, so its stress sensitivity was higher. With increasing confining pressure, the sample enters the elastic stage, and the stiffness of the rock matrix and its constituent particles determines the increasing trend of the velocity [55]. It is difficult to deform relative to the fracture, so its stress sensitivity is significantly reduced.

4.2. Elastic Anisotropy and Influencing Factors of Shale

Shale generally exhibits anisotropy in acoustic, mechanical, and electrical properties [43,56,57]. Studying the velocity anisotropy of shale is conducive to the accurate interpretation and inversion of seismic characteristics and acoustic logging data [58,59,60,61]. It is of great significance to guide the exploration of seismic imaging, estimate the phase state of fluid inside the rock, and invert the fracture density and aspect ratio [12,16]. The main factors leading to shale velocity anisotropy are as follows: (1) TOC content and maturity, (2) directional arrangement of clay minerals, and (3) fracture development degree [62,63]. Nevertheless, shale anisotropy is strongly affected by the preferred aligned microfractures [17,18,64], which is also the focus of this study.

The calculation results of the anisotropy parameters of the sample and the variation with the confining pressure are shown in

Table 7. Calculation results and change rate of anisotropy parameters of samples.

Sample Number	Calculated Value (5 MPa)		Change Rate after Loading (5~65 MPa)		Change Rate after Unloading (65~5 MPa)
	ε	γ	ε	γ	ε	γ
L1	0.30	0.20	−32.93	−16.67	−6.06	−7.20
L2	0.09	0.05	−21.09	−12.67	−7.76	7.84
L3	0.10	0.06	−23.71	−16.59	−7.34	−13.89
L4	0.09	0.05	−14.20	−8.78	−5.77	4.58

and Figure 10, which show the variation in the P-wave anisotropy parameter ε and the S-wave anisotropy parameter γ after pressure relief. The P-wave anisotropy parameter ε of the four groups of samples was 0.09~0.30 at the initial pressure point, and the S-wave anisotropy parameter γ was 0.05~0.20. The anisotropy parameters of the P-wave and S-wave decreased with increasing confining pressure. When the confining pressure was less than 40 MPa, the increase in pressure had a relatively significant influence on the acoustic velocity anisotropy, which was caused by the fracture closure mechanism. The influence gradually tended to become stable after exceeding 40 MPa. ε was always higher than γ; when the confining pressure increased from 5 MPa to 65 MPa, the decrease in the P-wave anisotropy parameter relative to the initial value was 14.20~32.93%, the decrease in the S-wave anisotropy parameter was 8.78~16.67%, and the decrease in ε was more obvious. This demonstrates that the stress sensitivity of the P-wave anisotropy parameter was stronger. When the confining pressure was reduced from 65 MPa to 5 MPa, except for a small increase in the S-wave anisotropy parameters of rock samples L2 and L4, the anisotropy parameters of the P-wave and S-wave of the remaining samples decreased. The P-wave anisotropy parameter was reduced by 5.77~7.76%, and the S-wave anisotropy parameter decreased by 7.20~13.89%. During the process of pressure relief, the anisotropy degree was reduced due to permanent damage to some microfracture-closed rocks.

Domestic and foreign scholars have done significant research on the correlation between the anisotropy parameters of shale P-waves and S-waves. Johnston et al. [30] established a linear correlation between the anisotropy parameters of P-waves and S-waves in shale under different effective stresses. Vernik et al. [12] established a correlation between the anisotropy parameters of dry and saturated shale. Deng et al. [29] established a non-differential intersection between the ε and γ of the foreign Devonian–Mississippian shale and the deep shale of the Longmaxi formation in the Sichuan Basin under an effective stress of 50 MPa and fitted the linear correlation. The establishment of this empirical formula has practical value for the mutual prediction of ε and γ. The ε and γ data obtained in this experiment were fitted, as shown in Figure 11. Although the data points of a single pressure point were few, they still had a high linear correlation (R² is 0.8939~0.9453). When the confining pressure increased from 5 MPa to 65 MPa, the slope of the fitting formula increased from 0.4896 to 0.8497, reflecting that ε and γ decrease with increasing pressure, and the difference between the two decreased. Establishing correlation formulas under different pressures can provide a reference for the calculation of in situ stress and the mutual prediction of anisotropic parameters in shale reservoirs.

The curve of the wave velocity with stress can be divided into compacted and elastic stages. The stress sensitivity of the wave velocity in the compacted stage is higher than that of the elastic stage, and the stress sensitivity of the rock sample perpendicular to the bedding direction is higher than that parallel to the bedding direction. The microfractures parallel to the bedding will close rapidly at this stage with stress loading. In contrast, the fracture development degree in the vertical bedding direction is low; therefore, the stress sensitivity of wave velocity is low, leading to the different stress sensitivities to wave velocities in different directions and a rapid reduction in anisotropy [51].

Since microfractures in shale are usually parallel to the bedding development, four groups of 90° samples in the experiment were selected for micro-CT scanning, which can be used to reconstruct the three-dimensional structure and quantitatively characterize the pore throat structure inside the rock sample. CT imaging uses grayscale to distinguish mineral components and pore space. High-density minerals such as pyrite are bright white in scanning images, low-density minerals are dark gray, medium-density minerals such as carbonate and clay are light gray, and pores and fractures are black [65,66]. The three-dimensional distribution characteristics and profile of rock sample micro-CT scanning fractures are shown in Figure 12. Taking rock samples L1 and L4 as examples, irregular microfractures can be seen in the gray body of the L1 rock sample and the top view in the XY direction. In contrast, the fractures in L4 are mostly filled with medium- and high-density minerals. With further extraction of the number of microfractures and the equivalent diameter distribution (Figure 13), the equal diameter of the sample microfractures were mainly distributed in 10~65 μm, among which the equivalent diameter of the sample L1 microfractures were distributed in 15~20 μm and 20~25 μm, the frequency was the highest, which was 21 and 28, respectively, and the total number of extracted microfractures was 91. The main distribution range of the equivalent diameter of microfractures in samples L2 and L3 was the same as that of L1. Still, the frequency was significantly reduced, and the total number of microfractures was 57 and 62, respectively. The frequency of the equivalent diameter distribution of the microfractures in the L4 sample was up to 7 in the range of 40~45μm, and the total number of microfractures was only 29. In general, the microfractures in sample L1 were the most developed, followed by L2 and L3, and the microfractures in sample L4 were the least developed. Therefore, the stress sensitivity coefficient of the wave velocity, the anisotropy parameter of the P- and S-waves, and the rate of change after loading were the highest in the L1 compacted stage (Table 6 and Table 7). At the same time, the wave velocity stress sensitivity, the anisotropy parameters of the P- and S-waves, and the rate of change after loading were the lowest in sample L4 due to the low degree of fracture development (see Table 6 and Table 7), which further verified that the microfracture was the main reason for the strong wave velocity sensitivity and anisotropy. The anisotropy of the rock sample wave velocity decreased rapidly under confining pressure loading during the compacted stage, also known as stress anisotropy [51,67].

After the wave velocity enters the elastic stage with the change in stress, the stress sensitivity coefficient of the wave velocity gradually decreases and tends to become stable. The difference in the stress sensitivity coefficient of the rock samples with different bedding angles was small, and the anisotropy parameter tended to become stable with the change in stress loading. This is because the microfractures were almost completely closed at the initial pressure loading stage, and the directional arrangement of the mineral particles caused a wave velocity change in the elastic stage. The arrangement of clay minerals is the key factor causing the high anisotropy of shale [5,6]. The higher clay content in the mineral composition renders a greater anisotropy [68]. This anisotropy changes little with stress, which is called material anisotropy [67].

Figure 14 compares the anisotropy parameters of the four experimental sample groups with the confining pressure change. For the L1 sample with the highest clay mineral content (V_clay = 52.5%), the anisotropy parameters of the P-wave and S-wave were the highest, and the change under stress was the most obvious. The 0°and 90° sample sections of L1 were observed using scanning electron microscopy (Figure 15). It can be seen that the 0°sample had a directional arrangement of clay minerals. In contrast, the 90° sample had a disordered arrangement of clay and no orientation. This parallel bedding plane preferentially arranges clay, which is the root cause of shale anisotropy. In general, the directional arrangement of microfractures and clay minerals was along the direction of the shale bedding, so the anisotropic stress sensitivity of the clay content of the samples was also relatively strong. Secondly, the anisotropy parameters of the low clay content sample L4 (Vclay = 15.9%) were the smallest, and the stress sensitivity was also low. In general, the anisotropy degree of the rock samples had a significant positive correlation with the content of clay minerals, which confirms that the content of clay minerals is a crucial factor in determining the anisotropy degree of shale.

5. Conclusions

The acoustic wave velocity of shale samples with different angles was V90° > V45° > V0°. The wave velocity and wave velocity ratio were positively correlated with the confining pressure. The change rate of the P-wave velocity after compression was significantly higher than that of the S-wave. The change rate of the wave velocity and the wave velocity ratio of samples with different angles was 0° > 45° > 90°, which shows that the stress response of the P-wave was higher than that of the S-wave. The stress response of the samples with acoustic waves propagating perpendicular to the bedding was higher than that of acoustic wave propagating horizontally. The elastic stiffness coefficient was anisotropic and increased monotonically with increasing confining pressure. Among them, C₃₃ had the highest sensitivity to confining pressure, and C₆₆ has the lowest. After unloading the confining pressure, some P-wave and S-wave velocities and wave velocity ratios increased relative to their initial values, and the elastic stiffness coefficients decrease to varying degrees. This was due to the irreversible plastic deformation caused by the permanent closure of some microfractures in the rock sample during loading.

The stress sensitivity coefficient of the P-wave (the mean was 3.00 m·s⁻¹·MPa⁻¹) was higher than that of the S-wave (the mean was 1.23 m·s⁻¹·MPa⁻¹), and the stress sensitivity coefficient of the low confining pressure compacted stage (the mean is 3.02 m·s⁻¹·MPa⁻¹) was higher than that of the high confining pressure elastic stage (the mean was 1.21 m·s⁻¹·MPa⁻¹). The stress sensitivity coefficient of the rock sample at the compacted stage with different angles was 0° > 45° > 90°. The fracture development degree and matrix stiffness determined the stress sensitivity of the compacted and elastic stages, respectively. It was convenient to quantitatively evaluate the sensitivity of the acoustic wave velocity stress using the acoustic stress sensitivity coefficient of acoustic wave velocity, which has practical guiding significance for engineering operations.

The degree of anisotropy of the P-wave and S-wave decreased with increasing confining pressure. When the confining pressure was loaded to 65 MPa, the change rate of the P-wave anisotropy coefficient was 23%, and its stress sensitivity was higher than that of the S-wave anisotropy coefficient (the change rate was 13.7%). After unloading the confining pressure, the anisotropy degree decreased due to the closure of some microcracks. The linear regression method establishes the empirical expressions of ε and γ under different pressures, which can provide a reference for the calculation of in situ stress and the mutual prediction of the P-wave and S-wave anisotropy parameters in shale reservoirs.

The variation degree of shale anisotropy with stress loading can be divided into stress anisotropy and material anisotropy, which are related to aligned microfractures and clay minerals, respectively. The distribution characteristics of microfractures and the content of clay minerals were verified through CT scanning. The quantitative characterization of shale anisotropy can be realized by evaluating the development degree of the reservoir fractures and mineral components.