Effects of Background Porosity on Seismic Anisotropy in Fractured Rocks: An Experimental Study

1. Introduction

Fractures are found throughout crustal rocks, and fracture detection is of critical importance for hydrocarbon exploration, water resource management, CO₂ capture and storage, and geo-thermal exploitation (Fan et al., 2018; Li et al., 2022). Aligned fractures can cause seismic anisotropy, making it possible to obtain fracture information from seismic data. Fracture-induced anisotropy is an integrated result of fracture parameters (fracture density, scale, aspect ratio, etc.), fluid properties (fluid modulus and viscosities), and background properties (modulus, porosity, and permeability). Background matrix porosity has been thought to greatly influence P wave velocity anisotropy (Gurevich et al., 2003; Sil et al., 2011). Because fractures are present in low, mid, and high porosity hydrocarbon reservoirs, understanding the effect of background porosity on seismic anisotropy could further improve seismic fracture characterization in reservoirs of different porosities.

To date, theoretical and experimental studies have been conducted to relate rock physical properties to seismic velocity and anisotropy, and various equivalent medium theories have been proposed (Eshelby et al., 1953; Schoenberg et al., 1988; Cheng et al., 1993; Gurevich et al., 2003). These theories were based on different assumptions and derived with different methods, resulting in different predictions. Thus, there is a desire to verify and calibrate these theoretical models through laboratory experiments. Laboratory experiments typically use synthetic samples with controlled and known rock physical properties and fracture parameters. A number of experimental studies have been conducted to observe elastic wave velocity and anisotropy in fractured rocks (Rathor et al., 1995; Ass’ad et al., 1996; Wei & Di, 2008; Tillotson et al., 2012, 2014; Amalokwu et al., 2014, 2015a, 2015b, 2016, 2017; Ding et al., 2014a, 2014b, 2016, 2017, 2018, 2020, 2021;). These studies have investigated how parameters (such as fracture density, scale, and aspect ratio) can influence elastic wave velocity and anisotropy, providing significant insights into the physical mechanisms for elastic wave propagation in fractured rocks. However, despite the variations in the background porosity of fractured reservoirs, few studies have focused on the effect of background porosity on seismic anisotropy in these reservoirs.

Gurevich (2003) discussed the P wave anisotropy variation with background porosity in water saturated rocks with aligned fractures in low frequency range using numerical simulations. When there is no background porosity, P wave anisotropy ε is very small. In the absence of background porosity, fluid is ‘isolated’ in fractures. The ‘isolated’ fluid acts to ‘stiffen’ the fractures and reduce fracture compliance. This results in a relatively low P wave anisotropy for the overall fractured rock. With the increase to characteristic porosity, P wave anisotropy exhibits a sharp increase. When surrounded by pores, the fluid in fractures can escape during compression, therefore, the fracture is more compliant than in the ‘isolated’ condition. As porosity continues to increase, P wave anisotropy flattens and decreases gradually. Other theoretical studies have also discussed the influence of porosity on anisotropic parameters (Sil et al., 2011; Huang et al., 2015). However, no experimental studies have been conducted to test and verify this owing to difficulties in manufacturing artificial fractured rocks with low porosity.

This study experimentally examined the influence of background porosity on P-wave anisotropy in fractured rocks under water saturation conditions. A set of artificial fractured sandstone samples with the same fracture density (0.06), but different porosities (ranging from 15.3% to 30.8%) were built. Laboratory ultrasonic measurements were conducted to obtain the P and S wave velocities and anisotropy. We analysed the P wave anisotropy and shear wave splitting (SWS) variation with background porosity. Results were compared with theoretical models and causes of the P wave anisotropy and shear wave splitting variation with porosity were discussed.

2. Equivalent Medium Theories

2.1. High frequency (unrelaxed) moduli

High frequency theories have modelled fractures as “isolated” inclusions; thus, there was no fluid exchange between inclusions and the background medium. They assumed that if the frequency was higher than the squirt characteristic frequency, f >> f_c, the fluid pressure does not have time to equilibrate between stiff and compliant pores during the half-wave cycle, or so-called “unrelaxed” state. Cheng (1978, 1993) derived the elastic stiffness tensor of a fractured medium at high frequency (unrelaxed) conditions based on Eshelby inclusion theory. Mathematically, the Eshelby—Cheng model is given by:

$$C_{\mathit{ij}}^{\mathit{eff}}=C_{\mathit{ij}}^{\left(0\right)}-{\phi }_f{C}_{\mathit{ij}}^{\left(1\right)}\hspace{1em}(1)$$

where $C_{\mathit{ij}}^{\left(0\right)}$ is the background matrix stiffness tensor, ${\phi }_f$ is fracture porosity, and:

$$C_{\mathit{ij}}^{\left(1\right)}={\mathrm{C}}_{\mathit{ik}}^{\left(0\right)}{A}_{\mathit{kj}}\hspace{1em}(2)$$

where $A_{\mathit{ij}}$ relates applied strain $e_i^A$ to “stress-free” strain $e_i^T$ of the inclusion:

$$e_i^T=A_{\mathit{ij}}{e}_j^A\hspace{1em}(3)$$

2.2. Low frequency (relaxed) moduli

In contrast to high frequency theories, low frequency theories considered the fluid exchange between fractures and the surrounding pores and micro-cracks. They assumed that if the frequency was lower than the squirt characteristic frequency, the fluid pressure can equilibrate locally between cracks and adjacent pores. By combining linear-slip theory and the anisotropic Gassmann equation, Gurevich (2003) derived expressions for a fluid saturated cracked porous medium in the low frequency condition in the form:

where $C_{ij}^0$ is the dry modulus of the fractured medium, and α_m is Biot’s coefficient:

for i=1,2, and 3, α₄ = α₅ = α₆ = 0. Scalar M is the direct analog of Gassmann’s pore space modulus:

where $K_g$ and $K_f$ is the bulk modulus of mineral grains and pore fluid respectively.

φ is the overall porosity of the fractured rock. $K^{*}$ denotes the so-called generalized drained bulk modulus, which is defined as:

3. Sample Preparation

The method of Ding (2014) was adapted to build artificial porous sandstones which have similar mineral components, pore structure, and cementation to natural rocks. In this method, powders of silica sand, feldspar, and kaolinite, were mixed in a ball mill for 24h to ensure homogeneity. Then these mineral powders were mixed with sodium silicate solution. The mixture was then paved into a mold layer by layer. In order to make the “meso-scale” fractures, touch paper discs were spread on the surface of mineral powder mixtures when layering the mixture in the mold. Touch paper is a special kind of paper which has been soaked in saltpeter. At high temperature touch paper would be decomposed into gas and leave nearly no remain.

In order to pave the mixture in the mould in layers, a special tool was made (Figure 1). The tool is consisted of a frame, a glass plate and a threaded rod. By rotating the threaded rod, the glass plate can be moved vertically. The procedure of paving the mixture is shown is Figure 2. Afterwards a predetermined number of discs were carefully placed on top of each mixture layer. Then the next mixture layer is paved , and the same number of discs were placed on it(except the top layer).

After that, a uniaxial compressing pressure was applied in the mold for 10 minutes to make mineral grains packed compactly. The compressing direction is perpendicular to fractures and layers. During this process, the pressure was adjusted to control the porosity of samples. During the next 10 hours the sodium silicate in the mixture became consolidated. In this way a block was formed and consolidated sodium silicate can give the block initial mechanical strength. Next the block was demolded, and then sintered in a muffle oven at 900℃. During sintering，touch paper discs were decomposed, leaving penny-shaped voids as fractures. The final step is to polish the sample surface for the convenience of ultrasonic velocity measurement. The size of all cuboid-shaped samples are 69.9±0.1mm in the two directions parallel to fractures, and 49.9 ±0.1mm in the direction perpendicular to fractures.

In order to produce artificial rocks with the same fracture densities but different porosities, we made the following modifications to the method of Ding (2014a): (1) Different amounts of mineral material mixture were layered in the mould (70 mm × 70 mm) in 37 layers when building the different porosity samples. (2) The mixture of different samples was compressed into the same thickness (50±0.3 mm) by adjusting the applied pressure, thus to form blocks of the same volume. Therefore, blocks made using more mixture would experience higher loading pressure and have lower porosity. (3) 126 touch paper discs were spread on the surface of each layer (except the top) in all blocks. Fracture density is unitless and calculated by:

$$e=\frac{N{a}^3}{V}\hspace{1em}(8)$$

Where N is the number of fractures, a is the radius of fractures, V is the total volume of the sample. Thus, all blocks would have the same volume (70mm×70mm×50 mm) and contain the same number of fractures(36 layers×126 fractures/layer),making fracture density practically identical for different porosities.

For each fractured sample, a blank sample (containing no penny-shaped fractures) was made as a reference sample to provide background properties. The blank samples were created using the same method as the fractured samples, excluding the introduction of fractures.

Four groups of artificial rocks with different porosities were constructed (Figure 3a). Each group consisted of two samples, one fractured sample containing penny shaped fractures with fracture density of 0.06 and a corresponding unfractured sample. Fracture diameter was set at 3 mm, and fracture thickness was 0.06 mm. Table 1 lists main parameters of all samples. The porosity was measured using helium porosimetry. In fractured samples, the fracture porosity was calculated to be 0.75%, but the porosity difference between the fractured and blank (unfractured) samples within each group varied from 0.4% to 1.5%, larger or smaller than the fracture porosity, indicating a small difference in background porosity between fractured and blank samples. This was caused by tiny errors in the manufacturing process. A section obtained by cutting a fractured sample is shown in Figure 3b and parallelly distributed fractures can be observed. A scanning electron microscope (SEM) image of the sectoin is shown in Figure 3c-d, pore and fracture structures with controlled geometry are clearly shown.

4. Ultrasonic Measurement

P- and S-wave velocity measurements were conducted on an ultrasonic bench-top pulse transmission system at room temperature and atmosphere pressure for water saturated conditions. The measurement was based on test standard DZ/T0276.24-2015. The ultrasonic transducers, heavily damped with a central frequency of 500 kHz, were excited by a voltage spike, giving rise to broadband ultrasonic pulses. The transmitted P- and S-wave signals were recorded with a digital oscilloscope and a desktop computer, time sampling interval for all experiments was 0.04 μs for both P- and S-wave signals. The measurement error was about 0.5% for P wave velocity and 0.8% for S wave velocity. Water saturation was obtained by immersing samples in a water filled container, the container was placed in a sealed bin to extract the air. To ensure full water saturation, the saturation rate is calculated by:

Where m_sat and m_dry is the wet and dry mass of the sample, V and φ is the volume and helium porosimetry measured porosity of the sample, ρ_f is bulk density of water. P-wave velocities were measured both parallel (90°) and perpendicular (0°) to the fractures. S-wave velocities were measured parallel to the fractures, and for polarization parallel (S1 wave) and perpendicular (S2 wave) to the fractures by rotating the transducers. Then Thomsen parameter ε which represents the difference in P-wave velocity between the parallel and perpendicular directions can be calculated by:

Where V_p ‖ and V_p ⊥ are P wave velocities in parallel and perpendicular direction respectively. Shear wave splitting (SWS) in the parallel direction can be calculated as:

Where V_s1 and V_s2 are fast and slow shear wave velocities. The numerical modelling experiment of wave propagation in anisotropic media by Dellinger & Vernik (1994) showed that if the wave front is propagating parallel or perpendicular to the layering(or in this case fractures), a true phase velocity is measured in laboratory ultrasonic experiments. Therefore, in this study, P and S wave phase velocities were measured in both parallel and perpendicular direction. The transmitted signals of blank samples were analysed by Fourier transform to obtain the dominant frequency and calculate the wavelength (Figure 4). The ratio of wavelength to fracture diameter (λ/d ratio) was between 5.2 and 10.8 for the P wave and between 3.9 and 12 for the S wave. Based on experimental and numerical studies (Ding et al., 2020), for this λ/d ratio, both P and S wave may exhibit Rayleigh scattering. This is undesirable because the equivalent medium criteria are violated. This problem has been frequently encountered in previous experimental studies (Rathor, 1995; Tillotson, 2014). However, similar to previous studies, we also found some interesting correlations between experimental results and theoretical predictions.

5. Results and Discussions

5.1. Laboratory results

Table 2 gives the P and S wave velocities in the blank samples. The P wave velocity in the parallel direction is slightly higher than that in the perpendicular direction, and the S1 velocity (polarized parallel to fractures) is slightly higher than the S2 velocity (polarized perpendicular to fractures). This is due to layering and uniaxial loading during the manufacturing process.

Figure 5 shows the measured P and S wave velocities in water saturated fractured samples. Significant velocity anisotropy was evident in all fractured samples with the faster P wave being measured parallel and the slower one perpendicular. Measured velocity anisotropy is discussed further in the following paragraphs. When the porosity increased from 15.3% to 22.1%, the P wave velocities, both parallel and perpendicular, dropped sharply. However, as the porosity continued to increase, the decreasing trend in the perpendicular direction is slower than that in the parallel direction, narrowing the difference between the parallel and perpendicular directions and indicating a decrease in anisotropy. However, when porosity increased from 15.1% to 30.8%, the decrease of both S1 and S2 wave velocities were generally similar.

Figure 6 shows the P wave anisotropy parameter ε and SWS in all blank samples. Both ε parameter and SWS decreased with the increase of porosity in the blank samples. Because higher porosity samples suffered lower uniaxial loading (Table 1), thus exhibited weaker anisotropy.

We aimed to investigate the effect of background porosity on the fracture-induced P wave anisotropy and SWS. Because blank samples exhibited transverse isotropy, so we inferred that the background matrix in the fractured samples had the same anisotropy as blank samples. In fractured samples, the measured anisotropy was a superposition of background matrix anisotropy and fracture induced anisotropy. To model the anisotropy caused by fractures, the anisotropy of the background matrix have to be considered. Therefore, we approximately represented the fracture-induced P wave anisotropy and SWS by subtracting the measured ε parameter and SWS of the blank samples from that of the fractured ones. Figure 7 and Figure 8 depicted the measured results of fracture-induced ε parameter and SWS as well as theoretical predictions. When porosity increased from 15.3% to 22.1%, both measured ε parameter and SWS exhibited small fluctuations. In this porosity range, both ε parameter and SWS are not so sensitive to porosity change. However, when porosity increased to 26.1%, P wave anisotropy showed a dramatic decline, from 4.5% to 3.5%, while SWS remained fairly constant. Afterwards，as porosity increased to 30.8% ,the decline trend of P wave anisotropy continued , from 3.5%to 2.2%, but SWS fluctuated very slightly.

5.2. Modelling insights and discussions

Inputs for equivalent medium theories are given in Table 1 and Table 3. Eshelby-Cheng (1993) theory (high frequency) predicted P wave anisotropy to be close to the measured results. Although expressions of Eshelby-Cheng theory do not contain background porosity, its predicted P wave anisotropy showed a decreasing trend with increasing background porosity. Eshelby-Cheng theory can be expressed by Equations (1) and (2). Note that in Equation (1), ${\phi }_f$ refers to the fracture occupied porosity, not the background porosity, and in all fractured samples the fracture porosity are practically identical. Predicted P wave anisotropy fluctuation is mainly dominated by parameter $A_{\mathrm{ij}}$ in Equation (2), which relates the applied strain to the “stress-free” strain of the inclusion. $A_{\mathrm{ij}}$ is expressed by:

$$A_{ij}=B_{ij}^{-1}{L}_{kj}\hspace{1em}(12)$$

Where:

$$\begin{gathered} B_{ij}=G{S}_j+2H{S}_{ij}+{\lambda }_b+2{\mu }_b{\delta }_{ij}\hspace{1em}(13) \\ {L}_{ij}=-G-2H{\delta }_{ij}\hspace{1em}(14) \\ H=\mu \text{'}-{\mu }_b\hspace{1em}(15) \\ G=\lambda \text{'}-{\lambda }_b\hspace{1em}(16) \\ {S}_j={\displaystyle \sum _{i=1}^3{S}_{ij}}\hspace{1em}(17) \end{gathered}$$

Where $\lambda$_b and $\mu$_b is lame constant of the background medium, $\lambda \text{'}$ and $\mu \text{'}$ is the lame constant of the inclusion material (water). $S_{ij}$is the Eshelby S tensor, which relates the strain field due to an ellipsoidal inclusion to the “stress-free” strain in the inclusion. For very thin fractures (small aspect ratio), $S_{ij}$ can be simplified as:

$$\begin{gathered} S_{11}=S_{22}=\pi \alpha (9-5R)/16\hspace{1em}(18) \\ {S}_{33}=1-\pi \alpha R/2\hspace{1em}(19) \\ {S}_{12}=S_{21}=-\pi \alpha (3-7R)/16\hspace{1em}(20) \\ {S}_{13}=S_{23}=-\pi \alpha R/4\hspace{1em}(21) \\ {S}_{31}=S_{32}=1-2R-\pi \alpha (3-5R)/4\hspace{1em}(22) \\ {S}_{55}=S_{44}=2-\pi \alpha (3-2R)/2\hspace{1em}(23) \\ {S}_{66}=\pi \alpha (3+R)\hspace{1em}(24) \end{gathered}$$

Where $\alpha$ is the fracture aspect ratio and:

$$R=\frac{{\mu }_b}{{\lambda }_b+2{\mu }_b}\hspace{1em}(25)$$

In the fracture normal direction (direction 3), since $\alpha$ is so small($\alpha =0.018$), so $S_{33}\approx 1$ and $S_3\approx 1$. Because the shear modulus of inclusion material(water) $\mu \text{'}$ is zero, therefore $B_{33}\approx \lambda \text{'}$ and $L_{33}\approx {\lambda }_b+2{\mu }_b-\lambda \text{'}$. With the increase of porosity, $L_{33}$ decreases while $B_{33}$ stays nearly constant. This leads to the decrease of $A_{ij}$. As a result, the fracture induced perturbation in the fracture normal direction decreases with the increase of porosity. However, as the aspect ratio of fractures was extremely small in this study, the fracture induced perturbation was nearly negligible in the parallel direction of the fractured samples. As the fracture induced perturbations decreased with increasing porosity in the perpendicular direction, the stiffness difference between the parallel and perpendicular directions also decreased. This leads to the declining trend of the predicted P wave anisotropy with increasing porosity.

In Figure 8, we can see that the Eshelby-Cheng predicted SWS is higher than measured results. The shear moduli of the S₁ ($C_{44}$) and S₂ ($C_{66}$) in Eshelby-Cheng theory are:

As ${\mu }^{\mathrm{\text{'}}}=0$ in our study, therefore SWS can be approximately expressed by:

where $e$ is fracture density, and φ_f = 2παe. As the value of R in four rocks is close (between 0.25 and 0.3), SWS prediction is similar and not sensitive to porosity change.

In theoretical predictions of Gurevich (2003), with the increase of porosity, P wave anisotropy showed a decreasing trend while SWS remained practically constant. This trend is similar to the measured results. However, theoretical predictions are higher than the measured results. For water saturated conditions, the theoretically predicted P wave anisotropy of Gurevich (2003) was:

where direction 3 is the fracture normal direction. $C_{11}^{sat}$ and $C_{33}^{sat}$ are the P-wave modulus of the parallel and perpendicular directions in water saturated condition, respectively. Meanwhile $C_{11}^{dry}$ and $C_{33}^{dry}$ are the P-wave modulus of the parallel and perpendicular directions in dry condition, and:

$$C_{33}^{dry}=\frac{3\mu (\lambda +\mu )}{3\mu (\lambda +\mu )+4(4{\mu }^2+4\lambda \mu )e}{C}_{11}^{dry}\hspace{1em}(30)$$

Where λ and μ are Lame parameters of the dry background matrix.So Equation (29) can be rewritten as :

$$\epsilon =\frac{1}{2}(\frac{C_{11}^{dry}+{\alpha }_1^{2}M}{\frac{3\mu (\lambda +\mu )}{3\mu (\lambda +\mu )+4(4{\mu }^2+4\lambda \mu )e}{C}_{11}^{dry}+{\alpha }_3^{2}M}-1)\hspace{1em}(31)$$

Firstly, with the increase of porosity, the value of $\frac{3\mu (\lambda +\mu )}{3\mu (\lambda +\mu )+4(4{\mu }^2+4\lambda \mu )e}$ increases from 0.62 to 0.75. Meanwhile, according to Equation (6), when the porosity increases, scalar M decreases greatly. Since α₁ is larger than α₃ (Equation (5)), so ${\alpha }_1^{2}M$ term decreases more than ${\alpha }_3^{2}M$ term. As a result, with the increase of porosity these two factors make the $C_{11}^{sat}$ decreases more than $C_{33}^{sat}$, making P wave anisotropy decrease. Predicted P wave anisotropy is higher than the measured results. Potential explanations for this phenomenon include: (i) Under experiment frequency (0.5 MHz), fluid didn’t have enough time to flow between pores and fractures, thus the system didn’t reach pressure equilibrium in pore space，and the equilibrium is assumed by Gurevich (2003) theory. This can cause a “stiffening” effect on fractures, reducing P wave anisotropy and making measured results lower than theoretical predictions. (ii) Interactions between fractures were not included in the theory of Gurevich (2003), which has been shown to reduce the overall anisotropy in a fractured medium in theoretical studies (Cheng, 1993).

In Gurevich (2003), the expression of SWS is:

$$SWS=\frac{C_{55}^{sat}-C_{44}^{sat}}{2C_{44}^{sat}}\hspace{1em}(32)$$

Where:

$$\begin{gathered} C_{55}^{sat}=\mu \hspace{1em}(33) \\ {C}_{44}^{sat}=\mu (1-{\Delta }_T)\hspace{1em}(34) \\ {\Delta }_T=\frac{\mu Z_T}{1+\mu Z_T}\hspace{1em}(35) \\ {Z}_T=\frac{16(\lambda +2\mu )e}{3\mu (3\lambda +4\mu )}\hspace{1em}(36) \end{gathered}$$

Combining equations from (32) to (36) can obtain:

$$SWS=\frac{4e}{3}(1-\frac{\lambda }{3\lambda +4\mu })\hspace{1em}(37)$$

where $e$ is fracture density. The term $\frac{\lambda }{3\lambda +4\mu }$ was extremely small in the four fractured samples (less than 0.077), therefore, . Compared with P wave anisotropy, SWS prediction is relatively close to the measured results.