超高应变率下岩石的破坏机理

赵高峰 孙建华 李世金 张奔

摘 要：通過对大理石的轻气炮试验进行三维数值仿真重现，探索了岩石在超高应变率下的破坏机理。大理石试样中两个高速应力计的实测应力信号被用作数值模拟的匹配目标，结果表明，在高速冲击压缩载荷下，细观静水压破坏是大理石的主要破坏机制，而不是通常认为的细观剪切或拉伸破坏。此外，为再现不同冲击速度下大理石的物理实验观测数据，提出了新的考虑时间非局部效应的细观速率相关岩石破坏模型，并且发现对岩石轻气炮的数值模拟需采用准静态测量的弹性力学参数而不是超声测量的弹性参数。最后对具有不同孔隙度的虚拟岩石样本进行进一步的轻气炮数值试验，结果表明，岩石的孔隙度越大，在冲击压缩载荷下的动态强度和动力学率效应越小。

关键词：动态破坏;轻气枪试验;冲击压缩荷载;格元模型

1 Introduction

The light-gas gun， first developed for hypervelocity research in aeronautics， was later used in the field of material science. The light-gas gun test is mainly used in material science for： 1） obtaining the Hugoniot curve of the material[1-3]， 2） measuring dynamic failure strength[4]， 3） studying the high-pressure phase transition[5]， and 4） investigating impact-induced chemical reactions[6]. The light-gas gun test has also been applied to rock and rock-like materials. Grote et al.[7] conducted a flat plate impact test on cement mortar and concrete. They reported that the average flow stress of cement mortar and concrete increased to 1.3 GPa and 1.7 GPa， respectively， compared with the unconfined quasi-static compressive strength of 46 MPa and 30 MPa. The dynamic effect of reinforced concrete under shock compressive loading was also explored using the light-gas gun test[8]. The bearing capacity of the reinforced concrete was found to be improved with the increase of impact velocity and reinforcement ratio[8]. Recently， similar phenomena were observed by Zhang et al.[1] with marble and gabbro.

The dynamic elasticity theory for solid materials was adopted to interpret experimental data of the light-gas gun test[1，7-8]. The constitutive model based on the dynamic elasticity theory can be used with the finite element method （FEM） for numerical simulation of the light-gas gun test. Liu et al.[9] simulated the dynamic fracture mode of tungsten alloy by LS-DYNA （an explicit FEM code） and captured the characteristics in their experimental results. Lopatnikov et al.[10] simulated the light-gas gun test of foamed aluminum plates using LS-DYNA and demonstrated that the dynamic deformation-time relationship and kinetic energy changes of the impact plates at different impact velocities predicted by the FEM were consistent with their theoretical model. Nevertheless， as pointed out by Riedel et al.[11]， the intrinsic heterogeneity of rock-like material might result in challenges to the determination of the equation of state parameters and the need for a mesomechanical model. Duan et al.[12] incorporated a statistical isotropic elastic micro-crack model into the FEM to simulate the propagation of plane shock waves in soda-lime glass and reproduced their experimental data. Recently， besides FEM with the mesomechanical model[13-14]， a number of discontinuum-based numerical methods[15-20] have also been applied in the study of the dynamic failure of rock. The molecular dynamics （MD）[15] and the lattice type model[18-20] were successfully applied to study the dynamic failure of rock. However， most of these studies only reproduced the failure mode of the rock using a mesoscopic constitutive model with tensile failure. The implementation and practicability of the model to simulate the light-gas gun test is still unclear.

In this work， we used a lattice type model to explore the light-gas gun test conducted on marble by Zhang et al.[1]. Our main purposes were to answer the following questions：

1） What kind of mesoscopic damage information can be obtained from the experimental data of the light-gas gun test？

2） What kind of elastic parameters， the dynamic elastic ones or the quasi-static elastic ones， control the shock loading propagation within the rock specimen？

3） How to describe the rate dependency of rock under shock compressive loading in the light-gas gun test using a mesoscopic constitutive model？

4） What is the role of the rock's mesostructure in the dynamic strength and dynamic effects （strain/loading rate dependency） of rock under shock compressive loading？

In this work， the distinct lattice spring model （DLSM） was adopted as the numerical tool due to its advantages in modeling rock failure based on simple mechanical elements， e.g.， Newton's second law and springs with simple constitutive models. These characteristics of the DLSM make it a suitable numerical tool for failure mechanism study. In addition， the DLSM is computationally appropriate for full 3D simulation. The following sections of this paper include a brief introduction of the physical test and numerical model， a brief description of the development of rate-dependent constitutive models and a comparison of model response with published data of the light-gas gun test conducted by Zhang et al.[1]， based on which the failure mechanism of rock under shock compressive loading is explored. Finally， a few conclusions are derived to answer the aforementioned questions.

2 Methods

2.1 The light-gas test

The light-gas test data reported by Zhang et al.[1]， obtained using the experimental facilities at the Cavendish Laboratory， Cambridge， UK[21]， are adopted in this work. Fig.1 shows the basic working principle of the experimental facilities as well as the components of the marble specimen and the copper flyer. Different from the traditional rock mechanics tests， two stress gauges were placed in the composite specimen made of marble and PMMA （see Fig.1（b））. The two stress gauges were used to obtain the stress history curve of the corresponding location inside the specimen during the test， that is， the waveform of the shock wave. Since there is a certain distance between the two stress gauges， the velocity of the shock wave in the rock specimen can be obtained based on these signals. After the light-gas gun test， both the composite rock specimen and the flyer became powdery and dissipated， as a result， no morphological data was recorded for the failure process， or the failure pattern of the rock specimen. Another important parameter of the light-gas gun test is the impact velocity of the copper flyer， which is the only controllable input of the light-gas gun test. A number of different impact velocities were adopted by Zhang et al.[1] to obtain the Hugoniot parameters of the marble. Basic material properties of the marble， PMMA and copper flyer were also provided in reference [1]. The longitudinal and shear wave velocities were obtained using an ultrasonic transducer， which can be further used to obtain the dynamic elastic parameters. The quasi-static elastic parameters of the marble were obtained by the standard uniaxial compression test， the parameters of which are listed in Table 1. The dynamic elastic parameters of PMMA and copper were close to their quasi-static counterparts， therefore， no quasi-static tests were conducted on the copper and PMMA.

2.2 Distinct Lattice Spring Model （DLSM）

The lattice spring model（LSM） was first developed by Hrennikoff in 1941[22]. Its basic principle is to represent the mechanical responses of a solid through a group of spring-like interactions. The LSM is regarded as the ancestor of both the FEM and the discrete element method （DEM）. Due to its simplicity， the LSM has been widely used in the study of many fundamental mechanical phenomena of solids[23-25]. The distinct lattice spring model （DLSM） was developed by Zhao et al.[26] to overcome the Poisson's limitation in the classical LSM. As shown in Fig.2（a）， the basic principle of the DLSM is to represent the solid as a group of particles linked through spring bonds， which consist of a normal spring and a shear spring. Defining the normal unit vector： n=（nx ny nz）T is directed from particle i to particle j， which are connected by a normal spring， and the normal deformation of the spring is defined as

in which ks is the stiffness of the shear spring and u*s is the ultimate shear deformation.

Eqs.（2） and （4） are the constitutive models adopted in the DLSM for brittle solids. Fig. 3（a） and （b） show the constitutive model for the normal spring and the constitutive model for the shear spring， respectively. Since a solid is represented as a spring network in the DLSM （as shown in Fig. 2（a））， the fracturing process is manifested by the gradual breakage of these springs （mesoscopic failure events）. It has the following advantages in simulating the fracturing of brittle materials. First， these constitutive models are easily implemented， and the simulation results are easily explained， as fewer parameters are required compared to other， discontinuous models， e.g.， the bonded DEM. Third， the straight forward parallelization of the DLSM[27-28] gives the model relatively high computational performance.

Another distinctive feature of the DLSM is that the relationship between the spring parameters（kn and ks） and the macroscopic material properties （E and v） is derived from the hyper-elastic theory. No calibration is required for the mesoscopic elastic parameters， as they can be calculated using the following equations.

where li is the initial length of the ith spring and Vm is the volume representing the computational model. More details of the DLSM and its recent development can be found in reference [25-27].

2.3 Computational model for the light-gas gun test

Fig.2（b） shows the computational model for the light-gas gun test conducted by Zhang et al.[1] on a marble specimen. Here， the marble specimen is represented as a group of mass points linked by springs. The dynamic contact between the particles was considered in the DLSM， so the impacts on the light-gas gun test could be directly simulated. The influences of thermal and chemical reactions were not considered in this work. Therefore， its mechanical response is characterized by a mass-spring system based on Newton′s second law and Hooke's law[26]. As shown in Fig.2（b）， the computational model consists of a copper flyer， a composite marble specimen and a PMMA plate. This three-dimensional numerical model can be viewed as a numerical test parallel to Zhang et al[1]， as shown in Fig.1（b）. The specific geometric parameters are consistent with the physical experiment. The only boundary condition input is the corresponding impact velocity of the copper flyer. For a closer comparison with the outputs of the actual experiment in Zhang et al.[1]， combined with the characteristics of the DLSM， a numerical stress gauge scheme was developed to record the stress data in the simulation. As shown in Fig.2（c）， a rectangular plane was placed in the model， and its normal vector is nA. During the calculation， assume the normal unit direction vector of the spring bond cut by this plane is nbi， and the total spring bond force accounting for the normal and shear springs is Fbi（t）=Fni（t）+Fsi（t）， then the corresponding stress of the numerical stress gauge is given as

Eq. （8） gives the time history of the normal stress on the stress gauge.

The particle size used in this work is 1 mm and there are about 100 000 particles. The input parameters of the light-gas gun test are based on experimental tests. The output of the numerical test is controlled by the selection of the meso-mechanical constitutive model and the corresponding constitutive parameters. A key point of this work is to minimize the difference between the numerical prediction and the physical experimental result. Thus， the failure mechanism of the rock under shock compressive loading may be interpreted from the searching process of the meso-mechanical constitutive model and its parameters.

2.4 Tri-linear hydro-compressive model

Many researchers suggest that macroscopic compression and shear failure is essentially tensile failure at the mesoscopic level[28-30]. However， most of those studies are limited to tensile failure of the normal spring[18， 20， 30]. Fig.3（a） shows the simple brittle constitutive model of the normal spring considering only tensile failure. More comprehensive forms， e.g.， considering the damage and plastic deformation， can be developed. Jiang and Zhao[31] recently developed a coupled damage plasticity model to describe the dynamic crack propagation of a gypsum-like 3D printing material. When considering the interaction between two particles， shear failure is a natural logical extension. However， there are very few studies on shear failure in classical LSMs due to the absence of shear interaction. In the DLSM， shear failure can be considered due to the introduction of the multi-body shear spring. The corresponding constitutive model is illustrated in Fig.3（b）. Shear failure was widely adopted in the bonded DEM， which usually takes into account both the influence of the cohesion of the contact bond between two particles and the frictional angle. The shear constitutive model presented in Fig.3（b） can be viewed as a special case when the frictional angle is zero. Therefore， the prediction of shear fracturing of the DLSM using the constitutive model shown in Fig.3（b） is not conservative.

When considering the interaction between two particles， it can be suggested that there is another possibility of damage， that is， the compression failure of the normal spring， which is related to the failure of the material under hydrostatic compression. If there is no hydrostatic failure involved， there is no need to introduce the hydro-compressive failure model into the DLSM. In this work， a trilinear constitutive model （see Fig.4） is introduced to describe the hydrostatic failure of the rock. Four parameters are required to determine the constitutive model without considering the loading and unloading. The first parameter， uc1， represents the first hydro-compression yield of the normal spring. The hardening or softening stage starts， followed by the normal spring yields. The second parameter， βredc1， corresponds to the reduction factor of the spring stiffness of the first yield. When it is greater than zero， it represents hardening， and when it is smaller than zero， it represents softening. The other two parameters， uc2 and βredc2， correspond to the secondary hydro-compression yield deformation and the stiffness reduction coefficient， respectively. The mathematic equation of this tri-linear constitutive model is as follows.

The response of the constitutive model is controlled by these four parameters，uc1，uc2，βredc1and βredc2

To introduce the effects of plasticity and damage， the same approach as that adopted by Jiang and Zhao[31] is used. First， rewrite Eq.（10） as follows.

where S（u） is a state parameter of the normal spring， which is defined as follows.

The interpretation of Eq.（10） can be expressed as the following formula if it is based on damage mechanics.

Then， the damage function can be expressed as follows.

In the calculation of the spring force， the maximum damage variable corresponding to the loading history is recorded as D*（u）， then the damage constitutive equation considering the loading and unloading can be written as

The loading and unloading response of this constitutive model is as shown in Fig.4（a）. If Eq.（10） is regarded as a plastic response， it can be expressed as

Compare Eq.（14）， we will have the plastic deformation up as

Similarly， to record the maximum plastic deformation up* experienced by the spring at the moment of loading， the pure plastic constitutive model considering loading and unloading can be written as

The corresponding pure plastic constitutive model is shown in Fig.4（b）. Combined with Eq.（15） and （18）， a coupled damage-plastic model can be constructed by introducing the damage-plastic coupling coefficient as follows.

When λdp=1， it represents the pure damage constitutive model， and when λdp=0， it represents the pure plastic constitutive model. In this work， we consider the copper behavior as pure plastic material， while the rock has a pure damage response.

Under dynamic loading， the material has obvious dynamic rate effects. To capture this phenomenon， it is convenient to introduce a rate-dependent model. The approach developed by Zhao et al.[18] is adopted in this work to bring rate dependency into the tri-linear constitutive model. The idea is to change the compressive strength parameters of the trilinear constitutive model as a function of the spring deformation rate. It was found that the instantaneous value of the spring deformation rate may not be able to reproduce the correct macroscopic dynamic effect[18]. To solve this problem， Zhao et al.[18] proposed the concept of time non-localization， which uses the average value of the spring deformation rate from the start time to the current time. In this work， following a similar idea， the local average deformation rate is given as

The average deformation rate of the spring is only calculated when the spring is deformed beyond 0.99uc1. Based on this concept， the strength parameters of the corresponding dynamic tri-linear constitutive model are given by the following formula.

where udync1 is the first yield deformation considering the dynamic effect， udync2 is the second yield deformation considering the dynamic effect， vref is the reference deformation rate， and α is a dimensionless coefficient. Fig.5（a） shows three different loading curves; each represents a different loading rate （deformation rate）. The loading rate has alternating positive and negative values， which is common in dynamic loads. Fig.5（b） shows the dynamic response based on the proposed dynamic tri-linear constitutive model， which shows good robustness for dynamic/cyclic loading. The dynamic constitutive model contains a total of six parameters with clear geometric and physical meanings. The selection of the parameters in the actual simulation process will be further illustrated in the following section.

3 Numerical modeling and discussion

3.1 Elastic parameters selection

A proper selection of parameters is fundamental to numerical simulation. For the light-gas gun test， Zhang et al.[1] performed two kinds of test to obtain the elastic parameters of their marble specimens. They conducted ultrasonic tests and obtained the wave velocities， which is a main variable of the dynamic elastic modulus and Poisson's ratio， as shown in Table 1. The quasi-static elastic modulus and Poisson's ratio of the marble specimens were also obtained through the traditional uniaxial compression test （see Table 1）. The quasi-static elastic parameters of the copper and PMMA are close to the dynamic ones， therefore， they were not tested and reported in reference [1]. Now， the challenge is the selection of the elastic parameters， dynamic or quasi-static ones， for the numerical modeling of the light-gas gun test. Here， the light-gas gun test with an impact velocity of 385 m/s is used as our modeling target. During the test， stress signals were recorded by two stress gauges， which were used to obtain the characteristics of the shock wave propagation within the specimen. It is reasonable to use the elastic parameters for ultrasonic wave propagation within the specimen to simulate the shock wave propagation. Therefore， dynamic elastic parameters are assigned to the marble as the first step， i.e.， E=55 GPa， v=0.31， for the numerical simulation. Here， a pure elastic dynamic simulation is performed by setting all failure parameters extremely large to prevent any damage from occurring， e.g. u*n=1 mm and u*s=1 mm. In order to achieve a quantitative comparison between the experimental data and the numerical data， an adjustment of the initial recording time is performed to make the initial waveform of the first numerical stress gauge match the corresponding physical signal. The time shift value was also applied to the signal of the second numerical stress gauge. As shown in Fig.6（a）， the initial stage of the numerical prediction of the second numerical stress gauge is significantly earlier than the experimental observation. This indicates that the numerically predicted shock wave velocity is markedly larger than that of the physical observation when the dynamic elastic parameters are used for the light-gas gun test. The quasi-static elastic parameters for the marble specimen， i.e.， E=40 GPa， v=0.20， are also adopted. The simulation results are shown in Fig.7（b）. It is noted that a better fitting is observed. Therefore， in the subsequent numerical modeling， for the marble specimen， the quasi-static elastic parameters are selected as the input values. In terms of the stress peaks， the numerical results are also much larger than the experimental observation， revealing that damage/failure occurred within the specimen during the light-gas gun test， otherwise a much higher peak value would be detected in the physical test. In the following section， we will further explore the corresponding failure mechanism of the specimen under shock compressive loading.

3.2 Failure mechanism of the marble specimen in the light-gas gun test

In this section， the mesoscopic failure mechanism of the marble specimen during the light-gas gun test conducted by Zhang et al.[1] is explored. First， the failure parameters were set as u*n=0.000 1 mm， and u*s=1 mm， which means only tensile failure was considered. According to the empirical equation proposed by Zhao[32]， the represented macro tensile strength is about u*nE/d-=4 MPa， where E is the elastic modulus and d- is the mean particle size of the computational model. Fig.7（a） shows the corresponding simulation results. Fig.8 shows the shock wave propagation and damage evolution of the marble specimen when considering only the mesoscopic tensile failure. Initially， the copper sheet collides with the rock sample at a high speed to generate a shock wave， and contact occurs at 3.00 microseconds. As shown in Fig.8， the speed of the copper sheet is immediately transmitted to the rock specimen and the part in contact with the copper sheet immediately breaks down. With the propagation of the shock wave， the damaged area spreads continuously around the copper piece （Fig.8（b））. The damage first propagates in the direction of impact， then it expands radially， and finally， the entire specimen is damaged. However， according to the observation presented in Fig.7（a）， the influence on the stress waveforms at the last stage is minor. The amplitude of the stress wave recorded at the first stress gauge is not reduced at all， even if the marble specimen is fully damaged in tensile failure. Considering the mesoscopic tensile failure has little effect on the experimental observations， it is not attributable to the failure mechanism of the marble specimen in the light-gas gun test. It should be mentioned that using the DLSM the uniaxial compressive failure can be reproduced with only mesoscopic tensile failure[32]. It is speculated that the failure mechanism of the marble specimen may differ from the compressive failure observed in the classical uniaxial compression test. To further investigate the influence of shear failure between particles， the failure parameters are set as u*n=1 mm， and u*s=0.000 1 mm. The corresponding numerical calculation results are shown in Fig.7（b）. In the DLSM， the shear constitutive model corresponds to the case where the friction coefficient between particles is zero， that is， the shear failure is fully considered. From the observations of Fig.7（b）， it is concluded that the mesoscopic shear failure is also not able to account for the failure mechanism.

In the light-gas gun test， the copper flyer was completely destroyed. Whether the recorded waveform can be affected by the damage to the copper flyer needs to be investigated. Here one may consider the compressive failure of the normal spring of the copper using the hydro-compressive constitutive model， as shown in Fig.4. The corresponding failure parameters were selected as uc1=0.001 mm， uc2=0.200 mm， βredc1=0， βredc2=0 and λdp=0 to represent the ideal elastoplastic response of the copper. Fig.9（a） shows the numerical simulation results. Compared with the experimental results， the peak is still unable to meet the experimental value. The compressive failure parameters of the marble then were taken as uc1=0.003 mm， uc2=0.200 mm， βredc1=0， βredc2=0 and λdp=1. The compression failure parameter uc1=0.003 mm， corresponds to the macroscopic failure strength of about 120 MPa， which is the typical value of the uniaxial compression strength of marble. Fig.9（b） shows the corresponding numerical simulation results. By comparing the numerical and physical results， it is concluded that the compression failure between particles in the light-gas gun test is caused by the hydro-static compression failure of the rock rather than uniaxial compression failure. In the following section， the compressive constitutive parameters will be adjusted to make the computational model reproduce the experimental results.

The tri-linear constitutive model developed in this work consists of four parameters. First， the values of these four parameters are taken as uc2=0.200 mm， βredc1=0，βredc2=0 and λdp=1. By continuously adjusting uc1and performing numerical modeling， it is found that when uc1 is 0.030 mm， the peak value of the numerical shock wave of stress gauge 1 is closest to the experimental results （see Fig.10（a））. However， the reduction of the stress wave of the second stress gauge was obvious. This attenuation is caused by the relatively coarse particle size adopted in the simulation. For a full three-dimensional model， there are already 110 000 particles at a particle size of 1 mm. Considering the computing capacity limitation of the available hardware and that the main purpose is to explore the failure mechanism rather than precisely capture the attenuation of the shock wave， the computational model with the particle size of 1 mm is preferred in the following simulation with a focus on the comparison of the stress waveform at stress gauge 1. Besides the attenuation， the platform segment of the stress wave at stress gauge 1 was also not well reproduced. By further adjusting other parameters， it is found that by setting βredc1=0.1， uc2=0.1 mm， and βredc2=0， the numerical simulation results give a closer fitting to the physical test （see Fig.10（b））， which indicates that a light hardening has occurred for the hydro-compression failure of the marble under shock compressive loading.

3.3 Dynamic effect of the mesoscopic hydro-compressive failure

So far，in the experimental data of Zhang et al.[1] while the impact velocity of 385 m/s is well simulated， one needs to generalize it for other impact velocities adopted by Zhang et al. The strength Sx can be determined as the peak value of the stress waveform recorded at numerical stress gauge 1. Fig.11 shows the experimental data of the dynamic strength of the marble Sx versus different particle velocities. As suggested by Zhang et al.[1]， the relationship between impact velocity vi and particle velocity vp may be estimated as vp=0.67vi， which is adopted in the data process of this simulation. Additional numerical simulations using different impact velocities， i.e.， 550 m/s， 650 m/s， 800 m/s and 1 000 m/s， were conducted using the same computational model and constitutive parameters. The rate-independent （RI） constitutive model was first used. As shown in Fig. 11， the numerical prediction using the DLSM with the RI constitutive model failed to reproduce the dynamic effect of the strength value observed in the light-gas gun test. Then， the rate dependent （RD） constitutive model is used， in which the strength parameters of the tri-linear constitutive model are described by Eq.（21） and （22）. Fig.11 shows the numerical results when v0ref=255 mm/s and α=0.43. From the numerical modeling results， it is concluded that the mesoscopic hydro-compressive constitutive model has to include the rate-dependent effect. In other words， the mesoscopic failure of the marble under shock compressive loading includes a dynamic effect.

3.4 Influence of mesostructure

In this section， the influence of the mesostructure on the dynamic failure of rock under shock loading is explored. Fig.12 shows the corresponding computational models for the rock specimens with different mesostructures. These computational models are formed by randomly removing a given number of particles from the marble specimen. The porosity， which was used to characterizs different computational models， is defined as n=number of particles being excavated / total number of particles in the original model. As shown in Fig.12， four pore models with a porosity of 5%， 10%， 15%， and 20% are constructed. All the constitutive parameters are the same as those of the previous section. Fig.13 shows the numerical results of the predicted strength of those different porosity models under different impact velocities. Overall， larger porosity results in lower dynamic strength， which seems logically correct. However， there are some variations， for example， when porosity is 5%， a higher strength is obtained. The reason might be that the inter-particle velocities become more violent under shock compressive loading when a small number of particles were removed from the original model. It might trigger the dynamic effect of the constitutive model and consequently， result in a higher strength. Nevertheless， when the porosity continuously increases， the inter-particle velocity increase may be released in the lateral direction due to the existence of a large mesoscopic free surface. This may result in the increase of strength due to the dynamic effect becoming inconspicuous （as shown in Fig.13）.

4 Conclusions

In this work，light-gas gun tests on marble were numerically investigated using the DLSM to study the failure mechanism of rock under compressive shock loading. Based on a detailed comparison of the numerical response with published experimental data， it was found that the elastic parameters controlling the compressive shock wave propagation in the rock specimen should be determined through quasi-static tests rather than ultrasonic tests. Moreover， mesoscopic tensile failure and shear failure between rock grains are not the mechanisms of rock failure under compressive shock loading. To reproduce the shock stress waveform observed in the light-gas gun test， a mesoscopic hydro-compressive constitutive model is needed. In other words， it is reasonable to conclude that the light-gas gun test may be able to provide information about the hydro-compressive failure of rock under dynamic loading and provide calibration data for the DLSM or other numerical methods to determine the constitutive equations under dynamic hydro-compression loading conditions. It is noted that a rate-dependent hydro-compressive constitutive model is also necessary to reproduce the experimental observation of the strength increase under different impact velocities. Finally， the influence of the mesostructure in terms of porosity was found to decrease both the dynamic strength and the dynamic effect. The findings in this work may provide a better understanding of rock failure under compressive shock loading， which might be useful for a more rational protective design of underground rock engineering structures under blasting loading[33-36].

Acknowledgements

This research is financially supported by the National Key R&D Program of China （Grant No. 2018YFC0406804） and the National Natural Science Foundation of China （Grant No. 11772221）.References：

[1] ZHANG Q B， BRAITHWAITE C H， ZHAO J. Hugoniot equation of state of rock materials under shock compression [J]. Philosophical Transactions of the Royal Society A： Mathematical， Physical and Engineering Sciences， 2017， 375（2085）： 20160169.

[2] ROSENBERG Z， YAZIV D， YESHURUN Y， et al. Shear strength of shock-loaded alumina as determined with longitudinal and transverse manganin gauges [J]. Journal of Applied Physics， 1987， 62（3）： 1120-1122.

[3] GEBBEKEN N， GREULICH S， PIETZSCH A. Hugoniot properties for concrete determined by full-scale detonation experiments and flyer-plate-impact tests [J]. International Journal of Impact Engineering， 2006， 32（12）： 2017-2031.

[4] GRADY D E. The spall strength of condensed matter [J]. Journal of the Mechanics and Physics of Solids， 1988， 36（3）： 353-384.

[5] DUNN J E， FOSDICK R， SLEMROD M. Shock induced transitions and phase structures in general media [M]. Springer， 1993.

[6] SEKINE T. Shock wave chemical synthesis [J]. European Journal of Solid State and Inorganic Chemistry， 1997， 34： 823-833.

[7] GROTE D L， PARK S W， ZHOU M. Dynamic behavior of concrete at high strain rates and pressures： I. experimental characterization [J]. International Journal of Impact Engineering， 2001， 25（9）： 869-886.

[8] LIU H F， FU J， NING J G. Experimental study on the dynamic mechanical properties of reinforced concrete under shock loading [J]. Acta Mechanica Solida Sinica， 2016， 29（1）： 22-30.

[9] LIU H Y， SONG W D， REN H L. Dynamic response of tungsten-nickel-iron composites under impact loadings [J]. International Journal of Nonlinear Sciences and Numerical Simulation， 2009， 10（8）： 993-1004.

[10] LOPATNIKOV S L， GAMA B A， HAQUE M J， et al. High-velocity plate impact of metal foams [J]. International Journal of Impact Engineering， 2004， 30（4）： 421-445.

[11] RIEDEL W， WICKLEIN M， THOMA K. Shock properties of conventional and high strength concrete： Experimental and mesomechanical analysis [J]. International Journal of Impact Engineering， 2008， 35（3）： 155-171.

[12] DUAN Z P， ZHANG Y G， ZHANG L S， et al. Multi-thickness target plate impact experimental approach to failure waves in soda-lime glass and its numerical simulation [J]. International Journal of Nonlinear Sciences and Numerical Simulation， 2012， 13（1）： 3-11.

[13] TANG C A. Numerical simulation of progressive rock failure and associated seismicity [J]. International Journal of Rock Mechanics and Mining Sciences， 1997， 34（2）： 249-261.

[14] ZHU W C， LIU J， YANG T H， et al. Effects of local rock heterogeneities on the hydromechanics of fractured rocks using a digital-image-based technique [J]. International Journal of Rock Mechanics and Mining Sciences， 2006， 43（8）： 1182-1199.

[15] KRIVTSOV A M. Relation between spall strength and mesoparticle velocity dispersion [J]. International Journal of Impact Engineering， 1999， 23（1）： 477-487.

[16] SILLING S A. Reformulation of elasticity theory for discontinuities and long-range forces [J]. Journal of the Mechanics and Physics of Solids， 2000， 48（1）： 175-209.

[17] ZHOU X P， WANG Y T， QIAN Q H. Numerical simulation of crack curving and branching in brittle materials under dynamic loads using the extended non-ordinary state-based peridynamics [J]. European Journal of Mechanics-A， 2016， 60： 277-299.

[18] ZHAO G F， RUSSELL A R， ZHAO X B， et al. Strain rate dependency of uniaxial tensile strength in Gosford sandstone by the Distinct Lattice Spring Model with X-ray micro CT [J]. International Journal of Solids and Structures， 2014， 51（7/8）： 1587-1600.

[19] ZHANG Z N， YAO Y， MAO X B. Modeling wave propagation induced fracture in rock with correlated lattice bond cell [J]. International Journal of Rock Mechanics and Mining Sciences， 2015， 78： 262-270.

[20] ZHAO G F， XIA K W. A study of mode-I self-similar dynamic crack propagation using a lattice spring model [J]. Computers and Geotechnics， 2018， 96： 215-225.

[21] BOURNE N K， ROSENBERG Z， JOHNSON D J， et al. Design and construction of the UK plate impact facility [J]. Measurement Science and Technology， 1995， 6（10）： 1462-1470.

[22] HRENNIKOFF A. Solution of problems of elasticity by the framework method [J]. Journal of Applied Mechanics， 1941， 8： A619-A715.

[23] FINEBERG J， GROSS S P， MARDER M， et al. Instability in dynamic fracture [J]. Physical Review Letters， 1991， 67（4）： 457.

[24] GUILLARD F， GOLSHAN P， SHEN L M， et al. Dynamic patterns of compaction in brittle porous media [J]. Nature Physics， 2015， 11（10）： 835-838.

[25] ZHAO G F. Developing a four-dimensional lattice spring model for mechanical responses of solids [J]. Computer Methods in Applied Mechanics and Engineering， 2017， 315： 881-895.

[26] ZHAO G F， FANG J N， ZHAO J. A 3D distinct lattice spring model for elasticity and dynamic failure [J]. International Journal for Numerical and Analytical Methods in Geomechanics， 2011， 35（8）： 859-885.

[27] ZHAO G F， FANG J N， SUN L， et al. Parallelization of the distinct lattice spring model [J]. International Journal for Numerical and Analytical Methods in Geomechanics， 2013， 37（1）： 51-74.

[28] ZHAO G F， CHEN H. Performance of the parallel four-dimensional lattice spring model using Alibaba cloud [J]. Journal of Civil and Environmental Engineering， 2019， 41（3）： 1-11.

[29] POTYONDY D O， CUNDALL P A. A bonded-particle model for rock [J]. International Journal of Rock Mechanics and Mining Sciences， 2004， 41（8）： 1329-1364.

[30] JIANG C， ZHAO G F. On crack propagation in brittle material using the distinct lattice spring model [J]. International Journal of Solids and Structures， 2017， 118/119： 41-57.

[31] JIANG C， ZHAO G F. Implementation of a coupled plastic damage distinct lattice spring model for dynamic crack propagation in geomaterials [J]. International Journal for Numerical and Analytical Methods in Geomechanics， 2018， 42（4）： 674-693.

[32] ZHAO G F. Development of micro-macro continuum-discontinuum coupled numerical method [D]. EPFL， Switzerland， 2010.

[33] XIE H， ZHU J， ZHOU T， et al. Conceptualization and preliminary study of engineering disturbed rock dynamics [J]. Geomechanics and Geophysics for Geo-Energy and Geo-Resources， 2020. https：//doi.org/10.1007/s40948-020-00157-x.

[34] ZHOU X P， YANG H Q. Micromechanical modeling of dynamic compressive responses of mesoscopic heterogenous brittle rock [J]. Theoretical and Applied Fracture Mechanics， 2007， 48（1）： 1-20.

[35] LIU P F， ZHOU X P， QIAN Q H， et al. Dynamic splitting tensile properties of concrete and cement mortar [J]. Fatigue & Fracture of Engineering Materials & Structures， 2020， 43（4）： 757-770.

[36] ZHANG B， ZHAO G F. Dynamic failure and energy dissipation of rock ring specimen based on 4D lattice spring model [J]. Journal of Civil and Environmental Engineering， 2019， 41（2）： 20-28.

（編辑 章润红）