# Sensitivity Analysis of Factors Affecting Time-Dependent Slope Stability under Freeze-Thaw Cycles

Authors: Zhiguo Chang; Qingxiang Cai; Li Ma; Liu Han
Journal: Mathematical Problems in Engineering
 (2018)

---

## Abstract

With open-pit mines excavated, the slopes will be exposed to the natural environment for a long time. Affected by factors like temperature, seepage, mining, freeze-thaw, etc., slope structural integrity and strength will gradually decline as slope exposure time extends. Besides, the development of defect structure within the rocks is closely correlated with time. In this paper, freeze-thaw cycle tests were conducted on the saturated sandstones collected from a certain open-pit mine. According to the test results, the mass density and longitudinal wave velocity gradually increased with more times of freeze-thaw cycles while mechanical properties such as internal friction angle, cohesion, elastic modulus, and uniaxial compressive strength decreased instead. The constitutive model of saturated rock deterioration was established by taking the volume of phase transition of water in microcracks as a variable. Based on the tests results and theoretical analysis, the sensitivity of the factors affecting slope stability under freeze-thaw damage was studied by Control Variable Method (CVM) and Orthogonal Design Method (ODM). It was determined that the internal friction angle and cohesion had a highly significant effect on the test results, while the mass density had a significant effect. The conclusions may play a certain role in guiding slope construction and protection.


---

## Body

## 1. Introduction

Rock is typical of highly heterogeneous material and there exist various microdefects in different scales within it. And their formation, expansion, and perforation have a great influence on rock deformation and even failure. The development of internal defect structure in rock is closely correlated with time. During the whole formation and evolution of rock engineering, a variety of environmental factors may lead to the gradual deterioration of mechanical properties. For years, many scholars have devoted themselves to the study of time-dependent characteristics of rock strength gradually weakening over time. A. Dragon et al. (1979) studied the brittle and plastic failure of rock according to the fracture surface and established corresponding time-dependent continuum model[1]. M. A. Kaehanov (1982) studied the failure mechanism of rock materials according to the structure and established corresponding damage constitutive model [2]. J. Kemeny (2003) derived a fracture mechanics model by utilizing subcritical crack growth, which resulted in a closed-form solution for joint cohesion as a function of time [3]. In [4], the creep properties of rock were simulated from the point of view that the mechanical properties of materials deteriorate with time. Having considered slope exposure time and service cycle, Q. X. Cai et al. (2008) proposed a theory of time-dependent slope, which has been successfully applied to steep slope mining practice [5]. This basic research shows that the time-dependent characteristics, namely, internal rock structure gradually weakening over time, play a crucial role in transforming microslope damage into macroinstability.Open-pit coal mining is a complex and systematic engineering. Generally speaking, it takes years, even decades of service period, to complete the surface excavation, slope outcropping, and the inner dumping and reclamation. Exposed to natural environment, slopes are easily affected by temperature, seepage, mining, freeze-thaw, and other factors. And as the exposure time extends, their structural integrity and strength gradually decline [6]. For open-pit mines located in seasonally frozen regions, the exposed soil and rock mass freeze in winter and thaw in spring. And such alternation of freezing and thawing brings about radical changes to their physical and mechanical properties. As a result, the long-term stability of slope engineering is directly compromised.The slope stability of open-pit coal mines under normal temperature conditions has been studied by many scholars and lots of results have been obtained [7–10]. However, as to freeze-thaw related slope disaster and its prevention in seasonally frozen regions, quite limited research has been done so far. As a matter of fact, certain geotechnical engineering research has been conducted under the condition of freeze-thaw cycle and it is mainly for the sake of regional development and resources exploitation in permafrost regions of high elevation or high latitude or in polar regions. To study the slope stability of opencast coal mines under freeze-thaw cycle effect, mature theories and methods such as slope stability evaluation and numerical simulation were quoted herein to analyze the variation pattern of physical and mechanical characteristics of soil and rock mass. On this basis, slope parameters and support programs were further optimized. M. G. Ferrick et al. (2010) completed three groups of experiments, with six experiments in each group, to quantify the differences in terms of soil erosion and rill development following a single freeze-thaw cycle [11]. K. M. Skarzynska (1985) conducted a series of experiments on fine-grained soil under freeze-thaw cycle effect. The experiment results showed that the soil structure, density, and porosity change under such effect, which brings about water redistribution [12]. K. D. Eigenbrod (1995) and J. Graham (2011) found that the strength, compressibility, and pore water pressure of the fine-grained soil were affected in the case of soil melting [13, 14]. Besides, by simulating rock freeze-thaw process, L. Wen et al. (2015) studied the changes of mass, volume, uniaxial compressive strength, tensile strength, freeze-thaw coefficient, and weathering degree of saturated hard rock under freeze-thaw cycle effect. In addition, the factor of safety (Fos) before and after going through such cycles was calculated by using strength reduction method [15]. In combination with theoretical research and numerical analysis in [16], a three-field coupling of stress field, temperature field, and water field was carried out with FLAC3D to simulate the influence that freeze-thaw cycles exert on the slope stability of open-pit copper mines. As a matter of fact, the fracture procedure is a whole process of crack initiation, crack propagation, and final fracture. P. Valko et al. (1994) studied the mechanism of crack propagation resulting from frost heave in heterogeneous rocks [17] by using continuum damage mechanics.In this paper, freeze-thaw cycle tests were conducted on sandstone collected from the slope of a certain open-pit mine. And they are designed to study the variation pattern, the damage and deterioration mechanism of physical and mechanical properties. Specifically, on the basis of ODM and FLAC3D, sensitivity analysis was conducted on the main controlling factors of the landslide. Meanwhile, the factors were also compared to determine their order of importance and correlation so as to provide theoretical basis for the safe construction and landslide protection of freeze-thaw slopes.


## 2. Properties of Saturated Rock in the Condition of Freeze-Thaw Cycles

The drilling cores of the sandstone herein were taken from a slope of an open-pit coal mine in Xinjiang province (a seasonally frozen region in northwest China). The samples were processed into standard cylinders of 50 mm diameter and 100 mm height. Considering the local climate and the service period of the open-pit mine, a total of 24 rock samples were divided into four groups, with six samples in each group. The four groups of samples were set to go through 0, 10, 20, and 30 times of freeze-thaw cycles, respectively. At the beginning of the test, a group of samples was randomly selected and soaked for 48 hours at standard atmospheric pressure before measuring their mass density and geometry. The measurement was followed by triaxial compression tests at a constant displacement loading rate of2×10-3mm/s. The confining pressures of the two rock samples were 0 and the other four rock samples were provided with two confining pressure levels. The frictional stress circle was plotted according to the maximum principal stress under different confining pressures, and the corresponding internal friction angle and cohesion were obtained and recorded. At the same time, the remaining samples were placed in a closed system equipped with two sets of independent refrigeration systems, whose control precision and display precision were, respectively, 0.5°C and 0.1°C, as shown in Figures 1(a) and 1(b). A freeze-thaw cycle was set as 24 hours, of which 14 hours was at a low-temperature of -20°C so as to ensure the complete freezing of the samples. During the rest of 10 hours, however, they were put in pure water of +20°C to thaw out. Similarly, the properties tests and triaxial compression tests were performed on the other three groups of samples that had experienced 10, 20, and 30 times of freeze-thaw cycles, as mentioned above. And the testing results were recorded accordingly, as shown in Table 1.Table 1
Properties of sandstone with different times of freeze-thaw cycles.
Freeze-thaw cycles Samplemass/(g) Wave velocity/(m/s) Compressive strength/(Mpa) Young’s modulus/(GPa) Poisson’s ratio Cohesion/(MPa) Angle of internal friction/(°C) 0 times 440.79 5140 138.68 8.3 0.214 28.6 25.6 10 times 440.86 5256 126.6 7.8 0.215 25.07 25.2 20 times 440.95 5441 118.87 7.2 0.214 24.55 24.6 30 times 441.01 5492 115.67 6.7 0.212 23.56 23.6Figure 1
(a) Refrigeration systems; (b) temperature control system.
(a) (b)( 1 ) The tests showed that, with freeze-thaw cycle times increasing, certain sandstone properties began to change. For instance, the uniaxial compressive strength, Young’s modulus, and cohesion obviously decreased as a result. In addition, since the rock samples were saturated when their mass was measured, it is reasonable to conclude that the rock mass increase was caused by water mass increase amongst the cracks. In other words, the space in the microcracks enlarged and the microcracks rate picked up accordingly. It indicated that freeze-thaw cycles resulted in the expansion and coalescence of microcracks in rock samples, which in turn led to internal structural destruction. At the same time, the proportion of effective rock per unit area decreased with gradual microcrack expansion; i.e., the effective stress-bearing area decreased as well, thus lowering the ultimate bearing capacity of the rock.( 2 ) According to the tests, there existed no visible natural spalling after every 10 times of freeze-thaw cycles. However, as can be seen from Table 1, the average mass density of saturated rock samples gradually enlarged with more times of freeze-thaw cycles. Obviously, affected by frost heave action, water freezing led to the microcrack growth. However, in the process of ice melting, the expanded microcrack space was filled with new water, which resulted in the increase of both water content and mass density of rock samples. In addition, the test data showed that the longitudinal wave velocity increased with more times of freeze-thaw cycles. As explained earlier, the closed cracks or unsaturated cracks in the rock sample have to go through frost heave. Considering the fact that wave propagates faster in water than in air, it is reasonable to deduce that it is because the air in closed cracks is replaced by water that not only the mass of the rock samples increased, but also the longitudinal wave velocity increased. They are positively related to each other, as shown in Figure 2.Figure 2
Mass and longitudinal wave velocity under the condition of freeze-thaw cycles.( 3 ) It should be noted that after 30 times of freeze-thaw cycles, the most remarkable changes in properties of sandstone were Young’s modulus, cohesion, and uniaxial compressive strength, with a reduction of 19.3%, 17.6%, and 16.6%, respectively. The internal friction angle increased by 7.8% at maximum and the Poisson's ratio hardly changed.


## 3. Mechanism of Freeze-Thaw Damage for Saturated Rock

The freeze-thaw cycle is actually the phase transition process of water amongst cracks. To be specific, when crack water is frozen into ice, the volume expansion leads to the extension of microcracks. And once it melts into water again, it would in turn fill up the newly expanded cracks. In this sense, long-term freeze-thaw cycles are sure to bring about continuous crack expansion, porosity increase, rock damage accumulation, and, finally, its mechanical strength decrease. With such influence taken into account, a constitutive model of rock deterioration is established herein by taking the phase transition volume of the water in microcracks as a variable.When the temperature drops below freezing point, the water in cracks will be frozen into ice. The microcracks in the rock mass were assumed to be uniformly distributed and the constraints of crack walls on volume expansion were ignored. In this case, the frozen volume of water in microcracks can be expressed as follows:(1)ΔViF=βVwμTwhere ΔViF is the frozen volume of water in microcracks at standard atmospheric pressure, β is the frost heaving coefficient of water at standard atmospheric pressure, Vw is the volume of water in saturated microcracks, and μT is the freezing rate of water atT degrees centigrade.When the water in saturated microcracks is frozen into ice, the crack wall will exert inward compressive stress on the ice while being subjected to outward expansion force. At the same time, with temperature lowering, the rock contracts and the resulting frozen contraction force would impose on the ice as well. Therefore, in addition to the free expansion, the ice volume would experience certain compressive deformation under the action of composite forces mentioned above. It is assumed that rock and ice are isotropic elastic media and the following equation can be obtained:(2)ΔViP=εx+εy+εzVw=1-2viEiσx+σy+σzVwwhere ΔViP is the volume compression of ice due to crack wall load, εx, εy, εz are the strain of ice in the directions ofx-,y- andz-axis, vi is Poisson’s ratio of ice, Ei is Young’s modulus of ice, and σx, σy, σz are the stress of ice in the directions ofx-,y-, andz-axis. Since the expansion force of ice is equal to the frozen contraction force of rock mass in every direction, (2) can be hence expressed as(3)ΔViP=31-2viEiσriVwwhere σri is the composite force exerting on the ice amongst the microcracks in arbitrary directions, σri=σx=σy=σz.So, when the water in the cracks is frozen into ice, the actual volume expansionΔVi is(4)ΔVi=ΔViF-ΔViP=βuT-31-2viEiσriVw.When the ice expansion stress exceeds the rock fracture strengthRt, the volume of microcracks would expand. At the moment of rock cracking, the composite force σri is equal to the fracture strength Rt in terms of their numerical value. So the expanded volume can be expressed as(5)Vi=Vw+ΔVi=1+βuT-31-2viEiRtVwwhere Vi is the microcrack volume after one freeze-thaw cycle.K, the microcrack propagation factor, is defined as(6)K=Vw+ΔViVw=1+ΔViVw=1+βuT-31-2viEiRt.Then afterN times of freeze-thaw cycles,Kn could be expressed as(7)Kn=1+βuT-31-2viEiRtn.For saturated rock mass, as soon as the microcracks within extend, there would be water filling them up in time. Therefore, the volume of microcracks afterN times of freeze-thaw cycles could be expressed as(8)Vin=1+βuT-31-2viEiRtnVw=KnVwwhere Vin is the volume of microcracks within rock mass afterN times of freeze-thaw cycles.Meanwhile, the microcrack rate of rock mass can be expressed as(9)en=VinV0+ΔVin=KnVwVw/e0.+KnVw-Vwwhere V0 is the initial volume of rock, en is the microcrack rate of rock mass afterN times of freeze-thaw cycles, and e0 is the initial microcrack rate. Based on (9), the following can be obtained:(10)en=Kne0.Kn-1e0+1.(11)Kn=en1-e0e01-enOnce the rock mass goes through many times of frost heaving, its microcrack rate changes and such change would directly affect its mechanical strength. A microunit of rock mass was taken as a sample herein and its arbitrary section is composed of microcolumn and microcracks. It showed that after several times of free-thaw cycles, microcracks per unit area extended while the number of microcolumns reduced, which resulted in less bearing area of effective stress. The ultimate bearing capacity of per unit volume before and after freeze-thaw cycles can be expressed as(12)F0=RcΔS0(13)Fn=RcΔSnwhere F0, Fn are the ultimate bearing capacities of the initial state and afterN times of freeze-thaw cycles, respectively. ΔS0, ΔSn are the bearing areas of effective stress in micro-rock-unit and Rc is the compressive strength of microrock column per unit area in microunit. Based on these two equations, we can obtain(14)F0Fn=ΔS0ΔSn.It is assumed that the microrock column is rigid and there exists no extrusion deformation in the condition of frost heave. The microcracks expand and contract with the phase change of crack water. The section area of the micro-rock-unit afterN times of freeze-thaw cycles is expressed as follows:(15)Sn=ΔS0+ΔS01-e0e0Knwhere Sn is the section area of the micro-rock-unit afterN times of freeze-thaw cycles.According to the assumption mentioned above, there is no compression deformation for microrock column when the cracks extend. So the section area of the extended unit will be larger than that of the initial state. By calculation, the effective area of the microrock column per unit area after expansion can be obtained as follows:(16)ΔSn=ΔS0S0Snwhere S0 is the section area of the initial unit.Therefore, the ratio of the effective stress-bearing areas in these two cases is(17)ΔS0ΔSn=SnS0=1+e0Kn-1/1-e0ΔS0ΔS0/1-e0=1+e0Kn-1.In combination with (7), (14), and (17), the relationship of ultimate bearing capacity in micro-rock-unit between the initial state and afterN times of freeze-thaw cycles can be obtained as follows:(18)Fn=F01+e01+βuT-31-2vi/EiRtn.It shows that the ultimate bearing capacity of rock is closely related to the times of freeze-thaw cycles. The rock damage resulting from frost heave is by nature a process of gradual accumulation and deterioration.


## 4. Sensitivity Analysis of Slope Stability Factors

### 4.1. Establishment of Slope Model

Rock is typical of highly heterogeneous material, and its internal microcrack extension is closely correlated with time and the surroundings. When it comes to the rock mass in open-pit slope, its deterioration is a nonlinear cumulative process that is closely related to time. Hence, the slope stability has typical time-dependent feature. During the entire service period ranging from the open-pit excavation, slope exposure, to inner dumping and reclamation, all the external conditions, such as blasting vibration, mining disturbance, weathering, groundwater infiltration, and frost heave, may cause the destruction and accumulation of fine structures inside the rock mass. As time goes on, the mechanical strength of rock mass attenuates and gradually converges to a lower limit value, which results in a gradual reduction of slope stability. Therefore, for the sake of safe construction and slope protection, it is of great significance to conduct sensitivity analysis on factors affecting slope stability.The parameters of the computational model in this paper are 80m, 10m, and 14m for whole height, single bench height, and bench width, respectively. The artificial boundaries are set in places that are far enough from the area of the slope. Specifically, it is 80m at both slope top and bottom boundary and 50m at the lower boundary. The details of modeled geometry are shown in Figure3. Boundary conditions of the model are set that nodes on the boundaries parallel to the X–Y plane, at X=0, X=260m, Y=0, and Y=80m, are constrained in the horizontal direction along X-axis and Y-axis. The nodes are free to move in the horizontal direction along Z-axis while bottom base nodes, at Z=0, are fully constrained. A Mohr-Coulomb material model is initially defined for all zones and the initial stress state is applied with gravity acting in the positive z direction. The model has 5200 zones and 6642 grid-points. The parameters and dimensionless values for each rock mass are shown in Table 2.Table 2
Properties of rock mass.
Stratum Young’s modulus /(GPa) Poisson’s ratio Mass density/(kg/m3) Cohesion/(MPa) Angle of internal friction/(°) Topsoil 1.80×10-2 0.3 1786 13.4 21.3 Sandstone 8.30 0.214 2245 28.6 25.6 Coal 0.65 0.36 1540 13.0 20 Limestone 6.11 0.224 2130 27.0 22.5Figure 3
The details of modeled geometry.Parameters in Table2 were substituted into the model for numerical simulation by using Flac3D. The convergence criterion for calculation is that the unbalanced force ratio falls below the default limiting value of 1×10−5. In this model, the calculation stops at step 7704 and converges to an equilibrium state. The results of the simulation are shown in Figures 4–6.Figure 4
Plastic state and shear-strain increment of the slope model.
(a)
Plastic state (b)
Shear-strain incrementFigure 5
Horizontal and vertical stress of slope.
(a)
Horizontal stress (b)
Vertical stressFigure 6
Displacement of slope in X and Z direction.
(a)
Displacement in X direction (b)
Displacement in Z directionAs shown in Figure4(a), there exist plastic zones within the slope that are parallel to the slope surface, while the basement and back of the slope are generally nondestructive and stable. Figure 4(b) demonstrates the shear-strain increment of the model. It is easy to observe that stress concentration appears at the bottom of the slope, which is close to the boundary of the plastic region. It is also in this area that the shear-strain increment reaches its maximum value. What is more, the stress concentration area and the interface of shear failure interpenetrate, thus bringing about potential landslide damage.According to Figure5(a), on one hand, the horizontal stress is small on the slope surface and gradually increases inward. It reaches its maximum at the bottom of the slope surface, which is not conducive to the slope stability. On the other hand, the vertical stress shown in Figure 5(b) enhances as the depth of the formation increases, indicating that, in vertical direction, the slope body is mainly affected by its own gravity. And such results conform to the variation law of the ground stress.Due to self-gravity, the slope body gradually settles from top to bottom and shear failure occurs at low-intensity places, causing the entire slope to gradually deform toward the empty surface of the slope. From Figures6(a) and 6(b), the maximum displacement in the vertical direction at the top of the slope is 46.3cm, while the maximum horizontal displacement at the bottom of the slope is 64.9cm and the vertical displacement is 9.3cm, indicating that the ground heave occurs at the bottom of the slope due to the compression of the upper landslide body, but no large-scale landslide has yet occurred. The corresponding Fos is 1.6992.


### 4.2. Control Variable Method (CVM) for Single Factor Sensitivity Analysis of Slope Stability

Controlling for a variable is often used in statistics for sensitivity analysis of multivariate problems; that is, when experiments are conducted in an attempt to evaluate the effect of one independent variable on another or more dependent variables, to ensure the measured effect is not influenced by external factors, other variables must be held constant. In this paper, the CVM is used to analyze the effect of a single factor on the slope stability. According to the degree of rock degradation measured in the tests, mass density, Young’s modulus, cohesion, and internal friction angle are taken as independent variables and four levels are designed, respectively, in which the mass density increased by 5% per level, while Young’s modulus, cohesion, and angle of internal friction reduced by 5% per level. The corresponding values of Fos are shown in Table3.Table 3
The value of Fos for slope stability based on CVM.
Level Fos Variables Mass density Young’s modulus Cohesion Angle of internal friction Level 1 1.6602 1.7070 1.6602 1.6602 Level 2 1.6289 1.6992 1.6211 1.6133 Level 3 1.5836 1.6992 1.5742 1.5664 Level 4 1.5566 1.7070 1.5332 1.5215In order to accurately and quantitatively analyze the sensitivity of each factor, the sensitivity coefficients are obtained as shown in Figure7.(19)Si=Δfi/fiΔxi/xiwhere Si is the sensitivity coefficient of each factor, Δfi/fi is the change rate of safety factor, and Δxi/xi is the change rate of the evaluation index.Figure 7
Trend of Fos.It can be seen that the sensitivity of all factors is relatively stable within the variation range of the factors studied. When it comes to the effect on slope stability, they rank like this: internal friction angle > cohesion > mass density > Young’s modulus, as shown in Figure8Figure 8
Histogram of sensitivity coefficient.


### 4.3. Orthogonal Design Method (ODM) for Multifactor Sensitivity Analysis of Slope Stability

In reality, there exist interactions amongst the factors affecting slope stability. Hence, when external conditions change, these influencing factors would change altogether. As can be seen from Figure8, both the internal friction angle and cohesion have a great influence on the slope stability. Instead of probing into these two factors and their possible interactions, the CVM only takes the change of one single factor into consideration. If full factorials are carried out, 44=256 tests should be required according to the factors and levels of this experiment, whose scale will be too large to implement. Therefore, the orthogonal design is used for test design.To explore the influencing factors of slope stability, the orthogonal tests (L1645) were designed and carried out. There were five factors in the tests, each at four levels, thus generating sixteen different combinations. Factors A-D, respectively, represented the mass density, Poisson’s modulus, cohesion, and internal friction angle while factor C∗D stood for the interaction of factors C and D. The details and calculation results are shown in Table 4.Table 4
Orthogonal design and calculation results.
Test No. A B C D C ∗ D Fos 1 1(105%) 1(95%) 1(95%) 1(95%) 1 1.5742 2 1 2(90%) 2(90%) 2(90%) 2 1.4863 3 1 3(85%) 3(85%) 3(85%) 3 1.4043 4 1 4(80%) 4(80%) 4(80%) 4 1.3223 5 2(110%) 1 2 3 4 1.4121 6 2 2 1 4 3 1.4082 7 2 3 4 1 2 1.4199 8 2 4 3 2 1 1.416 9 3(115%) 1 3 4 2 1.2988 10 3 2 4 3 1 1.3066 11 3 3 1 2 4 1.4629 12 3 4 2 1 3 1.4707 13 4(120%) 1 4 2 3 1.3223 14 4 2 3 1 4 1.4043 15 4 3 2 4 1 1.3105 16 4 4 1 3 2 1.3887SetKij as the standard deviation for factorj at leveli, and the difference between the maximum and minimum average standard deviations for the factor can be expressed as(20)Rj=max⁡K1j,K2j,⋯,Kij-min⁡K1j,K2j,⋯,Kij.The calculation results are shown in Table5. Larger range value indicates that the influencing factor makes a greater impact on system performance [18]. The influence of the five factors on slope stability was ranked as follows: D>C>A>B>C∗D. In other words, the internal friction angle and cohesion exerted the greatest influence, followed by mass density. And it could also be seen that the deltaRj of B was obviously smaller than that of D and C, which indicated that the change of Poisson’s modulus had little effect on the slope stability. Besides, the impact resulting from the interaction between internal friction angle and cohesion was much less than that of a single factor.Table 5
Value of range analysis.
Parameter A B C D C ∗ D K 1 j 1.447 1.402 1.474 1.467 1.402 K2j 1.414 1.401 1.42 1.422 1.414 K 3 j 1.385 1.399 1.381 1.394 1.401 K 4 j 1.372 1.415 1.343 1.335 1.4 R j 0.075 0.016 0.131 0.132 0.014  Affect order D>C>A>B>C∗ DFurthermore, theF distribution of the factors was tested by variance analysis. The sums of squared deviations of each factor can be expressed as(21)Sj=rn∑i=1rKij2-1n∑i=1nyi2whereSj is the sum of squared deviations of each factor,yi is the calculation result of testi, andi=1,2,3...,n,r is the level number of one single factor. The degree of freedom ofSj can be calculated as fi=r-1. Since all the columns in the orthogonal table have factors, that is, there are no empty columns herein, the minimum value in theSj is taken as the sum of squared deviations for errors. The statistics can then be defined as follows:(22)F=Sj/fjSe/fe=Sj¯/Se¯.By theF distribution table, F1-0.01 (3,3) = 29.5, F1-0.05 (3,3) = 9.28, and F1-0.1 (3,3) = 5.39 can be obtained. Given that FC=SC¯/Se¯=(0.038/3)/(0.001/3)=38>29.5=F1-0.01 (3,3), FD=SD¯/Se¯=(0.037/3)/(0.001/3)=37>29.5=F1-0.01(3,3), it is considered that the internal friction angle and cohesion have an extremely significant influence on the test results. FA=SA¯/Se¯=(0.013/3)/(0.001/3)=13>9.28=F1-0.053,3; it shows that mass density has a significant influence. Finally, FC∗D=S¯C∗D/Se¯=(0.001/3)/(0.001/3)=1<5.39=F1-0.13,3; it indicates that the effect that the interaction between internal friction angle and cohesion exerts on the test results is not significant. This is in accordance with the previous conclusions obtained by value analysis. The results of the tests are observed in Table 6.Table 6
Variance analysis of various factors.
Source of variance S j f j F = S j / f j S e / f e Significance A 0.013 3 13.000 ∗ ∗ B Δ 0.001 3 1.000 C 0.038 3 38.000 ∗ ∗ ∗ D 0.037 3 37.000 ∗ ∗ ∗ C ∗ D 0.001 3 1.000 Error 0.001 3 ∗ ∗ ∗ means extremely significant influence (α<0.01); ∗∗ means significant influence (0.01<α<0.05); ∗ means general influence (0.05<α<0.1); other values are not significant (0.1<α)


## 4.1. Establishment of Slope Model

Rock is typical of highly heterogeneous material, and its internal microcrack extension is closely correlated with time and the surroundings. When it comes to the rock mass in open-pit slope, its deterioration is a nonlinear cumulative process that is closely related to time. Hence, the slope stability has typical time-dependent feature. During the entire service period ranging from the open-pit excavation, slope exposure, to inner dumping and reclamation, all the external conditions, such as blasting vibration, mining disturbance, weathering, groundwater infiltration, and frost heave, may cause the destruction and accumulation of fine structures inside the rock mass. As time goes on, the mechanical strength of rock mass attenuates and gradually converges to a lower limit value, which results in a gradual reduction of slope stability. Therefore, for the sake of safe construction and slope protection, it is of great significance to conduct sensitivity analysis on factors affecting slope stability.The parameters of the computational model in this paper are 80m, 10m, and 14m for whole height, single bench height, and bench width, respectively. The artificial boundaries are set in places that are far enough from the area of the slope. Specifically, it is 80m at both slope top and bottom boundary and 50m at the lower boundary. The details of modeled geometry are shown in Figure3. Boundary conditions of the model are set that nodes on the boundaries parallel to the X–Y plane, at X=0, X=260m, Y=0, and Y=80m, are constrained in the horizontal direction along X-axis and Y-axis. The nodes are free to move in the horizontal direction along Z-axis while bottom base nodes, at Z=0, are fully constrained. A Mohr-Coulomb material model is initially defined for all zones and the initial stress state is applied with gravity acting in the positive z direction. The model has 5200 zones and 6642 grid-points. The parameters and dimensionless values for each rock mass are shown in Table 2.Table 2
Properties of rock mass.
Stratum Young’s modulus /(GPa) Poisson’s ratio Mass density/(kg/m3) Cohesion/(MPa) Angle of internal friction/(°) Topsoil 1.80×10-2 0.3 1786 13.4 21.3 Sandstone 8.30 0.214 2245 28.6 25.6 Coal 0.65 0.36 1540 13.0 20 Limestone 6.11 0.224 2130 27.0 22.5Figure 3
The details of modeled geometry.Parameters in Table2 were substituted into the model for numerical simulation by using Flac3D. The convergence criterion for calculation is that the unbalanced force ratio falls below the default limiting value of 1×10−5. In this model, the calculation stops at step 7704 and converges to an equilibrium state. The results of the simulation are shown in Figures 4–6.Figure 4
Plastic state and shear-strain increment of the slope model.
(a)
Plastic state (b)
Shear-strain incrementFigure 5
Horizontal and vertical stress of slope.
(a)
Horizontal stress (b)
Vertical stressFigure 6
Displacement of slope in X and Z direction.
(a)
Displacement in X direction (b)
Displacement in Z directionAs shown in Figure4(a), there exist plastic zones within the slope that are parallel to the slope surface, while the basement and back of the slope are generally nondestructive and stable. Figure 4(b) demonstrates the shear-strain increment of the model. It is easy to observe that stress concentration appears at the bottom of the slope, which is close to the boundary of the plastic region. It is also in this area that the shear-strain increment reaches its maximum value. What is more, the stress concentration area and the interface of shear failure interpenetrate, thus bringing about potential landslide damage.According to Figure5(a), on one hand, the horizontal stress is small on the slope surface and gradually increases inward. It reaches its maximum at the bottom of the slope surface, which is not conducive to the slope stability. On the other hand, the vertical stress shown in Figure 5(b) enhances as the depth of the formation increases, indicating that, in vertical direction, the slope body is mainly affected by its own gravity. And such results conform to the variation law of the ground stress.Due to self-gravity, the slope body gradually settles from top to bottom and shear failure occurs at low-intensity places, causing the entire slope to gradually deform toward the empty surface of the slope. From Figures6(a) and 6(b), the maximum displacement in the vertical direction at the top of the slope is 46.3cm, while the maximum horizontal displacement at the bottom of the slope is 64.9cm and the vertical displacement is 9.3cm, indicating that the ground heave occurs at the bottom of the slope due to the compression of the upper landslide body, but no large-scale landslide has yet occurred. The corresponding Fos is 1.6992.


## 4.2. Control Variable Method (CVM) for Single Factor Sensitivity Analysis of Slope Stability

Controlling for a variable is often used in statistics for sensitivity analysis of multivariate problems; that is, when experiments are conducted in an attempt to evaluate the effect of one independent variable on another or more dependent variables, to ensure the measured effect is not influenced by external factors, other variables must be held constant. In this paper, the CVM is used to analyze the effect of a single factor on the slope stability. According to the degree of rock degradation measured in the tests, mass density, Young’s modulus, cohesion, and internal friction angle are taken as independent variables and four levels are designed, respectively, in which the mass density increased by 5% per level, while Young’s modulus, cohesion, and angle of internal friction reduced by 5% per level. The corresponding values of Fos are shown in Table3.Table 3
The value of Fos for slope stability based on CVM.
Level Fos Variables Mass density Young’s modulus Cohesion Angle of internal friction Level 1 1.6602 1.7070 1.6602 1.6602 Level 2 1.6289 1.6992 1.6211 1.6133 Level 3 1.5836 1.6992 1.5742 1.5664 Level 4 1.5566 1.7070 1.5332 1.5215In order to accurately and quantitatively analyze the sensitivity of each factor, the sensitivity coefficients are obtained as shown in Figure7.(19)Si=Δfi/fiΔxi/xiwhere Si is the sensitivity coefficient of each factor, Δfi/fi is the change rate of safety factor, and Δxi/xi is the change rate of the evaluation index.Figure 7
Trend of Fos.It can be seen that the sensitivity of all factors is relatively stable within the variation range of the factors studied. When it comes to the effect on slope stability, they rank like this: internal friction angle > cohesion > mass density > Young’s modulus, as shown in Figure8Figure 8
Histogram of sensitivity coefficient.


## 4.3. Orthogonal Design Method (ODM) for Multifactor Sensitivity Analysis of Slope Stability

In reality, there exist interactions amongst the factors affecting slope stability. Hence, when external conditions change, these influencing factors would change altogether. As can be seen from Figure8, both the internal friction angle and cohesion have a great influence on the slope stability. Instead of probing into these two factors and their possible interactions, the CVM only takes the change of one single factor into consideration. If full factorials are carried out, 44=256 tests should be required according to the factors and levels of this experiment, whose scale will be too large to implement. Therefore, the orthogonal design is used for test design.To explore the influencing factors of slope stability, the orthogonal tests (L1645) were designed and carried out. There were five factors in the tests, each at four levels, thus generating sixteen different combinations. Factors A-D, respectively, represented the mass density, Poisson’s modulus, cohesion, and internal friction angle while factor C∗D stood for the interaction of factors C and D. The details and calculation results are shown in Table 4.Table 4
Orthogonal design and calculation results.
Test No. A B C D C ∗ D Fos 1 1(105%) 1(95%) 1(95%) 1(95%) 1 1.5742 2 1 2(90%) 2(90%) 2(90%) 2 1.4863 3 1 3(85%) 3(85%) 3(85%) 3 1.4043 4 1 4(80%) 4(80%) 4(80%) 4 1.3223 5 2(110%) 1 2 3 4 1.4121 6 2 2 1 4 3 1.4082 7 2 3 4 1 2 1.4199 8 2 4 3 2 1 1.416 9 3(115%) 1 3 4 2 1.2988 10 3 2 4 3 1 1.3066 11 3 3 1 2 4 1.4629 12 3 4 2 1 3 1.4707 13 4(120%) 1 4 2 3 1.3223 14 4 2 3 1 4 1.4043 15 4 3 2 4 1 1.3105 16 4 4 1 3 2 1.3887SetKij as the standard deviation for factorj at leveli, and the difference between the maximum and minimum average standard deviations for the factor can be expressed as(20)Rj=max⁡K1j,K2j,⋯,Kij-min⁡K1j,K2j,⋯,Kij.The calculation results are shown in Table5. Larger range value indicates that the influencing factor makes a greater impact on system performance [18]. The influence of the five factors on slope stability was ranked as follows: D>C>A>B>C∗D. In other words, the internal friction angle and cohesion exerted the greatest influence, followed by mass density. And it could also be seen that the deltaRj of B was obviously smaller than that of D and C, which indicated that the change of Poisson’s modulus had little effect on the slope stability. Besides, the impact resulting from the interaction between internal friction angle and cohesion was much less than that of a single factor.Table 5
Value of range analysis.
Parameter A B C D C ∗ D K 1 j 1.447 1.402 1.474 1.467 1.402 K2j 1.414 1.401 1.42 1.422 1.414 K 3 j 1.385 1.399 1.381 1.394 1.401 K 4 j 1.372 1.415 1.343 1.335 1.4 R j 0.075 0.016 0.131 0.132 0.014  Affect order D>C>A>B>C∗ DFurthermore, theF distribution of the factors was tested by variance analysis. The sums of squared deviations of each factor can be expressed as(21)Sj=rn∑i=1rKij2-1n∑i=1nyi2whereSj is the sum of squared deviations of each factor,yi is the calculation result of testi, andi=1,2,3...,n,r is the level number of one single factor. The degree of freedom ofSj can be calculated as fi=r-1. Since all the columns in the orthogonal table have factors, that is, there are no empty columns herein, the minimum value in theSj is taken as the sum of squared deviations for errors. The statistics can then be defined as follows:(22)F=Sj/fjSe/fe=Sj¯/Se¯.By theF distribution table, F1-0.01 (3,3) = 29.5, F1-0.05 (3,3) = 9.28, and F1-0.1 (3,3) = 5.39 can be obtained. Given that FC=SC¯/Se¯=(0.038/3)/(0.001/3)=38>29.5=F1-0.01 (3,3), FD=SD¯/Se¯=(0.037/3)/(0.001/3)=37>29.5=F1-0.01(3,3), it is considered that the internal friction angle and cohesion have an extremely significant influence on the test results. FA=SA¯/Se¯=(0.013/3)/(0.001/3)=13>9.28=F1-0.053,3; it shows that mass density has a significant influence. Finally, FC∗D=S¯C∗D/Se¯=(0.001/3)/(0.001/3)=1<5.39=F1-0.13,3; it indicates that the effect that the interaction between internal friction angle and cohesion exerts on the test results is not significant. This is in accordance with the previous conclusions obtained by value analysis. The results of the tests are observed in Table 6.Table 6
Variance analysis of various factors.
Source of variance S j f j F = S j / f j S e / f e Significance A 0.013 3 13.000 ∗ ∗ B Δ 0.001 3 1.000 C 0.038 3 38.000 ∗ ∗ ∗ D 0.037 3 37.000 ∗ ∗ ∗ C ∗ D 0.001 3 1.000 Error 0.001 3 ∗ ∗ ∗ means extremely significant influence (α<0.01); ∗∗ means significant influence (0.01<α<0.05); ∗ means general influence (0.05<α<0.1); other values are not significant (0.1<α)


## 5. Conclusions

( 1 ) The freeze-thaw cycle tests of saturated sandstone showed that the properties of rock mass were highly sensitive to increasing freeze-thaw cycles. Of all the indexes, physical indexes such as mass density, microcrack rate, and longitudinal wave velocity increased, while mechanical indexes like internal friction angle, cohesion, Poisson’s modulus, and compressive strength showed a decreasing trend. It is reasonable to conclude that freeze-thaw cycles caused the expansion and coalescence of microcracks within rock samples and the destruction of the internal rock structure and finally reduced the ultimate bearing capacity of rock mass.( 2 ) Considering the microcrack propagation factor, the constitutive model of saturated rock deterioration was presented by taking the volume of phase transition of water in microcracks as a variable. It was verified that the ultimate bearing capacity of rock mass gradually decreased under the freeze-thaw cycles, and the damage of frost heave to rock was a process of gradual accumulation and deterioration.( 3 ) Slope stability deteriorates due to freeze-thaw cycles. To study the sensitivity to such deterioration of influencing factors of slope stability, a slope model was established by Flac3D. Through the analysis of CVM and ODM, it was obtained that the sensitivity of all factors was relatively stable within the variation range of factors studied. Through variance analysis, FC=38>29.5 and FD=37>29.5 when the test level α=0.01 is given, which shows that the internal friction angle and cohesion are the most important determinant affecting the slope stability. Meanwhile, FB=FC∗D=1<5.39 when α=0.1; it indicates that Poisson’s modulus as well as the interaction between internal friction angle and cohesion are not significant. When it comes to the effect on slope stability, internal friction angle ranks the first, followed by cohesion, mass density, and Poisson’s modulus. Therefore, it is advisable to pay more attention to the change of significant factors in actual slope construction and protection.



---

*Source: 7431465-2018-06-04.xml*