Abstract
In order to further reveal the failure mode and dynamic response law of dangerous rocks with different degree of weathering in the rock cavity under the action of earthquake and to provide early warning and forecast for steep slope of dangerous rocks in similar earthquake areas, a typical steep slope of dangerous rock in earthquake area of Sichuan, China, was taken as the research object, after detailed geological survey, and according to the chain development law of dangerous rock, the steep slope of dangerous rock before the earthquake was restored. Based on the 3D particle flow software PFC3D, the dangerous rock was divided into 3 modes according to the degree of weathering of the mudstone rock cavity, and the three-dimensional discrete element dangerous rock model under different modes was established. By introducing the horizontal and vertical two-way coupled seismic waves in Wenchuan, Sichuan, in 2008, the failure evolution process of steep slope of dangerous rock under the action of the horizontal and vertical coupled seismic waves was dynamically simulated, which proved the rationality of the simulation. The frequency spectrum of velocity-time history signal of each rock block in the dangerous rock model was analyzed by MATLAB programming, and the time-frequency characteristics of each dangerous rock model under the action of coupled seismic wave were studied. The research results have important scientific guiding significance and practical value for the dynamic stability evaluation and prediction of such steep slope of dangerous rocks under the combined action of rock cavity weathering and earthquake.
1. Introduction
Sichuan Province, China, is located in the Helan Mountain-Liupan Mountain-Longmen Mountain-Hengduan Mountain seismic zone, where earthquake disasters occur frequently. As one of the most important secondary geological disasters in the earthquake area, examples of rock slope instability induced by earthquake are common, often causing huge economic and property losses. The dangerous rock steep slope is the main carrier of seismic slope instability disaster, because the rock mass is cut by different structural planes, when the main control structural plane at the rear of the dangerous rock block gradually penetrates until it breaks under the action of multiple factors, instability and destruction of dangerous rocks will occur. Existing studies have mostly focused on the instability and failure of rock slopes in nonseismic areas under the action of self-weight. The dynamic response and failure mode of the dangerous rock slope under the action of the earthquake in the earthquake zone need to be solved urgently, especially the study on the process of collapse deduction, failure mechanism, and dynamic response characteristics of rock cavern-type dangerous rock slopes under different weathering degrees under the action of the earthquake gradually become one of the research hotspots and technical difficulties in this field [1, 2].
For example, Chen et al. [3] deduced the chain development process of dangerous rock in the Three Gorges reservoir area from the perspective of geomorphology, the influence law of rock cavity weathering on the stability of dangerous rock slope was revealed, and the chain rule of dangerous rock development in the Three Gorges reservoir area was put forward. Based on fracture mechanics and damage mechanics, Tang et al. [4, 5] studied the rockfall sequence on the steep cliff of weak base and established the calculation method of rockfall time. Based on fracture mechanics and material mechanics, Wang et al. [6] established a method for calculating the interlayer load and stability coefficient of complex slow-dipping rock slope and revealed its collapse sequence with the finite element method. Zhou et al. [7] used the QUIVER code to realize the response analysis of the gentle slope under the action of earthquake and revealed the slope response and sliding displacement law under the action of earthquake. Huang et al. [8] used the discontinuous deformation analysis (DDA) method to simulate the movement characteristics of dangerous rock slopes under different seismic conditions and revealed the influence of seismic loads on the movement characteristics of collapse block of dangerous rock mass. Zhao et al. [9] used the particle flow software PFC2D to dynamically simulate the instability process of the rock slope and revealed the influence of the distribution of the soft structural plane on the instability and failure of the rock slope. Based on the cantilever beam calculation theory, Yuan et al. [10] used the time-history curve of seismic acceleration to obtain the time-history curve of dangerous rock under the maximum tensile stress and established the calculation method for the stability of dangerous rock under the action of the earthquake. Bian et al. [11] used the technology of PFC2D artificially synthesized rock mass to reveal the failure mode and dynamic response law of rock slopes with bedding and intermittent joints with different combinations of rock bridge inclination and joint spacing under the action of the earthquake. Based on the discrete element method (DEM), Zhao et al. [12] deduced the progressive failure evolution process of a three-dimensional steep rock slope under the action of seismic transverse and longitudinal waves.
The above studies are mostly based on the single factor effect of rock cavity weathering or seismic action, and only two-dimensional perspective is used in the discrete element numerical simulation process, which cannot directly reveal the instability and failure process of dangerous rock slope. As a result, it has poor applicability to the rock cavity type dangerous steep slope commonly existed in earthquake area and insufficient understanding of the law of seismic dynamic response and failure process of dangerous rocks. Based on this, this paper takes the dangerous rock slope of rock cavity type in Sichuan Province, China, as the research object, the combined conditions of rock weathering and seismic action were considered, discrete element particle flow software (PFC3D) was used to simulate the dynamic failure process, the dynamic response rule and failure mechanism under the action of earthquake were revealed, and the energy distribution and frequency characteristics were discussed. This method can directly deduce the dynamic process of the collapse of each dangerous rock block under the action of earthquake, improve the chain law of dangerous rock, and reveal the damage degree of each rock block affected by the earthquake under different weathering degrees of rock caves. It is of great scientific significance and practical value to evaluate and forecast the dynamic stability of the dangerous rock slope under the combined action of rock cavity weathering and earthquake.
2. Establishment of Dangerous Rock Model
2.1. Model Instance
In order to better describe the elastic-plastic behavior of soil around the pile group of wharf foundation under wetting-drying cycles and lateral cyclic loading, the Mohr-Coulomb model is utilized to simulate the elastoplastic behavior of the bank slope. The yield function of the Mohr-Coulomb model in terms of stress invariants can be expressed as follows: this paper takes a typical steep slope of dangerous rocks in Beichuan Qiang Autonomous County, Mianyang City, Sichuan Province, China, as the research object; as shown in Figure 1, the 2008 earthquake in Wenchuan, Sichuan Province, in China caused the collapse and instability of the dangerous rocks. There are no residents living below the dangerous rocks point, and it is far away from Provincial Highway 205, no casualties were caused by the earthquake. Most of the collapsed rocks were scattered at the foot of the slope, but there is a rock block with a volume of about 1.76 m3 still rolling down beside the Provincial Highway 205, posing a huge threat to the safe operation of Provincial Highway 205 and the nearby residents.

After detailed geological survey, the steep slope of dangerous rocks was restored. The height of the steep slope is about 23.0 m, and the upper rock mass is composed of arkstone sandstone with good integrity, height of about 18.0 m, natural bulk density of 24.7 kN/m3, and fracture toughness of 26.0 Mpa·M1/2. The lower rock mass is mudstone with poor integrity, a height of about 5.0 m, and a natural bulk density of 25.3 kN/m3. According to the chain law of the dangerous rock [3] and the structure plane layout of the site, there are two macrochains A and B on the steep slope of the dangerous rock from the slope to the inside of the slope. Among them, the width of chain A is about 1.8 m and the width of chain B is about 1.7 m. Each macrochain contains 4 microchains, and the dangerous rock mass is numbered from bottom to top (A chain is A1, A2, A3, and A4, and B chain is B1, B2, B3, and B4), as shown in Figure 1. The size and fracture size of the rest rock mass are shown in Table 1.
2.2. Establishment of Numerical Model
2.2.1. Model Establishment
The evolution process of natural weathering of shale cavity in this case is considered; after the reasonable generalization, the discrete element software PFC3D was used to establish 3 groups of rock cavity type slope models with different weathering degrees. The cavity depth of Model 1 to Model 3 is gradually deepened, which is 0.9 m, 2.1 m and 3.5 m, respectively, and Model 3 corresponds to the dangerous rocks in the examples before the earthquake in this paper, as shown in Figures 2(a), 2(b), and 2(c). With reference to [13], the linear parallel bond model was selected as the internal contact model of the rock in this paper, and the smooth-joint contact model was selected as the contact model of the joint plane between the rocks. As shown in Figure 2(d), the model joint surface is generated by the DFN rock joint system built in the PFC, and the joint mesoscopic parameters are controlled by it.

(a)

(b)

(c)

(d)
In this paper, rock mesoscopic parameters were tested by rock triaxial test and PFC numerical triaxial servo test [13, 14]. The obtained rock mesoscopic parameters are shown in Table 2, and joint parameters were referenced [11], as shown in Table 3.
2.2.2. Setting Up of Monitoring Scheme
In order to obtain the failure law of rock block movement under the action of earthquake, this paper intends to set up 16 monitoring points from bottom to top and from left to right along the vertical direction of A and B macrochains. The velocity, acceleration, and other characteristic quantities of rock particles are monitored, and the rock blocks are grouped to monitor the displacement, velocity, acceleration, and other characteristic quantities of each rock block of Model A1 to B4, as shown in Figure 2(c).
2.2.3. Seismic Wave Input
Referring to reference [11] and reference [15], in this paper, the first 25s of the E-W and U-D seismic time-history records obtained at Wolong Station in the Wenchuan Earthquake of Sichuan Province in 2008 recorded by the China National Earthquake Data Center are used as the horizontal and vertical two-way coupled seismic wave as the input seismic time-history wave, and the time-history curve of acceleration is shown in Figure 3.

(a)

(b)
3. Failure Mode and Dynamic Response of Rock Vibration
3.1. Establishment of Numerical Model
The failure modes of Model 1 and Model 2 are shown in Figure 4 when the earthquake occurs for 25 s. In Model 1, the development crack on the left joint surface of Chain A is LX1, and the development crack in the left joint surface of chain B is LX2. In Model 2, the development crack on the left joint surface of A chain is LX3, and the development and change of each crack with time during the seismic period are shown in Figure 5. It can be seen from Figures 4 and 5, the failure modes of each model under the action of earthquake are different due to the weathering degree of mudstone cavity. The crack between rock blocks develops gradually with the extension of the seismic action time, and the slope of each surface reflects the rate of crack expansion, indicating that the crack expansion rate is synchronized with the time-history curve of the input seismic acceleration, and the degree of rock collapse is positively correlated with the degree of weathering of the mudstone cavity. In Model 1, the weathering depth of cavity is 0.9 m. After the earthquake, the crack LX2 expands first, and the expansion speed and depth are much greater than that of LX1, causing the B chain to tilt outward. However, since the weathering depth of the cavity does not touch the location of crack LX2, the rock blocks of chain A and B did not collapse during the period of earthquake action, but the rock blocks of chain B showed a tendency to dump outwards as a whole chain. It can be seen that the potential failure mode of Model 1 can be transformed into tearing-toppling failure. In Model 2, because the weathering of the rock cavity at the bottom of chain B deepened to 2.1 m, the upper rock block was suspended. 25 s after the earthquake, the B1 rock block slipped and failed, and rock blocks B2, B3, and B4 suffered toppling failure. The A chain did not lose stability due to the low weathering of the rock mass at the bottom; it can be seen that the basic failure mode of Model 2 is slip-rotation damage.


As the weathering depth of the rock cavity in Model 3 is the largest and the model produces overall instability collapse failure after 25 s of earthquake action, therefore, Model 3 is adopted for key analysis, and the failure process of instability and collapse during the period of earthquake action is shown in Figure 6. Combining the time-history curve of the seismic wave in Figure 3, it is found that the seismic wave starts to increase gradually at about 5 s, causing the joint planes around the rock blocks B1 and B2 at the bottom of the B chain to first locally fracture and then slip failure. When the seismic action occurred for 10 s, the seismic acceleration continued to rise and fluctuate, and the rock blocks B3 and B4 were collapsed successively. When the seismic action occurred for about 15 s, the B chain rock block completely collapsed and the rock blocks B3 and B4 collided with the rock blocks B1 and B2. At the same time, the left joint surface of the A chain expanded and the crack development direction gradually deepened with the vertical direction. When the seismic action occurred for 20 s, the left crack of the A chain was completely penetrated, and the four rock blocks in A chain showed signs of dumping and failure outward, and the rock block B4 continued to roll outward. When the seismic action occurred for 25 s, the entire chain of A chain was dumped outwards, the rock blocks A1 and A2 collide with the rock blocks B2 and B3, and the particle velocity decreased and in a relatively static state. The rock blocks A3 and A4 continue to dump outwards, and the top rock block A4 has the highest particle velocity, up to 0.43 m/s. In conclusion, it can be seen that the basic failure mode of Model 3 is slip-rotation-toppling failure, and its overall instability failure confirms the rationality of this simulation.

(a)

(b)

(c)

(d)

(e)
3.2. Evolution of the Number of Force Chains in Dangerous Rock Blocks
In the PFC model, due to the physical and mechanical properties of rocks formed by contact force chains between meso-particles, macroscopic cracks are formed when a large number of force chains break between particles. Therefore, the fracture force chain number of a rock block can reflect the degree of crack development of the rock block and directly reflect the failure development process of the rock block. Therefore, this paper monitors and analyzes the number of force chains of the block rock. In view of the fact that no collapse and instability of rock blocks occurred in Model 1 during the period of seismic action, the time evolution process of the number of force chain of rock blocks in Model 2 and Model 3 under the action of seismic longitudinal and longitudinal coupled waves was monitored, which was divided into 4 and 3 stages according to their fracture rates, as shown in Figure 7.

(a)

(b)
In Model 2, the crack initiation stage (Stage 1) is from 0 to 5 seconds after earthquake action, in which the number of force chains decays slowly, and the decay number is between 7 and 17. The fracture expansion stage (Stage 2) is from 5 to 10 seconds after earthquake action, in which the decay rate of force chain is accelerated, and the decay number is between 8 and 22. 10–15 s is the rock block caving stage (Stage 3), the attenuation rate of force chain is the fastest in this stage, and the attenuation number is between 9 and 30. 15–25 s is the rock block collision stage, the rock block has basically separated from the parent rock in this stage, and the attenuation of force chain number gradually tends to be stable with the slowest decay rate, and the attenuation number is between 2 and 15. In Model 3, due to the large degree of weathering in the mudstone cavity, the rock blocks collapse faster under the action of earthquake, resulting in a more violent evolution of the force chain number. Therefore, it can be divided into three stages, among which the earthquake period from 0 to 4 s is the crack expansion stage (Stage 1), the number of force chain decreases slowly as in Model 2 in this stage, and the attenuation number is 5 to 20. 4–15 s is the rock block caving stage (Stage 2), in which the number of force chains attenuates the fastest and the attenuation number is 45–76, indicating that the number of force chain fractures is the largest in this stage. 15 to 25 s is the rock block collision stage (Stage 3), the attenuation number of force chain gradually tends to be stable with the slowest attenuation rate, and the attenuation number is between 0 and 28. It is worth noting that the rock block B1 in Model 2 collapses preferentially in 5–10 s (Stage 2), resulting in the earlier and fastest attenuation of the force chain number in this stage. Similarly, the rock block B2 in Model 3 has partial fractures in the 0–4 s period (stage 1), the number of force chains decays earlier, and the attenuation rate is the fastest.
By comparing the failure process and the attenuation process of the number of force chains of A and B, it can be found that the B chain is more obviously affected by the earthquake, and the reason is that there is no pedestal support at the bottom of the B chain rock block. By comparing the two figures (a) and (b), it can be seen that the attenuation rate of the force chain number of the above dangerous rock model changes slowly within about 0–5 s at the initial stage of seismic wave action. The difference between the two is that in the middle period of seismic action, the weathering degree of rock cavity of Model 3 is greater than that of Model 2 within about 5–15 s, resulting in greater force chain fracture rate and fracture degree than Model 2. On the whole, the number of force chain fractures in Model 3 within 25 s of the seismic action is significantly higher than that in Model 2, indicating that the dynamic fracture response is positively correlated with the degree of weathering of rock cavity.
In conclusion, the dynamic stability and failure mode of the rock cavity slope under the action of earthquake are controlled by seismic wave action and the weathering degree of rock cavity. The failure process of dangerous rock is synchronized with the time-history curve of seismic acceleration, and with the deepening of weathering, the number of fracture of force chain and the degree of collapse of rock block increase gradually.
3.3. Dynamic Response of Dangerous Rock Blocks under the Action of Earthquake
According to the aforementioned rock particle monitoring points and rock block layout plan, time-history changes of displacement, peak velocity, and PGA amplification factor distribution of the dangerous rock block within 25 s of earthquake action are obtained, and the dynamic response of dangerous rock block under the action of earthquake is studied accordingly. Among them, with reference to [13], this paper defines the ratio of the dynamic response peak ground acceleration of the rock particle monitoring point or monitored rock to the average peak acceleration of the input seismic wave (9.42 m/s2) as the PGA amplification acceleration.
3.3.1. Displacement of Dangerous Rock Blocks
In view of the fact that no instability failure occurred in the dangerous rock block of Model 1 during the period of seismic action, the displacements of dangerous rock blocks of Model 2 and Model 3 were analyzed in this paper; the displacement time-history curve of the rock blocks in the vertical Z direction is shown in Figure 8. It can be seen that the displacement of the dangerous rock blocks in Model 2 and Model 3 can be divided into 4 and 3 stages, respectively, which are consistent with the segmented analysis of the force chain number of the dangerous rock blocks. In Model 2, due to the shallow weathering depth of the rock cavity, there is no significant displacement of the internal A chain rock blocks under the action of earthquake.

(a)

(b)
During the period of about 0∼5 s, the displacement of B chain rock block is in the initial stage (Stage 1), in which the displacement increases slowly and the displacement between rock blocks has no obvious change in this stage. During the period of 5–10 s, the displacement of B chain rock block increases slowly (Stage 2), and the displacement increases slowly, and the displacement of Block B1 increased dramatically and grew rapidly at about 7 s after the earthquake, indicating that at this time, it was the starting time for block B chain to collapse. Moreover, the rock block B1 has priority over other rock blocks to collapse, and the displacement increases by about 2.70 m, while the displacement of the other rock blocks B2∼B4 increases by 0.08–0.43 m. During the period of 10–15 s, it is the stage of abrupt increase of displacement (Stage 3), which has the fastest growth rate of displacement. Among them, the displacement of rock block B2∼B4 increases by 0.76–2.32 m, and the displacement growth of rock block B1 slows down at this stage, about 1.02 m. During the period of 15–25 s, it is the continuous displacement stage (Stage 4), the displacement of rock blocks all shows a slow growth in this stage, and the displacement of rock blocks B1∼B4 increases by 0.63–1.25 m. In Model 3, during the period of 0∼4 s, it is the initial stage of displacement (Stage 1), the displacement curve is the same as Model 2 in this stage, and the displacement has no obvious change. 4–15 s is the increase stage of displacement (Stage 2), in which the A chain rock block is still in a stable state. The displacement of B chain rocks increases sharply, and the rock blocks B1 and B2 have the priority to slip and collapse when the earthquake action is about 4 s, and their displacement increases rapidly at the initial stage, indicating that this time is the starting time of B chain block caving. About 10 s after the earthquake, the rock blocks B3 and B4 collapsed and the displacement increased sharply, and the displacement of the rock blocks B1∼B4 increased by 3.94–10.40 m in this stage. During the period of 15–25 s, it is the continuous displacement stage (Stage 3), the A chain rock blocks gradually collapses in this stage, and the displacement of rock blocks A3 and A4 increases sharply, reaching 7.59 m and 19.80 m, respectively; the displacement of the B chain rock blocks continued to increase, especially the rock block B4 has the largest displacement, reaching 6.2 m, while the remaining rock blocks increased by 0.94–2.21 m.
Comparing the displacement evolution process of the rock blocks of Model 2 and Model 3, it can be seen that the displacement of the dangerous rock block increases with the increase of the earthquake action time. Within 25 s, the peak displacement of dangerous rock blocks under the action of earthquake reaches 4.45 m and 20.6 m, respectively. The initiation time of caving of dangerous rock blocks is about 7 s and 4 s, respectively, but both of them show that the collapse of B chain is superior to the A chain, indicating that under the action of earthquake, the instability failure degree of the rock blocks of the gently dipping rock slope is positively correlated with the weathering degree of the rock cavity, which further verifies the conclusion of the analysis of the number of force chains. Moreover, the deeper the rock cavity is, the earlier the rock block collapses start, but the instability and failure of the rock blocks of external chain are earlier than that of the rock blocks of internal chain.
3.3.2. Peak Velocity of Dangerous Rock Blocks
Figure 9 shows the distribution of peak velocity of macroscopic chain and microscopic chain rock blocks in each model under the action of earthquake. It can be found that the peak velocity of Model 3 is greater than that of Model 2 on the whole, while Model 1 is the smallest, indicating that the peak velocity of dangerous rock blocks increases with the increase of weathering degree of rock cavity. The rock blocks in the macrochains A and B show that the peak velocity of rock blocks of the macrochain B is greater than that of rock blocks of the macrochain A, indicating that the rock blocks of external macrochain are more prone to instability and failure under the action of earthquake than the rock blocks of internal macrochain. From the microchain sequence, it can be found that the higher the microchain sequence or the higher the vertical elevation of the rock blocks, the peak velocity of the rock blocks tends to increase, indicating that there is an obvious velocity and elevation amplification effect in the microchain rock blocks of the dangerous rock under the action of earthquake.

3.3.3. PGA Amplification Coefficient of Dangerous Rock Blocks
The distribution of PGA amplification coefficient of the rock blocks of the macroscopic chain and microscopic chain in each model under the action of earthquake is shown in Figure 10. It can be found that, on the whole, the distribution of PGA amplification coefficient shows the rule that Model 3 is larger than Model 2, and Model 1 is the smallest, indicating that the PGA amplification coefficient of the dangerous rock blocks under the action of earthquake is positively correlated with the weathering degree of the rock cavity. The rock blocks of macrochains A and B show the rule that the PGA amplification coefficient of the rock blocks of macrochain B is larger than that of the rock blocks of macrochain A, indicating that the rock blocks of external macrochain responds more violently than the rock blocks of internal macrochain under the action of earthquake. It can be found that the PGA amplification coefficient of the rock blocks in the microchain sequence increases with the increase of the microchain sequence of the rock block or the increases of vertical elevation of the rock blocks, indicating that there is an obvious PGA elevation magnification effect on the rock blocks of microchain under the action of earthquake.

3.4. Dynamic Response of Monitoring Points on the Joint Plane under the Action of Earthquake
The peak velocity and PGA amplification coefficient of 16 monitoring points arranged along the vertical direction of the joint planes of the two macroscopic chains under the action of the earthquake are shown in Figure 11. It can be found from Figure 11 that in Models 1–3, the peak velocities of rock particle monitoring points of the A chain are between 0.16 and 0.22 m/s, 0.19–0.32 m/s, and 0.21–0.56 m/s, respectively. The PGA amplification coefficient is between 0.74 and 2.00, 1.50–2.20, and 1.05–3.40, respectively. Moreover, in Model 3, the peak velocity of the rock particles of A chain and the PGA amplification coefficient are distributed along the elevation in a W-shaped distribution, and the maximum is located at the bottom of the joint. On the whole, in Model 3, the peak value of rock particle velocity and the peak value of PGA amplification coefficient of A chain are the largest, followed by Model 2, and Model 1 is the smallest. It can be found from Figure 12 that in Models 1–3, the peak velocities of the rock particles monitoring points of the B chain are between 0.25 and 0.75 m/s, 0.56–2.34 m/s, and 0.82–4.39 m/s. The PGA amplification coefficients are between 2.23 and 6.05, 4.25–12.71, and 5.82–24.30, respectively. The maximum peak velocity and the maximum PGA amplification coefficient are both located at the top of the joint. Moreover, in Model 3, the peak velocity of rock particle of B chain and the PGA amplification coefficient are the largest, followed by Model 2, and Model 1 is the smallest.

(a)

(b)

(a)

(b)
It can be seen that when the weathering degree of rock cavity does not reach to two macro joint surfaces (Model 1), the seismic force only caused a partial break at the top of the external B chain joint, and there are no obvious fracture signs of A chain. At this point, the seismic dynamic response of the rock particles of external B chain presents an elevation amplification effect in the vertical direction. When the weathering degree of rock cavity is deep into the external macroscopic joint surface (Model 2), the seismic force causes the external B chain joint to break, and the internal A chain joint is to break locally. At this time, the seismic dynamic response in the vertical direction of the rock particles of external B chain also presents elevation amplification effect. After the weathering degree of the rock cavity is deep into the internal joint surface (Model 3), due to the propagation of seismic waves, the bottom of A chain inside the dangerous rock first receives greater seismic power, which is transmitted to the external B chain, making the joints at the bottom of the external B chain expand preferentially. When the seismic wave is transmitted to the dangerous rock at 12.5 m, the rock mass is intensified by the superposition of the vertical and horizontal seismic waves, resulting in a sharp increase in peak velocity and PGA amplification coefficient of rock particles. Firstly, the rock blocks B3 and B4 collapsed preferentially, and then the other rock blocks collapsed (Model 3). At this time, the seismic dynamic response of the rock particles of A chain is distributed along the elevation in a W-shaped distribution and the rock particles of B chain still present an elevation amplification effect.
In summary, under the combined action of earthquake and rock cavity weathering, the peak velocity and the PGA amplification coefficient of the rock particles monitoring points of internal A chain fluctuated greatly; the rock particle monitoring points of external B chain show an increasing trend with the increase of the vertical elevation of the rock particles, showing a significant elevation amplification effect. On the whole, the peak velocity and PGA magnification coefficient of the rock particles monitoring points of A and B chains are shown as the rule of the maximum in Model 3, the second in Model 2, and the minimum in Model 1, indicating that the seismic dynamic response of rock particles is positively correlated with the weathering depth of the rock cavity. By comparison, the peak velocity and PGA amplification coefficient of the rock particle monitoring points of external B chain are generally greater than those of the A chain, it can be seen that the rock particles of external chain are more affected by the earthquake. The results further perfect the chain development law of the dangerous rock under the action of earthquake and rock cavity weathering, and it has important scientific guiding significance and practical value for revealing the stability evolution process of the dangerous rock.
4. Spectrum Analysis of Vibration Signal
The last section mainly analyzed the failure mode and dynamic response of the dangerous rocks under the action of earthquake, which confirmed the rationality of the process of numerical simulation. In this section, from the perspective of structural dynamics, the time-frequency analysis of vibration velocity-time history signal of rock block within 25s of seismic action was carried out by using MATLAB programming. The obtained seismic time-frequency information can provide an important reference for the stability analysis and related prevention and control projects of these kinds of dangerous rocks under the action of earthquake.
4.1. The Basic Principle
HHT (Hilbert-Huang Transform) is an analysis method for processing nonlinear and nonstationary signals proposed by Huang E and Hilbert-Huang in 1998. Nowadays, it is a major breakthrough in linear and steady-state spectral analysis based on Fourier transform in the past hundred years and currently widely used in spectrum energy analysis of seismic wave [16]. It is based on empirical mode decomposition (EMD), which decomposes the signal into a series of intrinsic mode functions (IMF) and then performs Hilbert transform (HHT) on each IMF to obtain the instantaneous frequency [17].
The basic steps are as follows:(1)Firstly, EMD is used to decompose the complex signal into a finite number of intrinsic mode functions (IMF). Then, cubic spline function is used to interpolate the extreme points on the original signal x(t), connect the envelopes fitted by the upper and lower extreme points, and define the mean value of the two envelopes as m (t). The original signal is subtracted from the mean value to get where x (t) is the original signal, m (t) is the mean value of the upper and lower envelopes, and y (t) is the difference. Taking y (t) as the original signal and repeating the above steps until the stopping criterion is satisfied after n iterations, the criterion is Among them, the SD value is generally between 0.2 and 0.3 to stop the iteration. At this time, c1 = y1 (t) obtained c1 after n screening is the first IMF component, and then let the difference value between x (t) and c1r (t) = x (t)-c1 (t) be the new signal, and the IMF components c2, c3,..., cn can be obtained in turn. So far, the original signal can be decomposed into the sum of n IMF components and the margin rn: where c (t) is the modal component; r (t) is the residual margin.(2)Then, Hilbert transform is performed on each intrinsic mode function component to obtain the instantaneous frequency and amplitude of each IMF component c (t) with time. where H is the Hilbert spectrum; PV is the Cauchy principal value.
The analytical signal z (t) is constructed as follows:where a (t) is the amplitude function; φ (t) is the phase function.
Based on formula (7), the instantaneous frequency f (t) is obtained, which is a function of time. It reflects a measure of the frequency concentration of the signal energy at a certain moment, that is, the instantaneous frequency of the signal. Its expression is as follows:
Hilbert transform of each IMF component:
By expressing the above formula as a function of three-dimensional time domain and frequency domain, the Hilbert spectrum can be obtained, that is, the spectrum reflecting the joint distribution of time-frequency-energy.
Hilbert instantaneous energy is defined as follows:
Instantaneous energy reflects the distribution of signal energy over time. The Hilbert energy spectrum can be obtained if the square of the amplitude is integrated with time.
In order to study the distribution law of seismic vibration frequency band energy with respect to frequency band, the frequency band energy analysis of seismic vibration signal can be carried out.where E (ω) is the marginal energy, which describes the distribution of signal energy with frequency. TE is the total energy of the signal based on the Hilbert transform. According to the sampling theorem [18], the Nyquist frequency is 250 Hz, that is, the highest frequency can be identified is 250 Hz since the sampling frequency is 500. In the process of wavelet packet analysis, the number of wavelet packet layers is set to 4, so 16 frequency bands are obtained, which are band 1, band 2,..., band 16, and the band degree of each frequency is 15.625 Hz; the frequency range of the frequency band is 0–15.625 Hz, 15.625–31.25 Hz,...,234.375–250 Hz.
For the convenience of analysis, 125–250 Hz is defined as a high-frequency band, collectively referred to as band 9. Namely, e1, e2, e3,..., e9 are set to the proportion of energy of band 1, band 2,..., band 9, respectively. The calculation formula is as shown in the following formulas:
4.2. Spectrum Analysis
At present, the analysis object of vibration frequency is mainly the Fourier dominant frequency, which is based on the Fourier transform of a whole segment of the signal. The method commonly used to analyze the frequency spectrum of discretely sampled signals is the fast Fourier transform (FFT) algorithm. Now, FFT transformation is performed by MATLAB 2016a programming. Figures 13 and 14 respectively show the time-history curve and frequency spectrum of each monitoring block in the Z direction. At the same time, for the convenience of analysis, the instantaneous Fourier dominant frequency corresponding to the obtained maximum amplitude was extracted, as shown in Figure 14.


(a)

(b)

(c)

(d)

(e)

(f)

(g)

(h)
The energy proportion of the vibration frequency band of each monitoring block in Model 1 is calculated in Table 4.
Figure 13 and Table 4 show that the Fourier dominant frequency of the monitoring block of Model 1 is 9.705 Hz, and the frequency band of each monitoring block is mainly concentrated in the frequency band of 0–31.25 Hz under the action of earthquake, accounting for over 91%, especially the frequency band of 0–15.625 Hz (low-frequency band), accounting for over 77%, indicating that the frequency bands of each rock block under the action of earthquake are basically close, and the frequency band distribution of Model 1 is mainly concentrated in the low-frequency band, and the vibration of earthquake causes the possibility of resonance in Model 1 to be high. The remaining frequency bands (31.25–250 Hz) account for less, and each frequency band accounts for less than 1%, which further indicates that the influence of low-frequency bands on dangerous rocks in Model 1 should be considered under the action of earthquake.
Based on the HHT transformation, the HHT three-dimensional time-frequency diagram (time-frequency-amplitude 3D diagram) of each monitoring block is obtained, as shown in Figure 15.

(a)

(b)

(c)

(d)

(e)

(f)

(g)

(h)
Compare the peak value of time-history curve of vibration velocity of each monitoring block in Figure 13 and the peak value of amplitude of HHT three-dimensional time-frequency diagram of each monitoring block in Figure 15:
From Figure 16, it is found that the peak value of amplitude of each monitoring block is positively correlated with the peak vibration velocity of the monitoring block in Model 1, and the larger the microscopic sequence of the rock block or the higher the vertical elevation of the rock block, the peak value of amplitude of the rock block shows a trend of increase, and the peak value of amplitude of the external B macrochain rock block was larger than that of the internal A chain on the whole, which is consistent with the distribution law of velocity and PGA amplification factor of the monitoring block mentioned above. At the same time, comparing the frequency distribution at different times, it can be seen that the peak values of amplitude are all located within 8–13 s, which is consistent with the distribution law of vibration peak velocity mentioned above, indicating the superiority of HHT time-frequency transformation in dealing with the vibration signals of the rock block of the dangerous rock model.

Same as Model 1, the vibration velocity signals of Model 2 and Model 3 were extracted and FFT transformation was carried out, and the distribution results of Fourier dominant frequency and proportion of frequency band energy of rock blocks in each model are shown in Figure 17.

It can be seen that Model 1 has better integrity under the action of coupled seismic waves, and the Fourier dominant frequencies of each rock block are the same, the Fourier dominant frequency of each rock block on the single chain of Model 2 and Model 3 gradually decreases with the vertical elevation of the rock block, and the Fourier dominant frequency of the B4 rock block on the external B chain of Model 3 is as low as 0.3052 Hz. On the whole, the distribution of Fourier dominant frequency of the external B chain is lower than that of the internal A chain. By contrast, the Fourier dominant frequency of each rock block in Model 3 is smaller than that in Model 2, and Model 1 is the largest. It shows that as the degree of weathering of mudstone caverns deepens, the Fourier dominant frequency of each block decreases gradually under the action of earthquakes. The greater the degree of weathering of mudstone cavity, the faster the attenuation of Fourier dominant frequency of dangerous rock blocks is with the increase of vertical elevation.
From Tables 5 and 6, it can be seen that the distribution of frequency band energy of Model 1, Model 2, and Model 3 under the action of coupled seismic waves is also mainly concentrated in the frequency band of 0–31.25 Hz, and the low-frequency band (0–16.625 Hz) is the main. The energy proportion of low-frequency band (0–16.625 Hz) of each rock block on the single chain of dangerous rocks is positively correlated with the vertical elevation of the rock block, the highest energy proportion of low-frequency band of the B4 rock block of the external B chain in Model 2 and Model 3 is 91.3318% and 99.3673%, respectively, and the proportion of other frequency bands is extremely small. On the whole, the energy proportion of the low-frequency band of the external B chain on the macrochain of each model is higher than that of the internal A chain. Frequency band 2 to frequency band 4 can be collectively referred to as the high-frequency band and then take the average value of the energy distribution of each rock block in each model and draw it as shown in Figure 18. It can be seen that the energy proportion of low-frequency in Model 3 is greater than that in Model 2; Model 1 is the smallest. The distribution rule of high-frequency band is that Model 3 is smaller than that of Model 2, and that of Model 1 is the largest, indicating that under the action of coupled seismic waves, the energy proportion of low-frequency band as a whole is positively correlated with the weathering degree of mudstone cavities and vice versa for high-frequency bands.

5. Conclusions
In this paper, with the help of three-dimensional particle flow software PFC3D, through the establishment of three groups of rock cavity dangerous rock models with different degrees of weathering, the progressive failure process, dynamic response law, and energy distribution characteristics of each dangerous rock under the action of seismic transverse and longitudinal two-way coupled waves are studied. The main conclusions are as follows:(1)The basic (potential) failure modes of dangerous rocks under the action of seismic transverse and longitudinal two-way coupled waves can be divided into tearing-toppling failure (Model 1), slip-rotation failure (Model 2), and slip-rotation-toppling failure (Model 3). The crack propagation speed and the failure process of the dangerous rocks are synchronized with the time-history curve of the input seismic acceleration. The evolution process of the rock force chain in Model 2 and Model 3 can be divided into 4 and 3 stages, respectively. The failure response of the rock force chain and the failure degree of the rock blocks are positively correlated with the weathering degree of mudstone cavity.(2)The displacement of the rock block increases with the increase of the earthquake action time. In Model 2 and Model 3, the evolution process can be divided into 4 and 3 stages, respectively. The greater the weathering depth of the rock cavity, the earlier the rock block caving starts under the action of earthquake, and the larger the peak value of the rock block displacement, and the external chain rock block has priority over the internal chain rock block to instability and failure. On the whole, under the action of earthquake, the peak velocity and the PGA amplification coefficient of the rock block are positively correlated with the degree of weathering of the rock cavity. The peak velocity and the PGA amplification coefficient of the external B macrochain rock block are both greater than that of the internal A macrochain rock block. The peak velocity and the PGA magnification coefficient of the rock blocks both increase with the vertical elevation, showing obvious velocity and PGA elevation magnification effects.(3)Under the combined action of earthquake and rock cavity weathering, the peak velocity and PGA magnification coefficient of the rock particles monitoring point of internal A chain fluctuated greatly, while the rock particle monitoring points of external B chain show an increasing trend with the increase of the vertical elevation of the rock particles, showing an obvious elevation magnification effect. The peak velocity and PGA amplification coefficient of the rock particle monitoring points of external B chain are generally larger than those of the A chain, indicating that the rock particles of external B chain are more affected by the earthquake. On the whole, the seismic dynamic response of rock particles is positively correlated with the weathering depth of the rock cavity.(4)The frequency spectrum analysis of the time-history signal of the velocity of the dangerous rock block shows that, under the action of coupled seismic wave, the Fourier main frequency of each rock block on the single chain of the dangerous rocks in each model decreases gradually with the vertical elevation of the rock block, and the Fourier main frequency of rock block B4 of external B chain in model 3 is as low as 0.3052 Hz. The Fourier frequency distribution of the external B chain is lower than that of the internal A chain. As the weathering degree of mudstone cavity deepens, the Fourier frequency of each block decreases gradually. The greater the weathering degree of mudstone cavity, the faster the attenuation of Fourier frequency of the rock block is with the increase of vertical elevation. The frequency band energy distribution of dangerous rock blocks is mainly concentrated in the 0–31.25 Hz frequency band. The energy proportion of the low-frequency band (0–16.625 Hz) of each rock block on the single chain of the dangerous rock is positively correlated with the vertical elevation of the rock block. The peak value of energy proportion of the low-frequency band in Model 1, Model 2, and Model 3 reaches 77.2934%, 91.3318%, and 99.3673, respectively. The low-frequency energy of the external B chain on the macroscopic chain of each model is higher than that of the internal A chain. The energy proportion of low-frequency of dangerous rocks is positively correlated with the weathering degree of mudstone cavity and vice versa for high-frequency band.
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 they have no conflicts of interest.
Acknowledgments
The authors gratefully acknowledge the support by the National Natural Science Foundation of China(41472262 and 51678097), Chongqing University Innovation Research Group (CXQT19021), Key Project of Chongqing Natural Science Foundation (cstc2020jcyj-zdxmx0012), The First Batch of Chongqing Elite Innovative Leading Talents (CQYC201903026), Science and Technology Research Program of Chongqing Municipal Education Commission (KJQN202001343), and Social Science Program of Chongqing Federation of Social Science Circles (2020QNJJ15).