Abstract
This paper investigates the wave transmission and reflection of an elastic P-wave at a single joint for normal incidence. First, considering a coupled joint (correction parameter λ, ), a normal deformation constitutive model of the joint ( model) under static or quasi-static loading is introduced and then extended to dynamic loading. The nonlinearity of the joint stress-deformation curve increases with increasing λ. Second, the interaction between the P-wave and the joint is investigated by using the method of characteristics and the displacement discontinuity method to deduce the differential expression of the transmitted wave’s particle velocity. The approximate analytical expressions of the transmission and reflection coefficients are obtained according to the Lemaitre equivalent strain assumption. Third, parametric studies are conducted to evaluate the effects of λ on transmission characteristics for a normally incident P-wave at a single joint. The results show that the particle velocity of the transmitted wave depends on λ. When λ takes the limit values 0 and 1, the transmitted wave’s particle velocities are then consistent with the conclusions of the classical exponential model and the Barton–Bandis model. In addition, the transmission and reflection coefficients are discussed with respect to λ and also to the ratio of the joint closure to the maximum allowable joint closure.
1. Introduction
In Earth’s crust, discontinuities, such as faults, joints, bedding planes, and fractures, are ubiquitous and they often control the hydraulic and mechanical behavior of rock masses [1–4]. Elastic waves propagating through jointed rock masses will attenuate (or slow down) because of the existence of abundant rock joints [5–9]. The properties of the jointed rock mass including the joint matching coefficient (JMC), joint roughness coefficient (JRC), and the joint stiffness affect not only the dynamic behavior of joint but also the stress wave energy attenuation confirmed by the laboratory experiment (i.e., split Hopkinson pressure bar (SHPB) and ultrasonic test) [10–12]. Other factors, such as rock acoustic impedance and rock thermal effect, also affect the stress wave propagation across the jointed rock mass [13, 14]. Due to these facts, it is necessary to further study the interaction between the P-wave and the joint to prove serving the application of the elastic waves in rock engineering.
Different constitutive models have been presented to investigate the deformation characteristics of a jointed rock mass under static or quasi-static loading/unloading. Amongst them, the Barton–Bandis (BB) model, a modification of Goodman’s hyperbolic model [15, 16], is widely used to describe the nonlinear deformation behavior of a single joint [17, 18]. Although the BB model is established under the quasi-static or cyclic loading condition, it is reasonable for the dynamic condition because jointed rock masses have experienced multiple deformations in geological history [19]. Many researchers have investigated the P-wave propagating through the joint based on the BB model [5, 20–22]. However, since the BB model is based on experimental data of five kinds of jointed rock masses (slate, dolerite, limestone, siltstone, and sandstone), it may not be universally applicable to all other jointed rock masses. Some improved constitutive models for jointed rock masses are proposed or deduced from different experimental data or mathematical principles. For instance, the matrix of compliance components for a rock joint was adopted to describe the stress-deformation relations under compression and shear loads [23]. The relationship between constant normal stress and time-dependent closure was formulated as a function of time in terms of the aperture distribution of a fracture and the relaxation modulus of rock [24]. With introducing a half-closed stress and correction parameters, Malama and Kulatilake [25] presented a generalized semiempirical exponential model, which can better fit the results of fracture closure experiments of jointed diorite and granodiorite rock specimens under monotonically increasing normal compressive loading. Yu et al. [26] proposed a 3-parameter constitutive model, which overcame the mathematical defects of the BB model and the classical exponential model. Following Yu’s viewpoint, Rong et al. [27] established the model () and the model () under static or quasi-static loading. The results indicated that the model can be more suitable for the coupled () and uncoupled () joint. However, the model has been made up for the mathematical defects and been verified by different experimental data. The collection parameter λ can be adjusted at the medium stress level, and the stress-deformation curve of the joint was more applicable than the classical exponential model and the BB model. The classical exponential model and the BB model have been extended to dynamic loading [19, 20]. The model’s dynamic characteristic is yet neglected.
In this paper, considering a coupled joint (), a normal deformation constitutive model of the joint ( model) under static or quasi-static loading is introduced first and then extended to dynamic loading. Second, the finite difference formula suitable for transmitted wave’s particle velocity is obtained found on the method of characteristics (MC) [5, 19, 22, 28] and the displacement discontinuity method (DDM) [29–33]. The approximate analytic solutions of the transmission and reflection coefficients are deduced theoretically based on damage mechanics. Finally, parametric studies and discussions of a P-wave transmission across a single joint which is regarded as a nonwelded interface are conducted according to the MATLAB numerical program. The results show that the parameter λ has a significant effect on the dynamic characteristics of the reflected and transmitted waves.
2. Theoretical Formulations
2.1. Dynamic Expansion of the Model
The model is a normal deformation constitutive model of a single joint and can be applied to both the coupled and uncoupled joint [27]. One of the key features of the model is adding the correction parameter λ (dimensionless) to speed up the normal deformation. The model is formulated as follows:where and denote the joint closure and the maximum allowable closure of the joint, respectively, is the joint normal stiffness at initial stress, and is the normal effective stress. The values of and can be determined from the laboratory measurements on JRC, joint surface compressive strength (JCS), and average aperture thickness () [17, 18]. The parameter λ is in association with joint weathering, roughness, fluctuation degree, matching of the joint surface, and strength of the wall of rock joints. In this paper, we only consider the case of for the coupled joint. The joint is assumed to be planar, large in extent, and small in thickness compared to the wavelength, and the joint and the intact rock are linearly elastic. If , according to the L’Hospital’s rule, the exponential term in equation (1) can be obtained as follows:to giveThat is, the classical exponential model is recovered for . On the contrary, when , the exponential term in equation (1) becomes
Substituting equation (4) in (1), we getand thus, the BB model is a special case of the model.
The schematic curves of the model with different values of λ, the BB model, and the classical exponential model of normal deformational behavior of a single joint are shown in Figure 1. The tangent slope of each curve can be defined as the joint stiffness [19]. Compared to the linear model, the model, the BB model, and the classical exponential model involve the joint stiffness that increases with the increasing normal effective stress. When the normal stress increases nearly to an infinite value (), the joint closure of each model approaches the maximum allowable closure (). The curve of the model is between the BB model and the classical exponential model. In addition, the tangent slope of the model’s curve increases with the increase of λ, which means that the degree of nonlinearity increases.
In the meantime, the model has been validated to be more suitable than other models by fitting different experimental results. The normal deformation of the BB model is less than test data, and the normal deformation of the classical exponential model is greater than it under the intermediate stress level, as shown in Figure 2. Under the high stress level, the approach speed of normal deformation in the BB model to the maximum joint closure is lower than that of the test data. The model and the model are more consistent with experimental data. Of them, the model fits best, where [27].
A large body of literature has confirmed that continuous cyclic loading and unloading can make the stress-deformation relationship become the elastic relationship without a hysteresis loop because natural rock joints undergo many deformations in the long geological history [5, 19, 20, 34]. Therefore, the hysteresis phenomenon in the initial cyclic loading and unloading process and the influence of the loading rate on the joint deformation can be neglected. It is reasonable for the application of the BB model and the classical exponential model in the wave propagation. This viewpoint is borrowed in the present study to extend the model to the dynamic condition for computational analysis.
2.2. Method of Characteristics
It is assumed that a joint exists at in a half-space of a linear elastic, homogeneous and isotropic medium. When a normal incidence plane P-wave impinges upon the boundary of the half-space, a reflected wave and a transmitted wave will be generated at the joint. Through the linear elastic DDM which has been applied to wave propagation, the stresses on the two sides of the joint are assumed to be continuous while the difference in displacements equals the closure (opening) of joint, which can be, respectively, expressed as follows:where and are the normal stress before and after the P-wave through the joint, respectively. Based on equation (1), the displacement difference between the two sides of the joint iswhere and are the displacement before and after the P-wave across the joint, respectively.
The derivative expression of equation (7) is as follows:where and are the particle velocity before and after the P-wave through the joint, respectively.
The method of characteristics has been widely used for solving the problems of one-dimensional wave propagation in continuous linearly elastic media, as shown in Figure 3. The relations of the particle velocity and the stress σ can be derived aswhere z is the rock acoustic impedance and equal to the product of the rock density ρ and the rock P-wave velocity .
Since the half-space is undisturbed at t = 0, both the particle velocity and the stress are zero at every point of the x-axis. Therefore, along left-running characteristic line ab,
Similarly, along the right-running characteristic line ac,where is the particle velocity input to the boundary at time and is the stress for = 0 at time .
Along another left-running characteristic line cd,
Summing up equations (11) and (12) yields
Based on equations (10) and (13), the relation between the particle velocity before and after the joint can be expressed as follows:
Substituting equations (8) and (10) into (14), we can derive as follows:
Similarly, if , on repeating the above derivation process we can obtain as follows:When , the expression of is consistent with Zhao and Cai [19] as follows:
For numerical calculation, equation (15) is expressed in differential form aswhere is the time interval. Equation (18) is an iterative calculating equation for calculating . Following the same viewpoint, equations (16) and (17) also have a corresponding differential iterative formula, respectively. The energy of the transmitted wave is examined by the energy transmission coefficient :where and are the period of transmitted and incident waves, respectively; and are the initial time of transmitted and incident waves, respectively; and are the particle velocity of transmitted and incident waves at time and position , respectively.
2.3. Transmission Coefficient and Reflection Coefficient
According to the analytical solutions obtained by Schoenberg [29] and Pyrak-Nolte et al. [30], the transmission coefficient and the reflection coefficient of a linear deformation joint in identical rock can be calculated as follows (for compressional/shear wave incidence):where ω is the angular frequency. For compressional wave incidence, k denotes the joint normal stiffness and z has the same definition as in equation (9):where ω, k, and z have the same definition as in equation (20).
In the model, the flexibility of joint iswhere is the equivalent stiffness [27].
In order to obtain the transmission and reflection coefficients for the model, here, we draw on the basic idea of the Lemaitre equivalent strain assumption in damage mechanics [35] and assume to represent the nonlinear coefficient of joint stiffness. The joint normal stiffness k in equations (20) and (21) can be replaced by the equivalent stiffness in equation (22). Therefore, the approximate analytical solutions of the transmission coefficient and the reflection coefficient of the model can be represented as follows:where is the ratio of the joint closure to the maximum allowable joint closure. It is found that the obtained and are in agreement with the corresponding results when γ is equal to 0 [29, 30]. Note that although the transmission and reflection coefficients of a harmonic wave are conducted in the frequency domain, the technique of the inverse fast Fourier transform (IFFT) can be used for an incident wave with the arbitrary waveform to the time domain [36–38].
3. Parametric Study and Discussions
To demonstrate that the model can be used for wave propagation, different parameters are discussed by the MATLAB numerical program. We mainly analyze the wave transmission and reflection of a P-wave travelling through a nonwelded single joint at normal incidence, including (i) the behavior of the transmitted wave, (ii) the amplitude-dependence and frequency-dependence of the energy transmission, and (iii) the relationship between the transmission/reflection coefficient and γ.
3.1. Particle Velocity of the Transmitted Wave
The threshold values of wave amplitude, such as the peak particle velocity (PPV), are essential to regulate the damage criteria of rock structures in earthquake engineering. The dynamic compressive strength of rock masses is generally much greater than the dynamic tensile strength. The selection of a wavelet does not cause damage to the jointed rock mass. In other words, the dynamic stress generated by the P-wave should be less than the dynamic tensile strength of the jointed rock mass. Using the fast Fourier transform (FFT) and IFFT, any arbitrary incident wave can be expressed as the sum of a series of harmonic waves [36, 37]. Without losing generality, the incident wave is assumed to be a positive part of one-cycle sinusoidal wave as follows:where is the particle velocity of the incident wave at time t and position = 0 and and f denote the amplitude and the frequency of the incident wave, respectively [21, 39]. The change of the particle velocity for the transmitted wave is analyzed when the half-sine incident wave propagates through a single joint with four different cases of λ, including (i) , (ii) , (iii) , and (iv) . The key parameters of a jointed rock mass are listed in Table 1.
When an incident wave reaches a joint in the half-space, a reflected wave and a transmitted wave are created and propagate in two opposite directions. According to equation (18) and the difference forms of equations (16) and (17), Figure 4 shows the particle velocities of the transmitted wave on the basis of the model, the BB model, and the classical exponential model ( = 0.1 m/s and f = 50 Hz). It is observed that the waveform of the transmitted wave is similar to that of the incident wave, but the phase has a delay. These phenomena are also observed in [19, 30]. The maximum particle velocities of the transmitted wave are 0.0874 m/s, 0.0911 m/s, 0.0937 m/s, and 0.0947 m/s, as , λ = 0.4, λ = 0.8, and , respectively. The amplitude of the transmitted wave’s particle velocity increases gradually with the increase of λ. If , the transmitted wave for the model is closer to the value of the classical exponential model, and if , it is closer to that of the BB model. These two special cases calculated under the dynamic condition are consistent with those of the static condition. It shows that the extension of the model from the static loading to the dynamic loading is workable. As shown in Figure 4, more waves are transmitted with the increase of λ (the weaker the degree of weathering, the less the roughness and joint closure) during the transmission of a P-wave pulse. In other words, the increase of the parameter λ will increase the joint stiffness, resulting in more incident waves being transmitted. The values of the joint closure calculated by equation (23) corresponding to the amplitude of the transmitted wave are listed in Table 2. It can be seen that the joint closure decreases as the parameter λ increases, further indicating that the joint stiffness is gradually increasing.
3.2. Amplitude Dependence and Frequency Dependence of Energy Transmission
In rock engineering, the energy transfer is often concerned when elastic waves propagate through a jointed rock mass. For example, the blast-induced stress waves propagate into the jointed rock mass. Thus, we investigate the properties of amplitude-dependence and frequency-dependence of energy transmission with an incident wave in a form as given in equation (25).
The particle velocity amplitude is set to a range of 0.05 m/s to 0.5 m/s. The corresponding stress amplitude () of the incident wave is from 0.54 MPa to 5.4 MPa. The frequency f is still 50 Hz. According to equation (19), the variation of the energy transmission coefficient with respect to the amplitude of the incident wave is shown in Figure 5. If the parameter λ is constant, increases with the increase of the incident wave’s amplitude. also increases with the increase of λ when the incident wave’s amplitude is fixed. Due to the nonlinear behavior of the joint, the specific stiffness increases with the amplitude of the incident wave which leads to a smaller difference between the joint stiffness and intact rock and more energy transmission.
With different values of λ, the relations of the energy transmission coefficient versus the frequency of the incident wave are shown in Figure 6, where the amplitude is 0.1 m/s and the frequency f changes from 50 Hz to 1000 Hz. decreases with the increase of the incident wave’s frequency. These phenomena indicate that the high-frequency components of the incident wave are filtered out, and only the low-frequency waves pass through. Many literatures have reached this conclusion [19, 21, 30, 32]. The parameter λ has little effect on the energy transmission when the frequency is constant.
3.3. Numerical Calculation of Transmission and Reflection Coefficients
According to equation (23), the transmission coefficient for four different values of λ increases until 1 as the parameter γ (the ratio of the joint closure to the maximum allowable joint closure) increases until 1, 0.97, 0.94, and 0.92, respectively (Figure 7(a)). If γ is fixed, also increases with the increase of parameter λ. It not only indicates that the joint with a smaller aperture generates more wave transmission across the joint but also that the rapid closure of the joint makes the joint stiffness increase faster. The overall growth trends of in Figures 7(a) and 7(b) are similar. However, the curves of at low frequency display an upward convex trend, while at high frequency, the S-shapes appear. at low frequency is greater than that at high frequency when λ and γ are the same.
(a)
(b)
Similarly, the reflection coefficient for four different values of λ gradually decreases until 0 with the increase of γ by using equation (24), as shown in Figure 8. For a given γ, decreases with increasing λ. From the point of view of energy conservation and without considering energy loss, the energy of the transmitted wave passing through the joint increases, and the energy of the reflected wave decreases, and vice versa. The difference between Figures 8(a) and 8(b) is that at low frequency is less than the high frequency when λ and γ are the same. As λ increases, the degree of the joint weathering is weaker and the amount of the joint closure is less, which leads to more transmitted waves and less reflected waves.
(a)
(b)
In general, the transmission coefficient of the nonlinear joint model (represented by ) is defined as the ratio of the transmitted wave’s peak particle velocity to the incident wave’s amplitude [10, 19, 36, 38]. However, the approximate analytical solution of the transmission coefficient in the present study is based on the viewpoint of damage mechanics, and the equivalent stiffness is used to replace the stiffness constants in equation (20). The difference between the two methods is compared using the combination of parameters listed in Table 3. It can be observed that and have a slight deviation, as shown in Figure 9. The deviation gradually decreases as the joint stiffness increases and maximum allowable joint closure decreases. This phenomenon indicates that the joint stiffness and the joint closure have a significant effect on the transmission coefficient calculated by the two methods. Therefore, if the joint stiffness is large and the amount of joint closure is little, it is a simple and convenient method to deduce the transmission coefficient and the reflection coefficient , ignoring the lengthy mathematical derivation process.
4. Conclusions
The model (a normal deformation constitutive model of the joint) is adopted to analyze a P-wave propagating through a nonwelded single joint. The following main conclusions are drawn:(1)The maximum amplitude of the transmitted wave increases gradually with the increase of parameter λ when using a half-sine wave as the incident wave. The waveform of the transmitted wave is like the waveform of the incident wave, but the phase has a delay. If , the transmitted wave is closer to the conclusion of the classical exponential model, while , it is in agreement with that of the BB model. It shows that the extension of the model from the static loading to the dynamic loading is feasible.(2)The energy transmission increases with the increase of the incident wave’s amplitude and decreases with the increase of the incident wave’s frequency. It also increases with the increase of λ when the amplitude is fixed. The filtering effect of the joint is further confirmed. However, parameter λ has little effect on the energy transmission if the frequency is constant.(3)The magnitude of the transmission and reflection coefficients depends on the joint stiffness and the extent of nonlinearity of rock joint deformation behavior. The transmission coefficient increases and the reflection coefficient decreases with the increase of λ. Besides, two methods calculating the transmission coefficient have less deviation when the joint stiffness increases and the maximum allowable closure of joint decreases.
Further studies are required to verify the accuracy of the model with different types of rock loading tests and nonlinear deformation behavior of joint as the P-wave obliquely propagates through a jointed rock mass.
Data Availability
The data used to support the findings of this study are available from the corresponding author upon request.
Conflicts of Interest
The authors declare that there are no conflicts of interest.
Acknowledgments
This work was supported by the National Key Research and Development Program of China (2016YFB1200401) and the Western Construction Project of the Ministry of Transport (2015318J29040).