Elastic waves in fractured rocks under periodic compression

Abstract

Background

One of the current problems in studying the mechanical properties and behavior of structurally inhomogeneous media with cracks is the characterization of acoustic wave propagation. This is especially important in Geomechanics and prognosis of earthquakes.

Methods

In this work, the authors propose an approach that could simplify characterization of wave propagation in medium with cracks. It is based on homogenization procedure performed at a set of equations characterizing acoustic wave propagation in media weakened by fractures under condition of external distributed loading. Such kind of loading in most cases is close to the real one in case of consideration of Geomechanics problems.

Results

On the basis of the proposed homogenization technique, we performed characterization of elastic properties and plane acoustic waves propagation in a pre-loaded linear elastic medium weakened by a large amount of cracks. We have investigated two special cases of loading: uniaxial compression and complex compression. We have also studied how the wavespeeds depend on averaged concentration and distribution of craks.

Conclusions

Effective elastic properties were theoretically characterized for fractured media under external loading. The results revealed high dependency of the longitudinal wave propagation speed on the relation between stresses reasoned by an external loading.

Background

Studying of acoustic wave propagation in structurally inhomogeneous media weakened by a large amount of fractures nowadays is a current problem for characterization of Earth shell, earthquakes prognosis, and mining. A complex hierarchical structure of geomaterials’ inhomogeneities became a reason for consideration of simplified models of acoustic wave propagation in such materials.

In the current work, authors propose characterization of plane wave propagation in pre-loaded linear elastic inhomogeneous media weakened by a large amount of isolated fractures. The existence of cracks makes the material to be highly heterogeneous with multiscale structure. Theoretical study of wave propagation in bodies and media with fractures was presented in many research papers (Schoenberg and Sayers 1995; Stiller and Wagner 1979; Schubnel and Gueguen 2003; Kachanov 1980). For example, in (Stiller and Wagner 1979), in order to obtain elastic wave velocities in fractured pre-loaded medium, the effective elastic moduli were used obtained for the case of uniaxial compression, and alteration of wave characteristics of the aforementioned medium was explained by fracturing effects. In Schoenberg and Sayers (1995); Schubnel and Gueguen (2003), plane wave propagation was studied in elastic media weakened by open fractures; this included problems on effects of anisotropy of open crack distribution at wave propagation speeds and decay characteristics. Numerical studies were performed by Zhang (2005); Chung et al. (2016) and other. In this research, the authors applied a linear-slip displacement-discontinuity model where fracture is assumed to have a vanishing width across which the tractions are taken to be continuous; however, displacements can be discontinuous. The novelty of our approach is the fact that we take into account the friction between the crack faces. Analysis of deformation properties of brittle materials revealed (Talonov and Tulinov 1988) that the characteristics of media weakened by fractures with edges interacting under compression stress fields depend on the stress-strain state at the stage prior to crack growth initiation. Therefore, the elastic medium weakened by fractures with interacting edges becomes anisotropic even at homogeneous and isotropic fracture distribution. In accordance with all the mentioned above, here, we propose investigation of wave propagation features in media weakened by fractures under condition of external periodically distributed loading.

Methods

Let us consider a plane wave propagating through a pre-loaded elastic material weakened by closed fractures (see Fig. 1).

Fig. 1

Coordinate system. Elastic wave propagates through pre-loaded fractured medium

The governing equations of motion are

$$ \rho \frac{\partial^{2} u_{i}}{\partial t^{2}}=\frac{\partial}{\partial x_{k}}\sigma_{ik}. $$

(1)

Here, t is the time, x _k is a coordinate (k=1,2,3), u _i is a component of a displacement vector (i=1,2,3), and ρ is the density of medium. Also, we assume that the stress tensor components arising due to wave propagation in the medium are small in comparison with the external stress.

Equations. 1 have to be accompanied by components of stress tensor σ _ik and initial conditions

$$ u\vert_{t=0}=v_{i},\qquad\frac{\partial u_{i}}{\partial t}\vert_{t=0}=w_{i}. $$

(2)

Let us consider an element of volume so that its characteristic size d is large in comparison with size of a crack (d≫b). The averaged strain tensor ε _ij in elastically homogeneous volume element V=d ³ containing N cracks can be written as a sum of two terms (Kachanov 1980)

$$ \delta \varepsilon_{ik}=\delta \varepsilon_{ik}^{0}+\frac{1}{2}\int \left(n_{i} \delta U_{k}+n_{k} \delta U_{i} \right)F(Z)dZ, $$

(3)

where $\varepsilon _{ik}^{0}$ are the components of strain tensor for a solid elastic material without cracks, n _i are the components of vectors normal to the faces of cracks, F(Z) is the function of set of parameters Z such as characteristic size of crack R and its orientation, and U _i are the components of the jump of displacement averaged over the faces of cracks.

We consider a loaded elastic material with Young’s modulus E ₀ and Poisson’s ratio ν ₀ weakened by isolated closed penny-shaped cracks. Such approximation of penny-shaped fractures could be useful for the aims of modeling of real material deformation properties when such materials are weakened by a large amount of isolated microfractures (Nemat-Nasser and Hori 1999). For the case of an axially symmetric stress (σ ₃₃<σ ₂₂=σ ₁₁≤0), we place a coordinate system at the center of an arbitrary crack so that the discontinuity lies in the same plane with principal stresses σ ₁₁ and σ ₂₂. For this stressed state, the crack orientation is governed by angle ψ between the principal axis corresponding to σ ₃₃ and normal n to the face of the crack. We assume that the crack faces react according to the Coulomb rule τ _f =μ σ _n , where μ is the friction coefficient and σ _n is normal stress. The angle ψ is α ₁<ψ<α ₂. In general, both α ₁ and α ₂ are functions of principal stresses σ ₁₁, σ ₂₂, and σ ₃₃. For this case, authors (Talonov and Tulinov 1988) have shown that the increment in components of the jump of displacement may be represented as

$$ \begin{aligned} \delta U_{k}&=J_{klm}\delta \sigma_{im},\\ J_{klm}&=\frac{32\left(1-\nu_{0}^{2}\right)R^{3}}{3E_{0}\left(2-\nu_{0}\right)}\left(n_{m}\delta_{lk}-n_{i}n_{k}n_{m}-\mu n_{l}n_{m}\tau_{k}\right),\\ \tau_{k}&=\frac{F_{k}}{\vert \overrightarrow{F}\vert}, \quad F_{k}=\sigma_{lm}\left(n_{m}\delta_{lk}-n_{l}n_{k}n_{m}\right), \end{aligned} $$

(4)

where δ _lk is the Kronecker symbol.

Using Eqs. (3)–(4), the increments can be rewritten in components of macroscopic strain tensor as

$$ \delta\varepsilon_{ik}=S_{ikjl}\delta \sigma_{jl}+A_{ikjl}\delta\sigma_{jl}. $$

(5)

In Eq. (5), S _ikjl are the components of the tensor S which defines elastic properties of homogeneous volume element V and A _jlpq are the components of tensor A which defines nonlinear properties due to cracks (i,j,k,l=1,2,3). For the isotropic elastic material, the components of the tensor S can be expressed as

$$ S_{ikjl}=-\frac{\nu_{0}}{E_{0}}\delta_{ij}\delta_{kl}+\frac{1+\nu_{0}}{2E_{0}}\left(\delta_{il}\delta_{kj}+\delta_{ik}\delta_{jl}\right). $$

(6)

The components of tensor A for the case of an axillary symmetric stress and isolated randomly oriented closed penny-shaped cracks with Coulomb friction are presented in Appendix 1.

On the basis of the problem stated above, we present an approach that is useful for characterization of effective elastic properties and acoustic wave propagation in the considered media with randomly distributed fractures under the condition of external loading.

Results and discussion

Elastic property characterization

Unlike the components of tensor S, components of tensor A depend on the external stress and, as a result, under the condition of external load effects, the initially isotropic media may become anisotropic. Along with this, the deformational properties of the media with closed fractures may also be altered at change of the external stress field. Regarding these facts, we consider a case when, under the conditions of a complex stressed state (σ ₃₃<σ ₂₂=σ ₁₁≤0), one of the components of the external stress, for example, σ ₃₃, will be a periodical function with respect to coordinates. Then, in accordance with the relations presented in Appendix 1, the deformational characteristics of the media with isolated closed fractures are also periodical functions of coordinates. In order to characterize wave propagation in the media with periodically distributed elastic characteristics, a method could be used developed in works Sanchez-Palencia (1980); Chen and Fish (2000); Andrianov et al. (2008); Bakhvalov and Panasenko (2012). Authors of (Chen and Fish 2000; Andrianov et al. 2008) applied this technique to homogenization of wave propagation in periodic composite materials.

In the current work, authors propose to utilize a technique to homogenization of wave propagation for characterization of the inhomogeneous external stress field effects at wave properties of the medium weakened by a large amoung of distributed microfractures.

Now, when the problem is formulated by (1)–(6), we substitute (6) into (1) and look for the solution as

$$ u_{i}(x_{n},t)=\tilde{u}_{i}(x_{n})e^{-i\omega t}. $$

(7)

The equations of motion (1) will be rewritten as

$$ \begin{aligned} -\rho \omega^{2} \tilde{u}_{i}(x_{n})&=\frac{\partial}{\partial x_{k}}\left(C_{ikjl}(x_{n})\frac{\partial\tilde{u}_{j}(x_{n})}{\partial x_{l}}\right),\\ C_{ikjl}&=\left(S_{ikjl}+A_{ikjl}\right)^{-1}. \end{aligned} $$

(8)

The components of tensor C for the case of an axillary symmetric stress and isolated randomly oriented closed penny-shaped cracks are presented in Appendix 2.

Acoustic wave propagation characterization

Here, we discuss the case when the wavelength λ is much larger than the characteristic size of heterogeneities. In this case, λ∼L, where L is a characteristic size of the geomaterial. Assuming that l (characteristic size of the inhomogeneity of external stress field distribution) is smaller than the wavelength and geological structure sizes, we consider processes in geomaterial at two scales: macroscale characterized by L and λ and microscale characterized by sizes of inhomogeneous distribution of external stress field. Under the assumption of separation of scale condition, we may write α=l/L∼l/λ≪1 introducing new dimensionless small parameter. Then, we are in a position where we can apply two-scale asymptotic homogenization. Suppose that any point is described by two special variables:

• “Slow” coordinate x, giving the general location of the point

• “Fast” coordinate ξ∈Y, the location of the point within the periodic cell

Coordinates x={x ₁,x ₂,x ₃} and ξ={ξ ₁,ξ ₂,ξ ₃} are related by ξ=x/α (α≪1). All the quantities will, in general, be functions of coordinates x and ξ.

Assume that u _i can be expanded in a power series in terms of parameter α:

$$ u_{i}=u_{i}^{(0)}+\alpha u_{i}^{(1)}+\alpha^{2}u_{i}^{(2)}+\ldots,\;\alpha\ll 1 $$

(9)

Assume that, according to the aforementioned conclusions, the expansions (9) are Y-periodic in ξ.

We substitute series (9) into the problem (8) and collect all factors at identical powers of α getting recurrent chain of equations.

At α ⁻²:

$$ \frac{\partial}{\partial \xi_{k}}\left(C_{ikjl}\frac{\partial u_{j}^{(0)}}{\partial \xi_{l}}\right)=0. $$

(10)

At α ⁻¹

$$ \begin{aligned} \frac{\partial}{\partial x_{k}}\left(C_{ikjl}\frac{\partial u_{j}^{(0)}}{\partial \xi_{l}}\right)&+\frac{\partial}{\partial \xi_{k}}\left(C_{ikjl}\frac{\partial u_{j}^{(0)}}{\partial x_{l}}\right)\\ &+\frac{\partial}{\partial \xi_{k}}\left(C_{ikjl}\frac{\partial u_{j}^{(1)}}{\partial \xi_{l}}\right)=0. \end{aligned} $$

(11)

At α ⁰

$$ \begin{aligned} \frac{\partial}{\partial x_{k}}\left(C_{ikjl}\frac{\partial u_{j}^{(0)}}{\partial xi_{l}}\right)&+\frac{\partial}{\partial \xi_{k}}\left(C_{ikjl}\frac{\partial u_{j}^{(1)}}{\partial x_{l}}\right)\\ &+\frac{\partial}{\partial x_{k}}\left(C_{ikjl}\frac{\partial u_{j}^{(1)}}{\partial \xi_{l}}\right)+\frac{\partial}{\partial \xi_{k}}\left(C_{ikjl}\frac{\partial u_{j}^{(2)}}{\partial \xi_{l}}\right)=0. \end{aligned} $$

(12)

Averaging over periodic cell can be defined as

$$ \langle f \rangle =\frac{1}{\vert Y \vert}\int_{Y_{i}} f d\xi, $$

(13)

where |Y| is volume of a periodic cell Y.

Homogenization of the recurrent chain of equations gives us the following results:

Averaging of Eq. (10) in accordance with (13) leads us to the conclusion that the first term $u_{i}^{(0)}$ in a serial expansion does not depend on the variable ξ

$$ u_{i}^{(0)}=u_{i}^{(0)}(x). $$

(14)

Later, $u_{i}^{(0)}$ will be defined as a solution of the homogenized macroscopic equation.

The next term $u_{i}^{(1)}$ in expansion (11) is represented as

$$ u_{i}^{(1)}=N_{ijl}(\xi)\frac{\partial u_{j}^{(0)}}{\partial x_{j}}. $$

(15)

We substitute (15) into (11) and, keeping in mind the independence of $u_{i}^{(0)}$ on ξ, perform averaging of (11) in accordance with (13). As a result, we arrive to the periodic cell problem

$$ \begin{aligned} \frac{\partial}{\partial \xi_{k}}\left(C_{ikjl}+C_{ikqp}\frac{\partial N_{ljq}}{\partial \xi_{p}}\right)=0,\\ \langle N_{ijl}(\xi)\rangle=0. \end{aligned} $$

(16)

In Eq. (16), N _ijl (ξ) are ξ-periodic solutions of the cell problem.

Averaging of Eq. (12) in accordance with (13) gives us the following homogenized macroscopic equations:

$$\begin{array}{*{20}l} -\rho\omega^{2}u_{i}^{(0)}&=C_{ikjl}^{eff}\frac{\partial^{2}u_{j}^{(0)}}{\partial x_{k}\partial x_{l}}, \end{array} $$

(17)

$$\begin{array}{*{20}l} C_{ikjl}^{eff}&=\langle C_{ikjl}\rangle+\left\langle C_{ikqp}\frac{\partial N_{ljq}}{\partial \xi_{p}}\right\rangle. \end{array} $$

(18)

We will be looking for the solution of homogenized problem (17) in the following form

$$ \begin{aligned} u_{i}^{(0)}&=\overline{u}_{i} e^{i\overrightarrow{k}\overrightarrow{x}},\\ \overrightarrow{k}\overrightarrow{x}&=k\left(m_{1}x_{1}+m_{2}x_{2}+m_{3}x_{3}\right), \end{aligned} $$

(19)

where $m_{3}=\cos \theta, m_{1}=\sin \theta \cos \phi, m_{2}=\sin \theta \sin \phi$ (see Fig. 1 for angles ϕ and θ), $\overrightarrow {k}$ is the wave vector, and $k=\left \vert \overrightarrow {k}\right \vert =\frac {2\pi }{\lambda }$ is the wavenumber.

Substituting (19) into (17), we get

$$ a\overline{u}_{i}=C_{ikjl}^{eff}m_{k}m_{l}\overline{u}_{j},\; a=\frac{\rho\omega^{2}}{k^{2}}, \;(i,j,k,l=1,2,3). $$

(20)

$$ \left\| C_{ikjl}^{eff}m_{k}m_{l}-a\delta_{ij}\right\|=0. $$

(21)

System of Eqs. (20)–(21) allows (in general case) to determine the speeds of acoustical waves propagating in the media with closed fractures at arbitrary direction with respect to external stresses. But the given system of equations could be substantially simplified in case of wave propagation in plane perpendicular to σ ₃₃ action direction (see Fig. 2). In this case, the components of displacement vector will depend only on one coordinate (e.g., x), and under the assumption that the media fracturing index value is small $\left (\text {i.e.,}\ \Omega =NR^{3}\ll 1\right)$, the longitudinal wave velocity will be defined by the following expression:

$$ \begin{aligned} \frac{v_{p}}{v_{p0}}&=\left(1-\frac{1-2\nu_{0}}{1-\nu_{0}}\right.\\ &\quad\times\left(\langle A_{1111}\rangle+\frac{\nu_{0}}{1-2\nu_{0}}\left(3\langle A_{1111}\rangle + \langle A_{1133}\rangle\right)\right.\\ &\quad+\left.\left.\left(\frac{\nu_{0}}{1-2\nu_{0}}\right)^{2}2\left(2\langle A_{1111}\rangle+\langle A_{3311}\rangle\right)\right)\right)^{\frac{1}{2}}. \end{aligned} $$

(22)

Fig. 2

Graphical statement of a problem

where v _p0 is the speed of longitudinal wave propagating through elastic medium without cracks

$$ v_{p0}=\sqrt{\frac{E_{0}(1-\nu_{0})}{(1+\nu_{0})(1-2\nu_{0})\rho}}. $$

(23)

In order to investigate the inhomogeneity of external stress field distribution effects on the acoustic wave speed, we considered a case of σ ₃₃ distribution as shown in Fig. 2. Stress σ ₃₃ continuously varied in the interval between σ ^∗−σ ₀ and σ ^∗+σ ₀ with period l.

Figures 3 and 4 demonstrate the results of the longitudinal wave velocities in the media with fractures depending on the inhomogeneity degree for stress distribution.

Fig. 3

Longitudinal wave speed variation with respect to external loading. Case of σ ₃₃ median value variation

Fig. 4

Longitudinal wave speed variation with respect to external loading. Case of σ ₃₃ amplitude value variation

Results on variation of the longitudinal wave speed with respect to changing of values of σ ^∗ and unaltered σ ₀ are shown in Fig. 3. It is obvious that the average value of σ ₃₃ significantly effects the change in the longitudinal wave velocity. At the same time, at the change of the amplitude σ ₀ of variation of σ ₃₃, the longitudinal wave speed changes only in quite narrow range of σ ₁₁.

Conclusions

We have considered the propagation of plane waves in a pre-loaded linear elastic medium weakened by a large number of cracks with periodic function of fracturing degree. The cracks were considered to be closed penny-shaped, isolated, randomly oriented. The model incorporates Coulomb friction between the faces of cracks. A homogenization technique is used to obtain a macroscopic equation for the case of plane wave propagation in effectively elastic media weakened by cracks. We have investigated low-frequency case where wavelength exceeds the characteristic size of heterogeneities. For two special cases of loading (uniaxial compression and complex compression with σ ₃₃<σ ₁₁=σ ₂₂≤0), we have characterized (theoretically) the effective elastic properties of the media. In terms of plane wave propagation through the fractured media, we mostly concentrated at speeds of longitudinal waves. We have studied how wave speeds depend on averaged concentration and distribution of cracks and changing of external load.

Appendix 1: evaluating components of tensor A

As it was shown by (Talonov and Tulinov 1988), nonzero components of tensor A can be written as

$$ \begin{aligned} A_{1111}&=A_{2222}=-\frac{D}{2}\left(f_{2}(\alpha_{2})-f_{2}(\alpha_{1})\right),\\ A_{3333}&=D\left(f_{1}(\alpha_{2})-f_{1}(\alpha_{1})\right),\\ A_{1133}&=A_{22332}=-\frac{D}{2}\left(f_{1}(\alpha_{2})-f_{1}(\alpha_{1})\right),\\ A_{3311}&=A_{3322}=\frac{D}{2}\left(f_{2}(\alpha_{2})-f_{2}(\alpha_{1})\right), \end{aligned} $$

(24)

where

$$ \,D=\frac{32(1-\nu_{0}^{2})\Omega}{3E_{0}(2-\nu_{0})},\;\Omega=NR^{3}, $$

(25)

$$ \begin{aligned} f_{1}(\psi)&=\frac{\cos^{5}\psi}{5}-\frac{\cos^{3}\psi}{3}-\mu\left(\frac{\sin^{3}\psi}{3}-\frac{\sin^{5}\psi}{5}\right),\\ f_{2}(\psi)&=\frac{\cos^{3}\psi}{3}-\frac{\cos^{5}\psi}{5}-\mu\frac{\sin^{5}\psi}{5}. \end{aligned} $$

(26)

There are two special cases of loading that specify corresponding values of α ₁ and α ₂:

• Case 1: uniaxial compression with σ ₃₃≠0, σ ₁₁=σ ₂₂=0. In this case, α ₁ and α ₂ do not depend on stress and can be evaluated as (see Talonov and Tulinov 1988)

  $$ \alpha_{1}=\frac{\pi}{2},\;\alpha_{2}=\beta=\arctan\mu. $$

  (27)

• Case 2: a complex stress with σ ₃₃≠0, σ ₁₁=σ ₂₂≠0. In this case, α ₁ and α ₂ depend on stress and can be written as (see (Talonov and Tulinov 1988))

  $$ \alpha_{1}=\frac{\beta}{2}+\frac{1}{2}\arcsin\left(\sin \beta\cdot\frac{\gamma+1}{\gamma-1}\right), $$

  (28)
  $$ \alpha_{2}=\frac{\pi}{2}+\beta-\alpha_{1}=\frac{\pi}{2}+\frac{\beta}{2}-\frac{1}{2}\arcsin\left(\sin \beta\cdot\frac{\gamma+1}{\gamma-1}\right), $$

  (29)

  where

  $$ \beta=\arctan \mu,\; \gamma=\frac{\sigma_{33}}{\sigma_{11}},\;\gamma>\frac{1+\sin\beta}{1-\sin\beta}. $$

  (30)

Appendix 2: nonzero dimensionless averaged components of tensor C

$$\begin{array}{@{}rcl@{}} \begin{aligned} C_{1122}&=C_{2211}=\frac{E_{0}}{(1-\nu_{0})B}\left((1+A_{3333})(\nu_{0}-A_{1111})\right.\\ &\quad+\left.(\nu_{0}-A_{1133})(\nu_{0}-A_{3311})\right),\\ C_{1111}&=C_{2222}=\frac{E_{0}}{(1-\nu_{0})B}\left((1+A_{3333})(1+A_{1111})\right.\\ &\quad-\left.(\nu_{0}-A_{1133})(\nu_{0}-A_{3311})\right),\\ C_{1133}&=C_{2233}=\frac{E_{0}(\nu_{0}-A_{1133}}{(1-\nu_{0})B},\\ C_{3333}&=\frac{E_{0}}{1+A_{3333}}\left(1+\frac{2(\nu_{0}-A_{3311})}{B}\right),\\ C_{3322}&=\frac{E_{0}(\nu_{0}-A_{1133})}{(1+A_{3333})(1+\nu_{0})}\left(1+\frac{2}{B}\left((1+A_{3333})(\nu_{0}-A_{1111})\right.\right.\\ &\quad+\left.\left.(\nu_{0}-A_{1133})(\nu_{0}-A_{3311})\right)\right),\\ C_{3311}&=C_{3322},\\ C_{ijji}&=C_{ijij}=\frac{E_{0}}{2(1+\nu_{0})},\;i,j=1,2,3,\\ B&=(1-\nu_{0}+2A_{1111})(1+A_{3333})\\ &\quad-2(\nu_{0}-A_{1133})(\nu_{0}-A_{3311}). \end{aligned} \end{array} $$