Numerical Analysis of Cyclic Impact Damage Evolution of Rock Materials under Confining Pressure

Abstract

To study the effect of cyclic impacts and confining pressure on the damage evolution of rock materials, numerical simulations of cyclic impact tests on rock materials under confining pressure were carried out by LS-DYNA using dynamic relaxation and full restart analysis. The static confining pressure was applied by dynamic relaxation, and cyclic impacts were realized by full restart analysis. As the crack generation and propagation result in the failure of elements in the finite element model, the damage variable defined by the crack density method was characterized by volume reduction. Numerical simulation indicates that both the confining pressure and amplitude of incident stress waves significantly affect the damage evolution of rock materials. High incident stress waves lead to severe damage, while large confining pressure results in minor damage. Under confining pressure, the damage to rock materials is alleviated due to the constraint effect on crack propagation. The number of cyclic impacts before macroscopic fracture increases as the confining pressure increases and decreases when the amplitude of the incident stress waves increases. The cumulative damage of rock materials under confining pressure progressively increases with the number of cyclic impacts, and the damage evolution exhibits three distinct stages: rapid rising, steady development, and sharp rising.

1. Introduction

In recent years, due to the rapid economic development and increasing demand for mineral resources, rock engineering involved in tunneling, mining resources exploitation, and underground space utilization has extended into deep environments [1,2,3]. Moreover, during rock excavation in deep environments, the dynamic disturbance induced by blasting or rock burst has a significant influence on the stability of rock engineering [4,5,6]. In fact, the instability of surrounding rock in deep rock engineering is not typically attributed to a single impact event but cumulative damage induced by cyclic loads [7,8,9]. Cumulative damage induced by cyclic loads leads to crack generation and propagation, fracture, and eventually instability [10,11,12,13]. Therefore, investigating the dynamic damage evolution of rock materials under confining pressure holds significant engineering application value and is crucial for the protection of underground projects.

The damage evolution of rock materials under cyclic impact conditions has been investigated by analyzing elastic modulus, acoustic emission, energy dissipation, wave impedance, and other parameters. Lu [14] and Mei [15] identified three stages in the damage evolution of rocks under cyclic impact, namely the initial stage, low-speed stage, and accelerated stage. Rock specimens experience significant impact damage at the initial stage, then the cumulative damage gradually increases at the low-speed stage, and finally, the damage accelerates rapidly, resulting in rock failure. Based on the results of cyclic impact tests, Jin et al. [16,17] employed wave impedance to define the damage variable of rocks and revealed that the wave impedance method provides a more accurate depiction of the damage evolution characteristics of rocks than other commonly used methods for defining damage variables. Gan et al. [18] performed cyclic impact tests on sandstone using split Hopkinson pressure bar apparatus and analyzed the damage evolution based on the Weibull distribution statistical damage model. Liu et al. [19] conducted cyclic impact tests on sandstone specimens under various confining pressures and observed that the confining pressure restricts the accumulation of damage, namely a significant confining pressure effect. In the uniaxial cyclic impact test of rock, Wang et al. [20] employed the static crack initiation stress in a static stress-strain curve to estimate the dynamic crack initiation stress. By conducting constant amplitude cyclic impacts on granite, Tian et al. [21] found that the number of cyclic impacts initially increases and then decreases with the rise in axial pressure, which implies a strengthening effect when the pre-axial pressure is lower than the rock’s crack initiation stress. Under the combined action of pre-axial compression and cyclic impact, Tang et al. [22] observed that the stress-strain curve of rocks under various axial pressures exhibits two distinct features, namely rebound and non-rebound, and the fragmentation of rocks is intensified as the axial pressure rises.

Lv et al. [23] investigated the dynamic mechanical properties of rocks subjected to cyclic impacts under confining pressure and observed that the confining pressure can enhance the rock’s resistance to cyclic loading. Yu et al. [24] performed cyclic impact tests on rocks using a SHPB apparatus with a confining pressure device and discovered a power function relationship between rock damage variables and average strain rate. Through three-dimensional coupled static-cyclic impact tests, Jin et al. [25,26] found that axial pressure and confining pressure are negatively correlated and positively correlated with the total number of cyclic impacts, respectively. Therefore, reducing the axial pressure or increasing the confining pressure can improve the impact resistance of rock. By analyzing the energy dissipation and damage of rock under cyclic impact, Guo et al. [27] found that confining pressure can greatly increase the ability of basalt to resist external impact, and the cumulative specific energy absorption value at failure is 10 times higher than that without confining pressure.

Although numerous scholars have investigated the dynamic mechanical properties and damage characteristics of rocks under cyclic impacts, few scholars have utilized numerical simulations to analyze the cyclic impact damage evolution of rock under confining pressure. Therefore, numerical simulations of cyclic impact tests based on the SHPB apparatus were carried out for rock under confining pressure by LS-DYNA. Then the damage evolution was investigated by the volume change of the rock specimen after cyclic impacts, and the damage evolution model was established using crack density as the damage variable. The objective of this study is to provide a new method for analyzing the cyclic impact damage evolution of rock materials under confining pressure.

2. Numerical Modeling and Verification

2.1. SHPB Fundamentals

SHPB (split Hopkinson pressure bar) apparatus is a widely employed equipment to investigate the dynamic mechanical properties of rock materials. A SHPB apparatus comprises an impact bar, an incident bar, a transmission bar, and an energy absorption device. During a SHPB test, the rock specimen is sandwiched between the incident bar and transmission bar. Axial pressure and confining pressure loading systems are incorporated with the SHPB apparatus to achieve combined dynamic and static loading.

By collecting incident wave ε_i(t) and reflected wave ε_r(t) on the incident bar and transmitted wave ε_t(t) on the transmission bar, the stress-strain relationship of the rock specimen can be determined according to the classical two-wave method.

$$\epsilon \left(t\right)=\frac{c_0}{l_0}{\int }_0^t{[\epsilon }_i\left(t\right)-{\epsilon }_r\left(t\right)-{\epsilon }_t\left(t\right)]dt=-\frac{2c_0}{l_0}{\int }_0^t{\epsilon }_r\left(t\right)dt$$

(1)

$$\sigma \left(t\right)=E{\epsilon }_t\left(t\right)$$

(2)

where E and c₀ are the elastic modulus and wave velocity of steel pressure bars, respectively, and l₀ is the length of the rock specimen.

During cyclic impact tests, a half-sine incident stress wave can effectively maintain a constant loading strain rate and ensure stress uniformity within the specimen during deformation [28,29]. In numerical simulations, the half-sine incident stress wave can be achieved by using initial velocity through a spindle-type impact bar or applying half-sine stress waves at the front end of the incident bar [30]. In this study, the half-sine incident wave is directly applied at the front end of the incident bar. The fitting formula of the half-sine incident stress wave is shown in Equation (3).

$$P_d\left(t\right)=\left\{\begin{array}{ll}{P}_{max}\left[\frac{1}{2}-\frac{1}{2}\cos (2\pi \omega t)\right]& \left(t\le 1/\omega \right)\\ 0& \left(t>1/\omega \right)\end{array}\right.$$

(3)

where P_max, ω and t are amplitude, angular frequency, and duration of the half-sine incident stress wave, respectively.

In this study, considering the physical and mechanical properties of rock, the amplitude of half-sine incident stress wave is taken as 170 MPa, 180 MPa, and 190 MPa.

2.2. Numerical Modeling

Based on the actual size of a SHPB apparatus with a diameter of 50 mm, the finite element model of the incident bar, transmission bar, and rock specimen is established using the solid element in Hypermesh software. In the numerical model, both the incident bar and transmission bar have a diameter of 50 mm and a length of 2 m, and the rock specimen has a diameter of 50 mm and a length of 25 mm. To accurately represent the damaging effect and damage evolution of rock under cyclic impacts, the rock specimen is meshed into 20 grids along the axial direction and 30 grids along the radial direction, as shown in Figure 1.

During the cyclic impact test, both the incident bar and transmission bar are in a linear elastic state. Therefore, the isotropic linear elastic model, namely MAT_ELASTIC, is applied for both the incident bar and transmission bar. Considering the damage evolution induced by cyclic impacts, a dynamic damage constitutive model accounting for strain rate effects, namely Holmquist–Johnson–Cook (HJC) model, is applied to rock material [31,32]. The HJC model can effectively simulate the dynamic response of rock-like materials under high strain rates, large deformation, and high confining pressure [33].

Table 1. Physical and static mechanical parameters of steel pressure bar.

ρ/(g·cm−3)	E/GPa	μ
7.80	210	0.3

and

Table 2. Physical and static mechanical parameters of rock material.

ρ/(g·cm−3)	E/GPa	μ	fc/MPa	G/GPa	K/GPa
2.42	12.9	0.25	88.1	5.16	8.6

give the basic physical and mechanical parameters of the steel pressure bar and rock material. In Table 1 and Table 2, ρ is density, μ is Poisson’s ratio, E is elastic modulus, f_c is static uniaxial compressive strength, G is shear modulus, and K is bulk modulus.

Keyword MAT_ADD_EROSION defines the failure threshold and controls the damage incurred by the material, and the maximum failure principal strain is set as 0.006. To make stress wave propagation closer to reality, the contact between the steel pressure bar and the rock specimen is set as an erosion contact, namely CONTACT_ERODING_SURFACE_TO_SURFACE, with the rock face being dominant. The remaining parameters in both contact types are set as the default values.

2.3. Cyclic Impact Simulation of Rock Materials under Confining Pressure

As illustrated in Figure 2, cyclic impact tests for rock material under confining pressure are implemented by applying static circumferential confining pressure on rock specimens before impact load. In Figure 2, P_c and P_d represent the static circumferential confining pressure and dynamic impact load.

LS-DYNA offers two methods to apply static circumferential confining pressure on rock specimens before cyclic impact, namely the dyain file method and the dynamic relaxation method [34]. But the dyain file method cannot inherit the damage of rock specimens in cyclic impact simulations. Therefore, the dynamic relaxation method is employed for a full restart analysis to apply static circumferential confining pressure on rock specimens. The convergence of dynamic relaxation can be assessed through the deformation movement energy. When the dynamic relaxation process ends, the prestressed state is transferred to the initial conditions for transient analysis, allowing for a full restart. The key steps of cyclic impact simulation of rock materials under confining pressure are listed as follows.

(1)

Segment setting. The lateral surface of the cylindrical rock specimen and the end surface of the transmission bar are defined as different groups by using the keyword SEGMENT_SET.

(2)

Boundary setting. Firstly, constrain all translation and rotation of the incident bar and transmission bar except the horizontal motion by using the keyword SPC_SET. Secondly, set the non-reflective boundary at the end surface of the transmission bar by using the keyword NON_REFLECTING.

(3)

Confining pressure applying. Firstly, define a preload curve according to static circumferential confining pressure by using the keyword DEFINE_CURVE. Secondly, apply the preload curve on the lateral surface of the cylindrical rock specimen by using the keyword LODING_ SEGMENT_SET. Thirdly, run the dynamic relaxation analysis by keyword CONTROL_DYNAMIC_RELAXATION to make the static circumferential confining pressure maintain the designed value before impact load.

(4)

Cyclic impacts realization. Cyclic impact simulation is realized by using a full restart analysis.

Cyclic impact simulations are conducted on rock specimens under confining pressures of 0 MPa, 2 MPa, and 4 MPa, respectively. The amplitude of half-sine incident stress for cyclic impact simulations is 170 MPa, 180 MPa, and 190 MPa under 0 MPa confining pressure, while it is only 180 MPa under 2 MPa and 4 MPa confining pressure. The objective is to investigate the effect of the amplitude of incident stress waves and confining pressure on the damage degree and damage evolution of rock materials.

2.4. Verification of Numerical Model

SHPB apparatus relies on two basic assumptions, namely, one-dimensional stress wave assumption and stress uniformity assumption. It is crucial to achieve stress uniformity in rock specimens to ensure the validity of numerical simulation. To satisfy the stress uniformity assumption, the dynamic stresses at both ends of the rock specimen must be approximately equal [35,36]. Therefore, stress uniformity can be verified by comparing the stress-time histories at both ends of rock specimens during cyclic impacts. As illustrated in Figure 3 and Figure 4, stress balance verifications are conducted using the initial impact in different conditions. In Figure 3 and Figure 4, In, Re, and Tr represent incident stress waves, reflected stress waves, and transmitted stress waves, respectively. It can be observed that the sum of the incident stress wave and reflected stress wave aligns with the transmitted stress wave prior to reaching the peak stress, indicating the achievement of stress balance during the impact.

3. Results and Analysis

3.1. Dynamic Stress-Strain Curve Analysis

When the amplitude of half-sine incident stress waves is 180 MPa, the dynamic stress-strain curves of rock specimens after cyclic impacts under different confining pressures are presented in Figure 5. In Figure 5, n represents the number of cyclic impacts.

As shown in Figure 5, dynamic stress-strain curves exhibit a relatively similar trend, which can be roughly divided into two stages, namely the initial nearly straight loading stage and the nonlinear unloading stage with a gradually decreasing slope. When rock specimens remain integrity, the partial rebound at the end of the unloading stage indicates the release of stored elastic energy.

When the confining pressure is constant, with the number of cyclic impacts increasing, both the slope in the straight loading stage, namely elastic modulus, and peak stress of the dynamic stress-strain curve gradually decreases. Under cyclic impacts, internal crack initiation and propagation after each impact load lead to the decrease of mechanical properties of rock specimens. Consequently, the cumulative crack development of cyclic impacts intensifies the damage to the internal structure of the rock. The larger the number of cyclic impacts, the greater the damage, leading to a significant reduction in mechanical properties. With the increase of confining pressure, the number of cyclic impacts increases, which implies that the confining pressure can slow down the crack propagation and reduce damage accumulation.

3.2. Numerical Simulation of Rock Fracture

When the amplitude of half-sine incident stress waves is 180 MPa, the rock fracture process after cyclic impacts under 0 MPa confining pressure and 2 MPa confining pressure are illustrated in Figure 6 and Figure 7, respectively. In numerical simulation, the crack generation and propagation result in the failure of elements.

As shown in Figure 6, the rock failure process exhibits a fracture surface approximately along the axial direction without the end effect. At the initial stage of the fracture process, a crack first appears on the periphery of the rock specimen, as shown in Figure 6b. Under the impact loads, the rock specimen experiences axial compression and radial tension. When the transverse tensile strain exceeds the criterion, short cracks generate in the middle of the outer side of the rock specimen, as depicted in Figure 6c. As the transverse tension deformation continues, the short cracks extend parallel to the length of the rock specimen until penetrating the entire specimen. Eventually, strip-shaped blocks are observed at the periphery of the rock specimen, as depicted in Figure 6d. Under cyclic impacts, the fracture of the rock specimen extends from the periphery toward the interior.

As illustrated in Figure 7, the fracture on the rock specimen increases with the number of cyclic impacts. In the early stage of the fracture process, a crack also first generates on the periphery of the rock specimen, as shown in Figure 7b. Due to the constraint effect of confining pressure on crack propagation, the number of cyclic impacts before macroscopic fracture in Figure 7 is much larger than that in Figure 6, and the number of cracks in Figure 7b is much less than that in Figure 6b. The damage to the rock specimen accumulates with the crack propagation after each impact. With the number of cyclic impacts increasing, the cracks extend from the periphery towards the interior, as illustrated in Figure 7c–g. Finally, strip-shaped blocks are also observed at the periphery of the rock specimen, as depicted in Figure 7h.

3.3. Defining Damage Variable

Damage variables can be defined in various ways, such as elastic modulus, wave velocity, energy dissipation, and crack density. In the historical progression of damage definitions, crack density is one of the earliest methods to characterize damage variables. However, the practical measurement of crack density poses significant challenges. In numerical simulation, the crack generation and propagation result in the failure of elements. As failed elements are deleted in the finite element model, the volume reduction induced by element failure provides insight into the propagation of cracks. The volume of rock specimen after each impact can be recorded in numerical simulations with relatively low computational cost; then, the crack density method is employed to characterize the damage variable. As the crack density can be characterized by volume reduction after each impact in numerical simulation, the damage variable can be expressed as:

$$D=1-\frac{V_i}{V}$$

(4)

where D is the damage variable, and V_i and V are the volume of rock specimens after ith impact and before impact, respectively.

The volume of rock specimens can be derived from the volume–time history curve obtained by the post-processing software LS-PREPOST. The volume of the cylindrical rock specimen before impact is 49.035 × 10⁻⁶ m³. The smaller the volume of the rock specimen after impact, the greater the damage to the rock specimen.

3.4. Damage Evolution Analysis

Rock is a typical brittle material in nature, characterized by numerous microcracks in it. When subjected to an external load or energy, rock experiences closure of original fissures, initiation and propagation of new cracks, and eventual coalescence of these cracks. This process significantly impacts the macroscopic behavior of rock. By incorporating the volume of the rock specimen after cyclic impacts into Equation (4), the damage of the rock specimen under various impact conditions is obtained. The calculated results are presented in

Table 3. Damage of rock specimen after cyclic impacts.

Pmax/MPa	Pc/MPa	n	Vi/(10−6 m3)	D
170	0	1	47.024	0.041
		2	46.730	0.047
		3	46.632	0.049
		4	45.455	0.073
		5	44.426	0.094
		6	40.552	0.173
180	0	1	45.556	0.071
		2	43.73	0.108
		3	33.736	0.312
	2	1	47.851	0.024
		2	47.417	0.033
		3	47.221	0.037
		4	46.730	0.047
		5	45.848	0.065
		6	43.445	0.114
		7	39.571	0.193
	4	1	48.937	0.002
		2	48.839	0.004
		3	48.741	0.006
		4	48.643	0.008
		5	48.250	0.016
		6	47.858	0.024
		7	46.927	0.043
		8	46.338	0.055
		9	45.014	0.082
		10	42.857	0.126
		11	32.804	0.331
190	0	1	42.611	0.131
		2	22.066	0.550

, and the trends of damage evolution are depicted in Figure 8.

As illustrated in Figure 8, the cumulative damage of the rock specimen under confining pressure progressively increases with the number of cyclic impacts, exhibiting three distinct evolutionary stages: rapid rising, steady development, and sharp rising. In the numerical simulation, as there are no pre-existing cracks in the rock specimen, the initial impact generates the cracks in the rock specimen, and the subsequent impacts lead to the propagation of the previously generated cracks. As the constraint effect of confining pressure, damage accumulates at a low speed after each impact. Once the damage reaches a specific level, the damage accelerates to form a macroscopic fracture due to the rapid propagation of cracks.

Both confining pressure and amplitude of incident stress waves significantly affect the damage evolution and the number of cyclic impacts. As seen from Figure 8a, under 0 MPa confining pressure, the damage of the rock specimen at 190 MPa incident stress wave is more severe than that at 170 MPa. The high incident stress wave leads to many cracks generated after the first impact and a small number of cyclic impacts. The higher the amplitude of the incident stress wave, the greater the damage after each impact and the smaller the number of cyclic impacts before macroscopic fracture. As seen from Figure 8b, under the same number of cyclic impacts, the damage of the rock specimen under 4 MPa confining pressure is significantly less than that under 0 MPa confining pressure. Due to the existence of confining pressure, the development of transverse cracks in the rock specimen is constrained; namely, the confining pressure has a significant constraint effect on the propagation of cracks in the rock specimen. To a certain extent, the constraint effect offsets the generated transverse tension and reduces the axial tensile cracks. Also, the constraint effect increases with the increase of confining pressure. The larger the confining pressure, the less the damage, and the bigger the number of cyclic impacts.

4. Discussion

Numerical simulation provides an efficient, flexible, and economical means to simulate the dynamic response of rock under complex impact loads. In this study, the numerical simulation of cyclic impact tests on rock materials under confining pressure was successfully realized through LS-DYNA by first applying static circumferential confining pressure before impact load with dynamic relaxation, then conducting cyclic impact tests with full restart analysis. Subsequently, the cyclic impact damage evolution of rock materials under confining pressure was analyzed based on simulation results. As only the HJC model is employed for rock materials in numerical simulation, the applicability of this study may be limited. Follow-up studies may consider more constitutive models for comparative studies to improve the applicability of simulation results. Moreover, the parameters of the constitutive model may also vary to represent various kinds of rocks.

Based on numerical simulation, both confining pressure and amplitude of half-sine incident stress waves show a significant effect on the damage evolution of rock materials. As the dynamic response of rock material under an impact load is closely related to the waveform of impact load, follow-up studies may consider various waveforms of incident stress waves, such as rectangular stress waves, trapezoidal stress waves, and triangular stress waves. The damage variable is characterized by volume change induced by failure elements in the finite element model. Therefore, the grid of the rock specimen has an influence on the accuracy of calculated damage. However, small element size consumes much time for numerical simulation. Therefore, follow-up studies may use the discrete element method, such as DEM and PFC, to investigate crack generation and propagation during cyclic impacts under confining pressure.

5. Conclusions

Numerical simulation of cyclic impact tests on rock materials under confining pressure has been carried out by LS-DYNA using dynamic relaxation and full restart analysis. This analyzes the cyclic impact damage evolution of rock materials under confining pressure by defining damage variables with the crack density method. The conclusions are drawn as follows.

(1)

Cyclic impact tests of rock materials under confining pressure can be simulated first by applying static confining pressure with a dynamic relaxation process, then conducting cyclic impact tests with full restart analysis.

(2)

In numerical simulation, the crack generation and propagation result in the failure of elements in the finite element model. Thus, the damage variable defined by the crack density method can be characterized by volume reduction.

(3)

The number of cyclic impacts before macroscopic fracture increases with the increase of confining pressure and decreases with the increase of the amplitude of the incident stress wave.

(4)

The cumulative damage of rock under confining pressure progressively increases with the number of cyclic impacts, and damage evolution exhibits three distinct stages: rapid rising, steady development, and sharp rising. A high incident stress wave leads to serve damage after the first impact and a small number of cyclic impacts.

(5)

Under confining pressure, the damage of rock is alleviated due to the constraint effect on crack propagation. Large confining pressure results in minor damage and a large number of cyclic impacts.