Failure Prediction of Open-Pit Mine Landslides Containing Complex Geological Structures Using the Inverse Velocity Method
2.1. Study Area
The Fushun West open-pit mine is located in the western part of the Fushun coalfield, at the northern foot of Qiantai Mountain on the southern bank of the Hun River. The geographical coordinates of the mine range from approximately 41°38′0″ N to 42°14′0″ N and 123°39′12″ E to 124°28′0″ E, as shown in Figure 1a. The open pit has a length of approximately 6.6 km, a width of around 2.2 km, and a total area of about 14.52 km2. The mining depth reaches 400 m. However, the Fushun West open-pit mine faces serious landslide hazards due to factors including open-pit mining, underground excavation, faults, and weak layers. Specifically, the mine has experienced over 900 collapse events attributable to landslides. More than 50% of these incidents occurred from June to September when rainfall is relatively concentrated. These landslide events have resulted in a total damaged area of 635,000 m2 and have given rise to a series of safety and geological environmental issues concerning open-pit mining.
On the evening of 25 July 2016, the Fushun area was struck by a rainfall event with a return period of 50 years, resulting in nearly 200 mm of precipitation. At 5:00 a.m. the following day, a partial landslide occurred on the northern slope of the Fushun West open-pit mine, as shown in Figure 1b. The elevation of the landslide’s rear edge was approximately +75 m, while the shear location at the front edge was around −25 m. The landslide spanned a north–south width of approximately 300 m, with a height difference of 110 m, and with an east–west width of approximately 500 m. The total area affected by the landslide was approximately 150,000 m2. The landslide caused the burial of the bottom sections 12 and 14 of the mainline, with the sliding tongue extending. This resulted in the complete interruption of the internal electrical railway lines in the eastern section of the mine, as well as the disruption of the western slope’s transportation roads, including the Xingping Road and the car transport highway. These disruptions had a significant impact on the internal drainage of the eastern open pit and the upper soil removal in the western area of the mine, severely impeding normal production in the mining area.
2.3. Data Description
GPS (Global Positioning System, hereafter called “GPS”) is a huge satellite-based system with global coverage for radio navigation and positioning. GPS technology has been developed rapidly in the application fields of navigation, positioning, precision measurement, etc. In particular, the GPS real-time monitoring system has been widely used in the field of real-time monitoring of landslide deformation in open-pit mines for its real-time nature, and has achieved better results.
In order to reduce the threat of landslide disasters in open-pit mines to national property and people’s lives, the Fushun West open-pit mine introduced a GPS real-time monitoring system, whose framework is shown in Figure 4. The system realizes the 24 h uninterrupted monitoring of geological disaster bodies and the remote automatic transmission of landslide displacement monitoring data, providing effective technical guarantees for the early warning of disaster bodies and the activation of emergency plans.
The system monitoring point deployment and GPS monitoring data are presented in detail below. A total of 11 monitoring profiles were established in different directions on the northern end slope of the Fushun West open-pit mine, Namely, E200, E300, E400, E500, E600, E700, E800, E900, E1000, E1100, and E1200. The engineering geological plan of the landslide is illustrated in Figure 5. As of June 2015, a total of 12 GPS monitoring points were deployed along the monitoring profile, as shown in Figure 5b. Its monitoring technology is mainly based on the main radar sensor transmitting microwaves, using differential aperture radar remote prism monitoring technology. The radar nominal precision is 0.1 mm, the frequency signal is 1575.42 Hz, and the wavelength is about 30~50 cm. These monitoring points were named GN1, GN2, GN3, GN, GN5, GN6, GN7, GN8, GN9, GN10, GN11, and GN12. For this study, data were selected from the period of 14 March 2016 to 31 August 2016.
Let the three-dimensional coordinate information of the landslide monitoring point at a certain time point acquired by remote prism synthetic aperture radar monitoring be (x, y, z). The 3D coordinate value corresponding to the initial moment t0 is (x0, y0, z0), the 3D coordinate value corresponding to any moment tn is (xn, yn, zn), and the cumulative displacement (Δx, Δy, Δz) component of the monitoring point within the moment from t0 to tt is:
The change in total displacement at the monitoring point is:
The cumulative displacement of the landslide is divided into cumulative horizontal displacement and cumulative vertical displacement, where cumulative vertical displacement is the cumulative displacement component in the z-direction and cumulative horizontal displacement Δh is:
The direction of landslide sliding is indicated by the displacement azimuth, which is α:
From the above definition, when α is positive, the displacement direction of the monitoring point is upward; when α is negative, the displacement direction of the monitoring point is downward. From the above derivation process, it can be seen that the size as well as the direction of the deformation of the landslide body is jointly determined by the magnitude of the x, y, and z directions. The deformation at any monitoring point on a landslide can be expressed in terms of the displacement components in the three directions of the monitoring point, or by horizontal displacement, vertical displacement, and displacement azimuth.
2.4. The Basal INV Method
The failure mechanism of slopes is defined as a complete paroxysmal collapse of rock and soil material. By analyzing a multitude of triaxial compression laboratory tests and in situ monitoring research, researchers have discovered that the deformation process of most landslides complies to the progressive characteristics [34,35,36] and three-stage law [37,38,39,40], as shown in Figure 6a. The whole process, from the initial deformation to the eventual failure, representatively comprises three stages: decelerating (green proportion), steady-state deformation (blue proportion), and acceleration deformation (pink proportion). Although the described methods have occasionally been successfully applied to a variety of cases such as man-made walls [30,41], rock and soil specimens [42,43,44], volcanic eruptions [45,46], or tunnels [45,47], these methods are primarily applied to unstable slopes. Hence, landslides are regarded as the principal research objects.
Fukuzono [21] further elaborated the classic three-stage creep theorem by propounding a simpler diagrammatic method (Figure 6b), which could be the most used and simple approach to provide a reasonable estimate of failure time. This method is valid for the tertiary stage. It is noteworthy that the method detects an OOA (onset of acceleration) point, which approximately distinguishes the secondary stage and the tertiary stage. The curve (Figure 6a) is separated into two segments by a demarcation point in the tertiary stage during the acceleration evolution process. (I) After an initial acceleration, Figure 6b displays a dotted line that is approximately parallel to the time axis with the landslide reaching equilibrium state; (II). Meanwhile, Figure 6b likewise displays a line whose value is towards ∞ (i.e., v−1 → 0) as the velocity asymptotically increases.
Several authors successively supplied suggestions and guidelines for proficient usage based on a mathematical generalization of Fukuzono’s solution. Representatively, Voight [22,23] encompassed the prediction of failure behavior and proposed the following equation (Equation (5)).
where Ω is the displacement, dΩ/dt and represents the “velocity” and “acceleration” of Ω, respectively. A and α are two empirical constants that denote characteristics of slope failure; recent investigations revealed that A and α are not independent of each other, varying with several factors comprising kinematic motion patterns [48], versatile types of materia, and macro or micro scales [49]. Consequently, Fukuzono proposed the following equation (Equation (6)) for predicting the failure time by combining the aforementioned equation (Equation (5)) with time:
where tf is the time of failure. This method consists in depicting a tangent line to the curve at an arbitrary point Ʌ1 that tallies to moment t1. The tangent passes across the horizontal axis at moment tc1 (tc1,0). Afterward, the point P1 is plotted vertically above Ʌ1, on a line that passes through Ʌ1 and parallel to the Y axis. The segments of t1Ʌ1 and t1tc1 have an equal displacement from the perspective of geometric shapes. The abovementioned procedure is repeated for another random point Ʌ2. Then, the time of failure tf can be obtained as the abscissa of the intercept of a straight line that passes through P1 and P2 (Figure 6c).
The major drawback of Equation (5) is represented by the necessity of determining two constants A and α. According to closely controlled laboratory conditions and studies by several authors [50,51], α commonly spans over three orders of magnitude. For α = 2, 1 α α > 2, the curve of inverse-velocity has a linear, concave, or convex shape (Figure 6d), respectively. For this condition, Segalini et al. [52], who considered 26 emblematic pre-failure landslide cases, proposed that A inclines to take on extremely low or high values as α deviates from 2. While α appears as intermediate fluctuation, this attribute can be sufficient to sensibly influence prediction results.
To solve this issue, the value of α with the assumption that it is equal to 2 can be generally applied to evaluating the time of failure. In terms of guaranteeing production schedules and staff safety, the assumption of α = 2 is often integrated with the mining industry environment because of its demanding promotion of visual feedback. Thus, Equation (6) is simplified into the following equation (Equation (7)):
As a result, the failure time tf is presumably provided for the point of abscissa of the extrapolated linear inverse velocity trend with the time axis.
2.5. The Moving Average Filtering INV Method Architecture
As formerly stated, the most powerful aspect of the INV method is probably its simplicity, and it is a useful resource in different instances, bypassing the intrinsic restriction for knowing the slope size, state of activity, and types of material. In addition, the tool resource also provides great convenience under many other aspects (e.g., risk assessment and management), if users can count on agile and suitable methods of appraising the state of the monitored circumstance and establish the probability of impending disastrous accidents, a task which is not always achievable because of hardly compensating restrictions as a consequence of measurement errors and random instrumental noise. Correspondingly, we verify two of the utmost prevalent and foolproof smoothing algorithms, i.e., finite impulse response models (also called moving average filter models). Three types of filter models are described below.
1. Short-term simple box filter (SSBF). As each new velocity datum sampling occurs, users can extract the unweighted mean of the antecedent data points through SSBF algorithm processing. This ensures that the alterations of v in the mean are coordinated with the alterations in the data (v) rather than being shifted in time. An example of a v simple equally weighted running mean is the mean over the latest k entries of a data set involving t entries. Let those velocity data points be v1, v2, …, vt. The mean over the latest k velocity data points is represented as SSBFk () and calculated as follows:
when the new velocity datum (SSBFk, next, ) is collected with the invariable sampling width m, the scope from t – k + 2 to t + 1 is considered. A new value vt+1 comes into the sum and the earliest value vt+1 drops out. This simplifies the computations by proceeding with the antecedent mean SSBFk, antecedent: ()
where the moving average cycle (k) of the SSBF model was set to 2 days (k = 2).
2. Long-term simple box filter (LSBF), where the moving average cycle (k) was selected to be 6 days (k = 6).
3. Exponentially weighted moving average (EWMA). Whereas in the short-term simple box filter (SSBF) and long-term simple box filter (LSBF), the past signal processing is weighted equally, the EWMA model is used to assign exponentially decreasing weights over time. The EWMA for a series can be calculated as follows:
where coefficient ξ represents the scale of recursion, a constant smoothing factor between 0 and 1. The smaller ξ is, the stronger the real-time performance of the moving average () is. On the contrary, the larger ξ is, the stronger the ability to absorb instantaneous burst value is, and the better the stability of the prediction model is. Hence, the smoothing factor with the assumption that ξ = 0.5 can be generally used to express and balance the recursion and attenuation properties of the prediction model.
Because the measurement instrument can be easily controlled by the geologist, the time interval between adjoining measurements can be given over a constant time interval. The pattern of filtering as in Equation (8) is equivalent to the easy mathematical statement utilized by Osansan and Stacey [52].
In Equation (11), we set 1/Ω equal to 1/Ωi (1/Ωi is the reciprocal of displacement rate at ti) and t equal to ti (t0 is the most recent instant).
The SSBF and LSBF models are customizable because they can be calculated for different numbers of time cycles (also called the order of the moving average). The biggest distinction of the SSBF and LSBF models is over setting the length of time cycles (k). Significantly, there is no regularly precise guideline or standard definition to set up the boundary between short-term and long-term cycles. It is noted that the selection of a suitable time cycle (k) value is due to the following two major elements: (i) monitoring data accuracy/quality; and (ii) data sampling frequency. To observe what the trend-cycle estimate looks like under different orders of moving average, we plot it (Figure 7) along with a group of monitoring data from an anonymous open-pit mine. It should be noted that the trend cycle (in red) is smoother than the original data after processing several moving averages (n
In Carlà et al. [30], the high-frequency rates of data acquisition, in the area of landslide monitoring programs, are representative of state-of-the-art radar monitoring technology (e.g., GPS [53,54], ground-based radar [31,55], total stations [56], and laser scanning [57]), and commonly require researchers to carry out smoothing over the bulk of measurements. Contrariwise, the low-frequency rates of data acquisition will produce low acquisition rates and will hide much of the background noise, resulting in the inability to trace short-term movements and delaying the identification of eventual trend changes; in such instances, smoothing should be performed over relatively lower measurements, compared to data obtained at high acquisition rates. Short-term averages respond quickly to changes in the price of the underlying security, while long-term averages are slower to react. The order of the moving average determines the smoothness of the trend-cycle estimate. In general, a large order means a smoother curve. The role played by the features of data sampling frequency and quality for the selection of suitable time cycles (k) is notable.
