# Dynamic Stability Discrimination Method for Concrete Dam under Complex Geological Conditions

Authors: Liaojun Zhang; Tianxiao Ma; Hanyun Zhang; Dongsheng Chen
Journal: Mathematical Problems in Engineering
 (2019)

---

## Abstract

The instability of dams will bring immeasurable personal and property losses to the downstream, so it has always been a trendy topic worthy of investigation. Currently, the rigid body limit equilibrium method is the most commonly used method for the dynamic stability analysis of dams. However, under the action of earthquakes, the instability of the integral dam-foundation system threatens the safety of the dams and is of great concern. In this paper, a stability analysis method that can reflect the complex geological structural forms of dam foundations is proposed in this paper. The advantages are that this method deals with the difficulty in assuming sliding surfaces and the lack of quantitative criteria for the dynamic instability analysis of dams with complex geological structural forms of dam foundations. In addition, through the method, the sliding channels that may appear in the dam foundations can be automatically searched under random earthquake action, and the safety factors of the dynamic instability of dams be quantitatively obtained. Taking a high RCC gravity dam under construction in China as an example, the proposed method is applied to the three-dimensional finite element model of the dam-foundation system of this dam, and then the dynamic stability calculation is carried out. Through this method, the formation process of the dam foundation’s plastic zone and the failure of sliding channels with different strength reduction coefficients are studied on and analyzed detailedly, and the quantitative acquisition of the safety factors is realized. The results show that the method is reasonable and feasible, and helps provide a new idea and method for the dynamic stability analysis of dams.


---

## Body

## 1. Introduction

In the construction of high concrete dams in seismic areas, the dynamic instability has always been one of the key issues that must be resolved. The instability of high dams under earthquakes will bring immeasurable loss of life and property to the downstream. Therefore, it is of great significance to study on the seismic safety of high dams. The earthquake action is characterized by randomness, reciprocity, and persistence. The mechanism of a dynamic analysis problem is completely different from a static one, and there are many influencing factors that make the dynamic stability analysis of high concrete dams rather complicated [1–3].The rigid body limit equilibrium method and the finite element method (fem) are two commonly used methods for the dynamic stability analysis of dams [4–6]. To some extent, the rigid body limit equilibrium method are advantageous because of its simpleness, the rich experience accumulated with time, and the mature safety judgment criteria developed in the engineering application of this method. However, this method still has many limitations. To be more specific, the rigid body limit equilibrium method only considers the equilibrium conditions of the system, providing a static solution that cannot reflect the randomness and duration of the earthquake action. In addition, before the calculation and analysis procedures begin, it is necessary to assume sliding surfaces through this method based on the geological structural forms, which can be not so accountable. Therefore, for the high dam projects with complex geological conditions, the rigid body limit equilibrium method may not represent the actual situation [7, 8]. Tu et al. [9] and Chen [10] believes that the rigid body limit equilibrium method, which does not account for the dynamic deformation coupling of rock masses in the dam abutments or the seismic dynamic effects of the rock mass in the dam foundation, does not completely reflect the actual instability state of the high dams.The main methods for the stability study of dam projects include the overload method and the strength reduction method. At present, as for the dynamic stability analysis of hydraulic structures, the seismic overload method is mainly adopted, which is proposed in the seismic design code for hydraulic structures in China [11]. In recent years, some scholars have studied on the strength reduction method for dynamic stability analysis of engineering structures [12–14] and applied the method to the dynamic stability analysis of slopes [15, 16].Because of the differences in the engineering geological structures of dams, the potential stability failure modes of dams are quite different and the situations can be divided into two kinds. One is that there is a distinct structural weak surface in the dam foundation. Assuming that the weak surface is a sliding one, with the seismic overload method applied, the appearance of displacement mutation at the key points is regarded as the limit state of instability. Because of the simplicity and efficient implementation, it is currently the main method utilized for the stability analysis of dams. The other situation is that there is no obvious structural weak surface in the rock mass, which may not induce an obvious sliding surface in a strong earthquake. And in this situation, the possible instability mode is that a certain weak part of the dam foundation first yields and breaks, a cracking or crushing zone forms, and then a plastic zone through between the dam heel and dam toe gradually develops as the rupture area of the rock mass enlarges. Some scholars believe that the inflection point on the displacement response curve, with the local cracking and sliding of the dam body and foundation taken into consideration, should be used as the criterion for the limit state of instability of the overall system from quantitative to qualitative changes [17]. However, the inflection point needs to be determined intuitively by researchers after a large number of trials.Due to the heterogeneity of materials and the complexity of structural forms, the rock mass in the dam foundation has shown obvious nonlinearity during the deformation process and even the final failure process. The catastrophe theory, which is a branch of the modern nonlinear theory, provides strong support for researching the nonlinear process of instability of the rock mass in the dam foundation. It is mainly utilized to describe the transition process from a continuous gradual state to a sudden change state of a nonlinear system. When solving nonlinear problems, the theory has incomparable advantages in terms of basic principles and specific operational processes. Among the methods based on the catastrophe theory, the cusp catastrophe method is the most widely used because the cusp catastrophe method is one of the most mature and reliable methods for structural stability analysis and the cusp catastrophe model is one of the most commonly used models in practical engineering to be specific; compared with the fold catastrophe model, this model has two control variables, which means that its precision is greatly improved. And it has been successfully used in the stability analysis of slopes [18], of the gravity dam-foundation interface [19], and of the dam-bedrock system [20, 21]. In addition, based on the cusp catastrophe theory, Wu et al. [22] investigated the crack status of concrete dams, and Sun and He et al. [23] studied on the stability of arch dams. Petraitis and Dudgeon et al. [24] analyzed the multiple stable states in marine systems, Yu and Liu et al. [25] researched the stability of interbed for salt cavern gas storage, and Xia et al. [26] studied on the catastrophe instability mechanism of the pillar-roof system. The transition process of the entire system from a stable state to an unstable state can be reflected by the changes at specific key points. Thereby, the quantitative determination of the dynamic stability of dams can be realized, which can help reduce the errors caused by the subjective judgment.In this paper, an adaptive quantitative discrimination method is proposed for researching and analyzing the dynamic stability of dams. And the dynamic time-history analysis method is utilized to simulate the yielding and fracture process of the rock mass and to study on the formation law of the potential sliding surfaces and zones in the dam foundation during the earthquake process. This method can overcome the shortcomings of the rigid body limit equilibrium method, which requires artificial assumption of sliding surfaces and lacks the corresponding quantitative criteria. Besides, it can automatically search for sliding channels under complex foundation conditions and quantitatively give the dynamic stability safety factor of the whole dam-foundation system. The correctness and practicability of the method are verified by the classic soil slope example and the actual gravity dam engineering project. The research results in this paper help provide a new idea and method for the seismic design in related fields.


## 2. Adaptive Quantitative Discrimination Method

### 2.1. Yield Criterion and Strength Reduction Method for Rock Mass Material

#### 2.1.1. Yield Criterion for Rock Mass Material

Under strong earthquakes, the damage of a dam foundation is mainly characterized by the shear and tensile failures. The Mohr–Coulomb yield criterion is adopted for the dam foundation. The formula for the shear failure criterion is(1)fs=σ1−σ3Nφ+2c′Nφ=0,Nφ=1+sinφ′1−sinφ′,where c′ and φ′ are the cohesion and internal friction angle of the rock mass, respectively.The formula for the tensile failure criterion is(2)ft=σt−σ3=0,where σt is the tensile strength of the rock mass.


#### 2.1.2. FEM Combining with Strength Reduction Method

O. C. Zienkiewicz first proposed the concept of the reduction coefficients for the shear strength of geotechnical materials in 1975. The principle of the concept is to increase the safety coefficients of the materials’ strength through reducing the calculated strength gradually until the structural system eventually becomes unstable. The safety coefficient is the ratio of the actual strength of the rock and soil to the reduced strength in the occurrence of instability of the system. The reduction formula is as follows:(3)cf′=c′k,φf′=arctantanφ′k,where c′ and φ′ are the cohesion and internal friction angle of the rock mass before reduction, respectively; cf′, and φf′ are the cohesion and internal friction angle of the rock mass after reduction, respectively; and k is the reduction coefficient, which is greater than 1.0.With the development of the elastoplastic FEM, many scholars have considered combining it with the strength reduction method. Hence, the elastoplastic finite element method combing with the strength reduction method has been developed and greatly improved upon, and it is widely used in the stability studies of slopes and dams.


### 2.2. The Criterion of Dynamic Stability Based on Catastrophe Theory

Different from the static problem, the displacement at a key point of the dam-dam foundation system changes with time during the earthquake, and it is difficult to determine the stability of the dam by the displacement mutation at a certain moment. In this paper, displacement is taken as the main instability investigation target for the cusp catastrophe model. First, the typical control points of the dam-foundation system are selected as the key points for the stability analysis of the catastrophe model. Second, according to the finite element calculation result with a certain reduction coefficientk, the cumulative displacement y along the river direction at the key point is selected as the dynamic antisliding investigation target of the gravity dam. Based on this, a cusp catastrophe model is finally established and with multiple reductions of the strength parameters of the materials in the dam foundation, the dynamic reduction coefficient k and the instability investigation target y are fitted by a quartic polynomial equation as equation (4):(4)y=fk=a0+a1k+a2k2+a3k3+a4k4.The Tschirnhaus transformation is performed on equation (4), and q=a3/4a4,p=k+q are defined and substituted into equation (4) to obtain(5)y=b0+b1p+b2p2+b4p4,where(6)b0b1b2b4=q4−q3q2−q1−4q33q2−2q106q2−3q10010000a4a3a2a1a0.Let(7)p=14b44r,b4>0,−14b44r,b4<0,u=b2b4,v=b1b44,and substitute (7) into (5) to obtain(8)y=14r4+12ar2+br+c,c=b0,u=b2b4,v=b1b44,b4>0.Equation (8) takes r as the state variable and u and v as the control variables. Derive r in equation (8), and let the derivative be 0 to get the equilibrium surface equation:(9)r3+ur+v=0,(10)Δ=4u3+27v2.In the cusp catastrophe model, the equilibrium surface is divided into three regions: the upper, middle, and lower leaves (Figure1). In the upper and lower leaves of the model, if the equilibrium point of the system has a minimum value, it indicates that the state of the system is stable. However, when the equilibrium point is located in the middle part of the model, it corresponds to the critical instability state of the system.Figure 1
Equilibrium surface and bifurcation set of the cusp catastrophe model [27].Based on the cusp catastrophe theory, whenu and v are outside the bifurcation point set, namely, when Δ > 0, the dam-foundation system is stable; when u and v are inside the bifurcation set, namely, when Δ < 0, the system is unstable, and when Δ = 0, the system is in a critical state.


### 2.3. Solution Method of the Adaptive Quantitative Criteria

#### 2.3.1. Solution Process of Adaptive Quantitative Criteria

On the basis of the dynamic finite element reduction method and the catastrophe theory mentioned above, a three-dimensional adaptive quantitative criterion for high concrete dams is proposed in this paper. Assume that the strength reduction is performed for n times, and the strength reduction parameterki is used for the dynamic calculation of the dam-foundation system using the FEM. According to equations (4)∼(10), a cusp catastrophe model is established to obtain a catastrophe coefficient Δi. The Bessel formula is used to fit the relationship between the strength reduction coefficient ki and the Δi value, and the stability safety coefficient is obtained at the critical state. The steps are as follows:(1)
Establish a three-dimensional finite element model of the dam-foundation system that can simulate the complex geological structures and reflect the strength characteristics of the rock mass and the possible weak surfaces in the dam foundation.(2)
The analysis selects a limited number of reduction coefficientski and uses a high-accuracy dynamic time-history analysis to calculate the dynamic response of the dam-foundation system corresponding to each reduction coefficient ki, which is intended to study the yield and fracture laws of the rock masses in the dam foundation during earthquakes as well as the possible opening and closing of sliding surfaces or the formation of potential sliding channels.(3)
Select the cumulative displacementy along the river at the key points such as those near the dam heel and the dam toe with the reduction coefficient ki as the target of instability investigation. Then establish a cusp catastrophe model for calculating the catastrophe coefficient Δi with the reduction coefficient ki, and then plot the k-y curve and the k-Δ curve. The stability safety coefficient k∗ with Δ = 0 is eventually obtained by the Bessel function.(4)
Substitute the calculated stability safety coefficientk∗ for the dynamic finite element calculation and check the stability of the dam-foundation system.Based on the strength reduction method and the cusp catastrophe theory, the adaptive discrimination method can quantitatively determine the antisliding stability safety coefficientk of the dam-foundation system.


#### 2.3.2. Bessel Function and Determination of Critical Safety Coefficient

The change of state of the dam-foundation system from stable to unstable is a process in which the value of Δ gradually changes from greater than zero to zero and finally to less than zero. However, in order to obtain the safety coefficientk∗ at the critical state corresponding to Δ = 0 or Δ ≈ 0, it is necessary to carry out trial calculation many times, which greatly increases the calculation workload.In order to quickly determine the safety coefficient at the critical state, the existing data points are fitted by the Bessel function. Because the Bessel curves have the advantages of high fitting accuracy and the ability to be formulated by equations, they are widely used in the establishment of smooth curve models.According to the known data points fitted to the cubic Bessel parametric formula, the relationship between the strength reduction coefficientk and the value of Δ can be obtained and plotted. The cubic Bessel parametric formula is shown as follows:(11)Bt=P01−t3+3P1t1−t2+3P2t21−t+P3t3,t∈0,1,where P0,P1,P2,andP3 are four control points (Figure 2), which define a cubic Bessel curve in a plane or a three-dimensional space. The curve starts at P0 and ends at P3, and P1 and P2 are two directional control points, which only provide direction reference for the curve; t is the parameter, and Bt is the position of the point on the curve when the parameter equals t.Figure 2
Cubic Bessel curve.Based on this, a custom function is programmed in this paper. The Bessel curve fitting between every two data points is used to obtain the Bessel parametric formula of the curve. And then the relationship curve between the dynamic strength reduction coefficientk and the value of Δ can be drawn. In order to get the position of any point on the drawn curve, substitute the value of Δ into the Bessel equation and determine whether the equation is solvable. And if the equation has a solution, the value of k corresponding to each Δ in the interval is calculated based on the solved t. In particular, the safety coefficient k∗ when Δ = 0 or Δ ≈ 0 can be conveniently determined. In the ADINA software, it is only necessary to check the catastrophe displacement when the reduction coefficient equals k∗, through which the critical safety coefficient can be verified and determined, thus reducing the workload of trial calculation.


### 2.4. Numerical Verification Example of Homogeneous Soil Slope

In order to verify the correctness of the method in the previous section, the homogeneous soil slope analyzed by Dawson et al. before is utilized in this section as an example because it has been verified by many scholars using a variety of methods. As shown in Figure3, it is a homogeneous soil slope with a height of 10.0 m, and the slope angle is 45°. The density of the soil is 20 kN/m3, the cohesion is 12.38 kPa, and the internal friction angle is 20°. According to the limit equilibrium method, the safety coefficient of stability of the soil slope in this example is 1.0.Figure 3
Schematic diagram of homogeneous soil slope.According to Figure3, a finite element analysis model is established, in which the soil material adopts the Mohr–Coulomb model and the stability of the slope is analyzed by the strength reduction method, with the reduction coefficient ranging from 0.5 to 2. The cohesion and internal friction angle of the soil are reduced according to equation (3), and the displacements at the point on the top of the soil slope with different reduction coefficients are calculated. Then, the relationship between the displacement at the point on the top of the slope and the strength reduction coefficient is plotted in Figure 4.Figure 4
Relationship curve between displacement and strength reduction coefficientk.It can be seen from Figure4 that when the strength reduction coefficient changes within the range of 0.5–0.9, the displacement on the top of the slope is relatively stable, indicating that the slope is stable. While the reduction coefficient varies from 0.9 to 1.0, the displacement encounters a sudden increase. Therefore, the catastrophe point lies in the interval when the reduction coefficient varies from 0.9 to 1.0. And the displacement in this interval is analyzed by the cusp catastrophe theory, Δs with different reduction coefficients are calculated (Table 1), and the k-Δ curve is plotted as in Figure 5.Table 1
Catastrophe discrimination results of the homogeneous soil slope with different reduction coefficients.
Reduction coefficient
Displacement at the point on the top of the slope (m)
Δ
Result of catastrophe discrimination
0.9
0.000389467
0.00321
Noncatastrophe
0.95
0.000349479
0.00295
Noncatastrophe
0.96
0.000327267
0.00289
Noncatastrophe
0.97
0.000296854
0.00268
Noncatastrophe
0.98
0.000272011
0.00094
Noncatastrophe
0.99
−0.052985201
−0.00150
Catastrophe
1
−0.255328012
−0.00164
CatastropheFigure 5
Relationship curve between Δ and strength reduction coefficientk.When the strength reduction coefficientk ranges from 0.9 to 0.99, the distributions of the plastic zones of the soil slope are shown in Figure 6. When k = 0.9, the soil at the foot of the slope first yields plastically. And with the increase of the reduction coefficient, the plastic zone gradually develops to the top of the slope. When k = 0.98, the plastic zone has penetrated through, and a sliding channel forms; when k = 0.99, the range of the plastic zone proceeds to expand. Through analyzing the changing process of the plastic zone, the safety coefficient of stability of the soil slope lies between 0.98 and 0.99. Combining the results of the analysis based on the cusp catastrophe theory and the analysis of the changing process of the plastic zone, the safety coefficient of stability calculated by the adaptive quantitative method can be determined and verified to be 0.984, which is very close to the theoretical solution with the limit equilibrium method applied. The results show that the adaptive quantitative method is feasible and reasonable.Figure 6
Distribution diagram of plastic zone with different reduction coefficients. (a)k = 0.9, (b) k = 0.97, (c) k = 0.98, and (d) k = 0.99.
(a)
(b)
(c)
(d)


## 2.1. Yield Criterion and Strength Reduction Method for Rock Mass Material

### 2.1.1. Yield Criterion for Rock Mass Material

Under strong earthquakes, the damage of a dam foundation is mainly characterized by the shear and tensile failures. The Mohr–Coulomb yield criterion is adopted for the dam foundation. The formula for the shear failure criterion is(1)fs=σ1−σ3Nφ+2c′Nφ=0,Nφ=1+sinφ′1−sinφ′,where c′ and φ′ are the cohesion and internal friction angle of the rock mass, respectively.The formula for the tensile failure criterion is(2)ft=σt−σ3=0,where σt is the tensile strength of the rock mass.


### 2.1.2. FEM Combining with Strength Reduction Method

O. C. Zienkiewicz first proposed the concept of the reduction coefficients for the shear strength of geotechnical materials in 1975. The principle of the concept is to increase the safety coefficients of the materials’ strength through reducing the calculated strength gradually until the structural system eventually becomes unstable. The safety coefficient is the ratio of the actual strength of the rock and soil to the reduced strength in the occurrence of instability of the system. The reduction formula is as follows:(3)cf′=c′k,φf′=arctantanφ′k,where c′ and φ′ are the cohesion and internal friction angle of the rock mass before reduction, respectively; cf′, and φf′ are the cohesion and internal friction angle of the rock mass after reduction, respectively; and k is the reduction coefficient, which is greater than 1.0.With the development of the elastoplastic FEM, many scholars have considered combining it with the strength reduction method. Hence, the elastoplastic finite element method combing with the strength reduction method has been developed and greatly improved upon, and it is widely used in the stability studies of slopes and dams.


## 2.1.1. Yield Criterion for Rock Mass Material

Under strong earthquakes, the damage of a dam foundation is mainly characterized by the shear and tensile failures. The Mohr–Coulomb yield criterion is adopted for the dam foundation. The formula for the shear failure criterion is(1)fs=σ1−σ3Nφ+2c′Nφ=0,Nφ=1+sinφ′1−sinφ′,where c′ and φ′ are the cohesion and internal friction angle of the rock mass, respectively.The formula for the tensile failure criterion is(2)ft=σt−σ3=0,where σt is the tensile strength of the rock mass.


## 2.1.2. FEM Combining with Strength Reduction Method

O. C. Zienkiewicz first proposed the concept of the reduction coefficients for the shear strength of geotechnical materials in 1975. The principle of the concept is to increase the safety coefficients of the materials’ strength through reducing the calculated strength gradually until the structural system eventually becomes unstable. The safety coefficient is the ratio of the actual strength of the rock and soil to the reduced strength in the occurrence of instability of the system. The reduction formula is as follows:(3)cf′=c′k,φf′=arctantanφ′k,where c′ and φ′ are the cohesion and internal friction angle of the rock mass before reduction, respectively; cf′, and φf′ are the cohesion and internal friction angle of the rock mass after reduction, respectively; and k is the reduction coefficient, which is greater than 1.0.With the development of the elastoplastic FEM, many scholars have considered combining it with the strength reduction method. Hence, the elastoplastic finite element method combing with the strength reduction method has been developed and greatly improved upon, and it is widely used in the stability studies of slopes and dams.


## 2.2. The Criterion of Dynamic Stability Based on Catastrophe Theory

Different from the static problem, the displacement at a key point of the dam-dam foundation system changes with time during the earthquake, and it is difficult to determine the stability of the dam by the displacement mutation at a certain moment. In this paper, displacement is taken as the main instability investigation target for the cusp catastrophe model. First, the typical control points of the dam-foundation system are selected as the key points for the stability analysis of the catastrophe model. Second, according to the finite element calculation result with a certain reduction coefficientk, the cumulative displacement y along the river direction at the key point is selected as the dynamic antisliding investigation target of the gravity dam. Based on this, a cusp catastrophe model is finally established and with multiple reductions of the strength parameters of the materials in the dam foundation, the dynamic reduction coefficient k and the instability investigation target y are fitted by a quartic polynomial equation as equation (4):(4)y=fk=a0+a1k+a2k2+a3k3+a4k4.The Tschirnhaus transformation is performed on equation (4), and q=a3/4a4,p=k+q are defined and substituted into equation (4) to obtain(5)y=b0+b1p+b2p2+b4p4,where(6)b0b1b2b4=q4−q3q2−q1−4q33q2−2q106q2−3q10010000a4a3a2a1a0.Let(7)p=14b44r,b4>0,−14b44r,b4<0,u=b2b4,v=b1b44,and substitute (7) into (5) to obtain(8)y=14r4+12ar2+br+c,c=b0,u=b2b4,v=b1b44,b4>0.Equation (8) takes r as the state variable and u and v as the control variables. Derive r in equation (8), and let the derivative be 0 to get the equilibrium surface equation:(9)r3+ur+v=0,(10)Δ=4u3+27v2.In the cusp catastrophe model, the equilibrium surface is divided into three regions: the upper, middle, and lower leaves (Figure1). In the upper and lower leaves of the model, if the equilibrium point of the system has a minimum value, it indicates that the state of the system is stable. However, when the equilibrium point is located in the middle part of the model, it corresponds to the critical instability state of the system.Figure 1
Equilibrium surface and bifurcation set of the cusp catastrophe model [27].Based on the cusp catastrophe theory, whenu and v are outside the bifurcation point set, namely, when Δ > 0, the dam-foundation system is stable; when u and v are inside the bifurcation set, namely, when Δ < 0, the system is unstable, and when Δ = 0, the system is in a critical state.


## 2.3. Solution Method of the Adaptive Quantitative Criteria

### 2.3.1. Solution Process of Adaptive Quantitative Criteria

On the basis of the dynamic finite element reduction method and the catastrophe theory mentioned above, a three-dimensional adaptive quantitative criterion for high concrete dams is proposed in this paper. Assume that the strength reduction is performed for n times, and the strength reduction parameterki is used for the dynamic calculation of the dam-foundation system using the FEM. According to equations (4)∼(10), a cusp catastrophe model is established to obtain a catastrophe coefficient Δi. The Bessel formula is used to fit the relationship between the strength reduction coefficient ki and the Δi value, and the stability safety coefficient is obtained at the critical state. The steps are as follows:(1)
Establish a three-dimensional finite element model of the dam-foundation system that can simulate the complex geological structures and reflect the strength characteristics of the rock mass and the possible weak surfaces in the dam foundation.(2)
The analysis selects a limited number of reduction coefficientski and uses a high-accuracy dynamic time-history analysis to calculate the dynamic response of the dam-foundation system corresponding to each reduction coefficient ki, which is intended to study the yield and fracture laws of the rock masses in the dam foundation during earthquakes as well as the possible opening and closing of sliding surfaces or the formation of potential sliding channels.(3)
Select the cumulative displacementy along the river at the key points such as those near the dam heel and the dam toe with the reduction coefficient ki as the target of instability investigation. Then establish a cusp catastrophe model for calculating the catastrophe coefficient Δi with the reduction coefficient ki, and then plot the k-y curve and the k-Δ curve. The stability safety coefficient k∗ with Δ = 0 is eventually obtained by the Bessel function.(4)
Substitute the calculated stability safety coefficientk∗ for the dynamic finite element calculation and check the stability of the dam-foundation system.Based on the strength reduction method and the cusp catastrophe theory, the adaptive discrimination method can quantitatively determine the antisliding stability safety coefficientk of the dam-foundation system.


### 2.3.2. Bessel Function and Determination of Critical Safety Coefficient

The change of state of the dam-foundation system from stable to unstable is a process in which the value of Δ gradually changes from greater than zero to zero and finally to less than zero. However, in order to obtain the safety coefficientk∗ at the critical state corresponding to Δ = 0 or Δ ≈ 0, it is necessary to carry out trial calculation many times, which greatly increases the calculation workload.In order to quickly determine the safety coefficient at the critical state, the existing data points are fitted by the Bessel function. Because the Bessel curves have the advantages of high fitting accuracy and the ability to be formulated by equations, they are widely used in the establishment of smooth curve models.According to the known data points fitted to the cubic Bessel parametric formula, the relationship between the strength reduction coefficientk and the value of Δ can be obtained and plotted. The cubic Bessel parametric formula is shown as follows:(11)Bt=P01−t3+3P1t1−t2+3P2t21−t+P3t3,t∈0,1,where P0,P1,P2,andP3 are four control points (Figure 2), which define a cubic Bessel curve in a plane or a three-dimensional space. The curve starts at P0 and ends at P3, and P1 and P2 are two directional control points, which only provide direction reference for the curve; t is the parameter, and Bt is the position of the point on the curve when the parameter equals t.Figure 2
Cubic Bessel curve.Based on this, a custom function is programmed in this paper. The Bessel curve fitting between every two data points is used to obtain the Bessel parametric formula of the curve. And then the relationship curve between the dynamic strength reduction coefficientk and the value of Δ can be drawn. In order to get the position of any point on the drawn curve, substitute the value of Δ into the Bessel equation and determine whether the equation is solvable. And if the equation has a solution, the value of k corresponding to each Δ in the interval is calculated based on the solved t. In particular, the safety coefficient k∗ when Δ = 0 or Δ ≈ 0 can be conveniently determined. In the ADINA software, it is only necessary to check the catastrophe displacement when the reduction coefficient equals k∗, through which the critical safety coefficient can be verified and determined, thus reducing the workload of trial calculation.


## 2.3.1. Solution Process of Adaptive Quantitative Criteria

On the basis of the dynamic finite element reduction method and the catastrophe theory mentioned above, a three-dimensional adaptive quantitative criterion for high concrete dams is proposed in this paper. Assume that the strength reduction is performed for n times, and the strength reduction parameterki is used for the dynamic calculation of the dam-foundation system using the FEM. According to equations (4)∼(10), a cusp catastrophe model is established to obtain a catastrophe coefficient Δi. The Bessel formula is used to fit the relationship between the strength reduction coefficient ki and the Δi value, and the stability safety coefficient is obtained at the critical state. The steps are as follows:(1)
Establish a three-dimensional finite element model of the dam-foundation system that can simulate the complex geological structures and reflect the strength characteristics of the rock mass and the possible weak surfaces in the dam foundation.(2)
The analysis selects a limited number of reduction coefficientski and uses a high-accuracy dynamic time-history analysis to calculate the dynamic response of the dam-foundation system corresponding to each reduction coefficient ki, which is intended to study the yield and fracture laws of the rock masses in the dam foundation during earthquakes as well as the possible opening and closing of sliding surfaces or the formation of potential sliding channels.(3)
Select the cumulative displacementy along the river at the key points such as those near the dam heel and the dam toe with the reduction coefficient ki as the target of instability investigation. Then establish a cusp catastrophe model for calculating the catastrophe coefficient Δi with the reduction coefficient ki, and then plot the k-y curve and the k-Δ curve. The stability safety coefficient k∗ with Δ = 0 is eventually obtained by the Bessel function.(4)
Substitute the calculated stability safety coefficientk∗ for the dynamic finite element calculation and check the stability of the dam-foundation system.Based on the strength reduction method and the cusp catastrophe theory, the adaptive discrimination method can quantitatively determine the antisliding stability safety coefficientk of the dam-foundation system.


## 2.3.2. Bessel Function and Determination of Critical Safety Coefficient

The change of state of the dam-foundation system from stable to unstable is a process in which the value of Δ gradually changes from greater than zero to zero and finally to less than zero. However, in order to obtain the safety coefficientk∗ at the critical state corresponding to Δ = 0 or Δ ≈ 0, it is necessary to carry out trial calculation many times, which greatly increases the calculation workload.In order to quickly determine the safety coefficient at the critical state, the existing data points are fitted by the Bessel function. Because the Bessel curves have the advantages of high fitting accuracy and the ability to be formulated by equations, they are widely used in the establishment of smooth curve models.According to the known data points fitted to the cubic Bessel parametric formula, the relationship between the strength reduction coefficientk and the value of Δ can be obtained and plotted. The cubic Bessel parametric formula is shown as follows:(11)Bt=P01−t3+3P1t1−t2+3P2t21−t+P3t3,t∈0,1,where P0,P1,P2,andP3 are four control points (Figure 2), which define a cubic Bessel curve in a plane or a three-dimensional space. The curve starts at P0 and ends at P3, and P1 and P2 are two directional control points, which only provide direction reference for the curve; t is the parameter, and Bt is the position of the point on the curve when the parameter equals t.Figure 2
Cubic Bessel curve.Based on this, a custom function is programmed in this paper. The Bessel curve fitting between every two data points is used to obtain the Bessel parametric formula of the curve. And then the relationship curve between the dynamic strength reduction coefficientk and the value of Δ can be drawn. In order to get the position of any point on the drawn curve, substitute the value of Δ into the Bessel equation and determine whether the equation is solvable. And if the equation has a solution, the value of k corresponding to each Δ in the interval is calculated based on the solved t. In particular, the safety coefficient k∗ when Δ = 0 or Δ ≈ 0 can be conveniently determined. In the ADINA software, it is only necessary to check the catastrophe displacement when the reduction coefficient equals k∗, through which the critical safety coefficient can be verified and determined, thus reducing the workload of trial calculation.


## 2.4. Numerical Verification Example of Homogeneous Soil Slope

In order to verify the correctness of the method in the previous section, the homogeneous soil slope analyzed by Dawson et al. before is utilized in this section as an example because it has been verified by many scholars using a variety of methods. As shown in Figure3, it is a homogeneous soil slope with a height of 10.0 m, and the slope angle is 45°. The density of the soil is 20 kN/m3, the cohesion is 12.38 kPa, and the internal friction angle is 20°. According to the limit equilibrium method, the safety coefficient of stability of the soil slope in this example is 1.0.Figure 3
Schematic diagram of homogeneous soil slope.According to Figure3, a finite element analysis model is established, in which the soil material adopts the Mohr–Coulomb model and the stability of the slope is analyzed by the strength reduction method, with the reduction coefficient ranging from 0.5 to 2. The cohesion and internal friction angle of the soil are reduced according to equation (3), and the displacements at the point on the top of the soil slope with different reduction coefficients are calculated. Then, the relationship between the displacement at the point on the top of the slope and the strength reduction coefficient is plotted in Figure 4.Figure 4
Relationship curve between displacement and strength reduction coefficientk.It can be seen from Figure4 that when the strength reduction coefficient changes within the range of 0.5–0.9, the displacement on the top of the slope is relatively stable, indicating that the slope is stable. While the reduction coefficient varies from 0.9 to 1.0, the displacement encounters a sudden increase. Therefore, the catastrophe point lies in the interval when the reduction coefficient varies from 0.9 to 1.0. And the displacement in this interval is analyzed by the cusp catastrophe theory, Δs with different reduction coefficients are calculated (Table 1), and the k-Δ curve is plotted as in Figure 5.Table 1
Catastrophe discrimination results of the homogeneous soil slope with different reduction coefficients.
Reduction coefficient
Displacement at the point on the top of the slope (m)
Δ
Result of catastrophe discrimination
0.9
0.000389467
0.00321
Noncatastrophe
0.95
0.000349479
0.00295
Noncatastrophe
0.96
0.000327267
0.00289
Noncatastrophe
0.97
0.000296854
0.00268
Noncatastrophe
0.98
0.000272011
0.00094
Noncatastrophe
0.99
−0.052985201
−0.00150
Catastrophe
1
−0.255328012
−0.00164
CatastropheFigure 5
Relationship curve between Δ and strength reduction coefficientk.When the strength reduction coefficientk ranges from 0.9 to 0.99, the distributions of the plastic zones of the soil slope are shown in Figure 6. When k = 0.9, the soil at the foot of the slope first yields plastically. And with the increase of the reduction coefficient, the plastic zone gradually develops to the top of the slope. When k = 0.98, the plastic zone has penetrated through, and a sliding channel forms; when k = 0.99, the range of the plastic zone proceeds to expand. Through analyzing the changing process of the plastic zone, the safety coefficient of stability of the soil slope lies between 0.98 and 0.99. Combining the results of the analysis based on the cusp catastrophe theory and the analysis of the changing process of the plastic zone, the safety coefficient of stability calculated by the adaptive quantitative method can be determined and verified to be 0.984, which is very close to the theoretical solution with the limit equilibrium method applied. The results show that the adaptive quantitative method is feasible and reasonable.Figure 6
Distribution diagram of plastic zone with different reduction coefficients. (a)k = 0.9, (b) k = 0.97, (c) k = 0.98, and (d) k = 0.99.
(a)
(b)
(c)
(d)


## 3. Engineering Application

### 3.1. Engineering Situation and Calculation Model

Based on an RCC gravity dam in southwestern China, the feasibility and effectiveness of the proposed method for studying the dynamic stability of dams are discussed in this paper. As shown in Figure7, this dam section is 162 m high, 137.2 m wide at the bottom, and 17 m thick. The dam foundation is mainly composed of granite and sandstone, which presents alternative distribution features, and there is a soft interlayer in the upstream of the dam foundation, making it difficult to determine the sliding surfaces by experience and to carry out the traditional rigid body limit equilibrium analysis. Based on the finite element software of ADINA, taking the dynamic interaction among the dam, foundation, and reservoir water into consideration, a three-dimensional finite element model is established in this paper, which can fully reflect the material characteristics of the rock mass in the dam foundation (shown in Figure 7). And the whole model is divided into 8,856 elements, with a total of 10,847 nodes. The coordinate system of the model takes the transverse direction (pointing to the right bank of the river) as the X axis, the longitudinal direction (pointing downstream) as the Y axis, and the vertical direction (pointing upward) as the Z axis. The Mohr–Coulomb yielding criterion is applied for the analysis of the rock mass in the dam foundation. According to the Code for Seismic Design of Hydraulic Structures of Hydropower Project [11], the standard value of the dynamic elastic modulus of concrete should be increased by 50% compared with that of the static elastic modulus. The material parameters are shown in Table 2.Figure 7
Three-dimensional finite element model of the dam-foundation system.Table 2
Physical and mechanical parameters of concrete and rock mass.
Material
Density
Young’s modulus
Poisson’s ratioν
Mohr–Coulomb parameters
ρ (kg/m3)
E (GPa)
φ
′
c
′ (MPa)
RCC (C15)
2500
23
0.18
—
—
RCC (C20)
2500
25
0.18
—
—
RCC (C25)
2500
26
0.18
—
—
Granite
2650
12
0.21
50.19°
1.50
Sandstone
2625
7
0.24
52.43°
1.70
Of fault
2625
3
0.33
19.29°
0.065The calculated static loads include the weight of the dam, the upstream and downstream hydrostatic pressure, the uplift pressure at the dam foundation, and the sediment pressure. The upstream and downstream water heads are, respectively, 113 m and 33.41 m. The uplift pressure is determined according to the Code for Seismic Design of Hydraulic Structures of Hydropower Project. The floating density of the silt is 6.0 kN/m3, and the internal friction angle is 12°.The basic earthquake intensity at the dam site is VIII degrees, and the seismic category of the water retaining structure is A. The peak value of the seismic horizontal acceleration with a probability exceeding 2% in 100 years is 0.316 g. The vertical peak acceleration takes 2/3 of the horizontal peak acceleration. Taking the site spectrum as the target spectrum to synthesize the horizontal and vertical seismic waves, the calculation time step is 0.01 s, and the ground motion duration is 20 s. The time-history curves of the accelerations are shown in Figure8.Figure 8
The artificial seismic wave. (a) Horizontal direction. (b) Vertical direction.
(a)
(b)The influence of the dynamic water pressure of the reservoir on the seismic response of the dam is considered by applying the Westergaard added mass method. In this paper, the Rayleigh damping method is utilized in order to account for the material damping, and the damping ratio of the dam is 10%.In the case of the foundation model, an appropriate boundary condition for the far end nodes is required to prevent reflection of the outgoing waves back into the system. And the consistent equivalent viscoelastic boundary is adopted in this paper to simulate the ground radiation damping [28–30].


### 3.2. Analysis of Results

Based on the above three-dimensional finite element model, the adaptive discrimination method introduced in Section2 is used for the dynamic stability analysis in this paper and the calculation results with different strength reduction coefficients, k = 2.2, 2.4, 2.6, 2.8, 3.0, 3.2, are selected and analyzed.


#### 3.2.1. Analysis for the Dynamic Stability of the Dam Foundation

Since the seismic load changes with time, the displacement at specific points suddenly becoming larger cannot be considered as a basis for judging whether the dam is stable or not, but the final cumulative displacement after the earthquake can be utilized as a criterion for instability discrimination. For this example, two key points, respectively, located at the dam toe and the dam heel, are selected and the displacement-time curves at the two points are drawn (shown in Figures9 and 10).Figure 9
Displacement-time curves of the point at the dam toe.Figure 10
Displacement-time curves of the point at the dam heel.According to the dynamic calculation results with different reduction coefficientsk, the value of the cumulative displacement, y, at the two points under seismic action is obtained. When the reduction coefficient is ki, take the first i values of k and the displacement y to do the quartic polynomial fitting. And then substitute ki into equations (5)–(10) and calculate the mutation factor Δ with different reduction coefficients.The cumulative displacement curves of the point at the dam toe with the specific reduction coefficients are shown in Figure11. When the reduction coefficient equals 2.2, a quartic polynomial fit is performed based on the results with ks that are smaller than 2.2 (Figures 11(a)), and the corresponding mutation factor Δ = 0.082235031 is calculated. Simiarly, the calculation of the mutation factors proceeds when k, respectively, equals to 2.4, 2.6, 2.8, 3.0, and 3.2. It can be seen in Figures 11(a)–11(e) that the goodness of fit is very close to 1, indicating that the accuracy of the quartic polynomial is sufficient and that the cusp catastrophe model can effectively reflect the stability state of the system. In addition, in order to determine the value of the safety factor at Δ = 0, a Bessel function is used to fit the relationship between the dynamic strength reduction factor k and the mutation factor Δ (see Figure 12). It can be seen from the curve in Figure 13 that as the dynamic strength reduction coefficient increases, the value of Δ at the dam toe gradually decreases. And when the value of Δ is less than 0, the system is considered to be unstable; when the curve is nearly close to the horizontal axis, it indicates that the system is on the verge of the critical state. Last but not the least, by writing a macro command, Δ = 0 is substituted into the Bessel relations of cumulative displacement of the point at the dam toe, and the corresponding k is calculated to be 2.58089.Figure 11
(a) The fitting curve of reduction coefficient against displacement of the point at the dam toe (k = 2.2). (b) The fitting curve of reduction coefficient against displacement of the point at the dam toe (k = 2.4). (c) The fitting curve of reduction coefficient against displacement of the point at the dam toe (k = 2.6). (d) The fitting curve of reduction coefficient displacement of the point at the dam toe (k = 2.8). (e) The fitting curve of reduction coefficient against displacement of the point at the dam toe (k = 3.2).
(a)
(b)
(c)
(d)
(e)Figure 12
Relationship between strength reduction.Figure 13
(a) The fitting curve of reduction coefficient against displacement of the point at the dam heel (k = 2.2). (b) The fitting curve of reduction coefficient against displacement of the point at the dam heel (k = 2.4). (c) The fitting curve of reduction coefficient against displacement of the point at the dam heel (k = 2.6). (d) The fitting curve of reduction coefficient against displacement of the point at the dam heel (k = 2.8). (e) The fitting curve of reduction coefficient against displacement of the point at dam heel (k = 3.2).
(a)
(b)
(c)
(d)
(e)Similarly, the curve of cumulative displacement against reduction coefficient of the point at the dam heel is shown in Figure13. It can be seen that the goodness of fit of the quartic polynomial is rather close to 1 as well, indicating that the cusp catastrophe model can well reflect the stability of the dam. And with the bessel function applied, δ = 0 is substituted into the bessel relations of cumulative displacement of the point at the dam heel, and the corresponding k is calculated to be 2.57759 (Figure 14), which is slightly less than that at the dam toe. Besides, it can be drawn from Figures 11–14 that due to the reciprocating effects of the earthquake, the dynamic instability mechanism of the dam is completely different from the static one.Figure 14
Relationship between strength reduction factor and mutation factor.The results of the catastrophe determination of the points, respectively, at the dam toe and dam heel with different reduction coefficients are shown in Table3.Table 3
Catastrophe discrimination results of the points at the dam toe and heel with different reduction coefficients.
Reduction coefficient
Key points
u
v
Δ
Result of catastrophe discrimination
2.2
Dam heel
−0.183078066
−0.0333325
0.005453211
Noncatastrophe
Dam toe
−0.490153037
−0.1431487
0.082235031
Noncatastrophe
2.4
Dam heel
−0.082455007
−0.0122254
0.001793047
Noncatastrophe
Dam toe
−0.339602824
−0.0862048
0.043978566
Noncatastrophe
2.6
Dam heel
−0.093978032
0.00895964
−0.001152577
Catastrophe
Dam toe
−0.160788464
−0.0219337
−0.003638031
Catastrophe
2.8
Dam heel
−0.209970046
−0.0229496
−0.022807694
Catastrophe
Dam toe
−0.214829545
−0.0217247
−0.026916081
Catastrophe
3.0
Dam heel
−0.118377439
−0.0031959
−0.006359630
Catastrophe
Dam toe
−0.146569898
−0.0122332
−0.008554322
Catastrophe
3.2
Dam heel
−0.185233084
−0.0167385
−0.017857578
Catastrophe
Dam toe
−0.19446602
−0.0162687
−0.022270416
Catastrophe


#### 3.2.2. Sensitivity Analysis to the Strength Parameters of the Dam Foundation

Due to the constraints of engineering survey technology, it is difficult to accurately obtain the material parameters of the dam foundation. Therefore, it is necessary to perform a sensitivity analysis to the strength parameters of the dam foundation.First, the strength parameters of the rock mass at the dam foundation are increased by 20%, and the dynamic stability analysis of the dam-foundation system is conducted to obtain the values of the mutation factor Δ corresponding to each different strength reduction coefficient (Table4). And then the Bessel function is used to determine the reduction coefficient when Δ equals 0, and the relationship curves between the mutation factor Δ and the strength reduction factor k are plotted (Figures 15 and 16). It can be seen from the figures that after increasing the strength parameters at the dam foundation by 20%, the dynamic stability safety factor at the dam toe is equal to 2.58368, which is slightly larger than that with the parameters not increased (2.58089). And the dynamic stability safety factor at the dam heel is 2.58790, which is slightly greater than that with the parameters not increased (2.57759).Table 4
Catastrophe discrimination results of the points at the dam toe and heel with different reduction coefficients (20% increase in the dam foundation strength).
Reduction coefficient
Key points
u
v
Δ
Result of catastrophe discrimination
2.2
Dam heel
0.036515
0.005371
0.000974
Noncatastrophe
Dam toe
0.072672
0.009195
0.003818
Noncatastrophe
2.4
Dam heel
−0.07127
−0.01014
0.001325
Noncatastrophe
Dam toe
0.05249
0.012883
0.005060
Noncatastrophe
2.6
Dam heel
−0.05448
0.001862
−0.000553
Catastrophe
Dam toe
−0.06151
0.000709
−0.000917
Catastrophe
2.8
Dam heel
−0.13847
−0.00122
−0.010581
Catastrophe
Dam toe
−0.16737
0.001011
−0.018727
Catastrophe
3.0
Dam heel
−0.11747
−0.00111
−0.006450
Catastrophe
Dam toe
−0.11705
−0.00115
−0.006379
Catastrophe
3.2
Dam heel
−0.18426
0.00634
−0.023939
Catastrophe
Dam toe
−0.15204
0.00152
−0.013995
CatastropheFigure 15
Relationship between strength reduction factor and mutation factor at the dam toe (20% increase in the dam foundation strength).Figure 16
Relationship between strength reduction factor and mutation factor at the dam heel (20% increase in the dam foundation strength).Through the above analysis, it can be drawn that when the strength parameters at the dam foundation are increased by 20%, the dynamic stability safety factors become slightly larger, but the change is not obvious.Similar to the above procedures, the strength parameters at the dam foundation are decreased by 20%, and the values of the mutation factor Δ corresponding to each different strength reduction coefficient are obtained (Table5). And then the Bessel function is used to determine the reduction coefficient when Δ equals 0, and the relationship curves between the mutation factor Δ and the strength reduction factor k are plotted (Figures 17 and 18). It can be seen from the figures that after decreasing the strength parameters at the dam foundation by 20%, the dynamic stability safety factor at the dam toe is equal to 2.56634, which is slightly smaller than that with the parameters not increased (2.58089). And the dynamic stability safety factor at the dam heel is 2.55017, which is slightly smaller than that with the parameters not increased (2.57759).Table 5
Catastrophe discrimination results of the points at the dam toe and heel with different reduction coefficients (20% decrease in the dam foundation strength).
Reduction coefficient
Key points
u
v
Δ
Result of catastrophe discrimination
2.2
Dam heel
0.023094
0.005313
0.000811
Noncatastrophe
Dam toe
0.047434
0.011929
0.004269
Noncatastrophe
2.4
Dam heel
0.05625
0.008574
0.002697
Noncatastrophe
Dam toe
0.069669
0.020231
0.012404
Noncatastrophe
2.6
Dam heel
−0.10016
0.009008
−0.001829
Catastrophe
Dam toe
−0.08765
−0.0005
−0.002687
Catastrophe
2.8
Dam heel
−0.15368
−0.00321
−0.014241
Catastrophe
Dam toe
−−0.18609
−0.00098
−0.025751
Catastrophe
3.0
Dam heel
−0.13073
0.00042
−0.008933
Catastrophe
Dam toe
−0.12937
−−0.00341
−0.008347
Catastrophe
3.2
Dam heel
−0.15189
−0.00204
−0.013905
Catastrophe
Dam toe
−0.15459
0.002175
−0.014649
CatastropheFigure 17
Relationship between strength reduction factor and mutation factor at the dam toe (20% decrease in the dam foundation strength).Figure 18
Relationship between strength reduction factor and mutation factor at the dam heel (20% decrease in the dam foundation strength).Through the above analysis, it can be drawn that when the strength parameters at the dam foundation are decreased by 20%, the dynamic stability safety factors become slightly smaller, but the change is not obvious.


#### 3.2.3. Determination of the Sliding Area in the Dam Foundation

As for the influence of the weak layer orientation on the stability analysis in this paper, the relevant engineering practice and experimental research results show that there may be two types of instability failure modes for gravity dams on rock foundations: one type is that there are inclined weak layers at the rock mass in the dam foundation. And with the weak layers being subject to external loads, sliding blocks may possibly form at the dam foundation, which will result in deep sliding (Figure19). The other is that there are no obvious weak layers in the dam foundation, making it difficult to assume the possible sliding direction of the dam foundation. And due to the reciprocating action of earthquakes, the cracking and crushing occur at the rock mass in the dam foundation, and the yielding area gradually expands until the final penetration. Finally, the dam body and some part of the foundation suffer instability and failure as a whole, as shown in Figure 20.Figure 19
Failure mode with obvious sliding surface.Figure 20
Failure mode without obvious sliding surface.As for the project in this paper, the rock mass in the dam foundation presents a layered distribution approximately perpendicular to the horizontal plane, with similar material parameters. Besides, there is no obvious weak layer. Therefore, the instability mode of the project belongs to the second type. According to the results through finite element calculation, the rock mass near the upstream dam heel and the downstream dam toe yielded first. And with the increase of the strength reduction factor, the yield area gradually expanded until it was completely penetrated, eventually causing the instability failure of the dam body and some part of the dam foundation as a whole.Due to the complex distribution features of the rock mass in the dam foundation, it is difficult to delineate the possible sliding area by experience. Therefore, the capability of searching for the sliding area for the dynamic stability analysis of the dam is of great concern. With the increase of the strength reduction coefficient, the rock mass in the dam foundation gradually yields, and the plastic zones develop with time until final penetration, which is an important sign of the dam instability. In other words, if the plastic zones in the dam foundation have not been run-through, then even if the displacement at the key point of the dam changes abruptly, the system is not necessarily considered to be unstable. Therefore, when applying the above method, the dynamic changes of the plastic zones with different reduction coefficients need to be considered comprehensively to analyze the dynamic stability of the dam.When the reduction coefficient varies from 1.0 to 2.0, the soft interlayer, which tends to the upstream of the dam, first enters the plastic yielding stage, while the plastic zones gradually and continuously develop at the dam heel, dam toe, and other weak parts. When the reduction coefficient equals 2.4, most of the rock mass in the dam foundation is in the plastic yielding state, with the plastic zones not being entirely run-through. And when the reduction coefficient is 2.6, the rock mass in the dam foundation completely enters the plastic yielding state, and the plastic zones are entirely run-through (shown in Figure21). As described in Section 3.2.1, the displacements at the dam heel and the dam toe are confronted with mutations, with the value of δ less than zero. And according to the results combining the analyses of cumulative displacement mutations and of the dynamic distribution of the plastic zone, it can be concluded that the dam-foundation system is in an instability state. The sliding area searched by the adaptive quantitative method is the plastic zone in the dam foundation in Figure 21(d). This method can solve the problem that the sliding area is difficult to determine by experience in the case of complex geological forms of the dam foundation.Figure 21
Distribution diagram of the plastic zone in the foundation with different reduction coefficients. (a)k = 1.0, (b) k = 2.0, (c) k = 2.4, and (d) k = 2.6.
(a)
(b)
(c)
(d)In addition, because the duration of the earthquake is much shorter than the change time of the seepage field and the earthquake action is reciprocating, the influence of the change of the seepage field is generally not considered during the earthquake action. However, after the earthquake, if only a small part of the dam heel yields and the yield area does not extend to the curtain position, the seepage field will have little effect on the stability of the dam. However, if the yield area at the dam foundation extends to the location of the curtain, possibly leading to the damage of the curtain, the osmotic pressure at the dam foundation will gradually increase, which is detrimental to the stability of the dam and requires danger-reinforcement treatment.


## 3.1. Engineering Situation and Calculation Model

Based on an RCC gravity dam in southwestern China, the feasibility and effectiveness of the proposed method for studying the dynamic stability of dams are discussed in this paper. As shown in Figure7, this dam section is 162 m high, 137.2 m wide at the bottom, and 17 m thick. The dam foundation is mainly composed of granite and sandstone, which presents alternative distribution features, and there is a soft interlayer in the upstream of the dam foundation, making it difficult to determine the sliding surfaces by experience and to carry out the traditional rigid body limit equilibrium analysis. Based on the finite element software of ADINA, taking the dynamic interaction among the dam, foundation, and reservoir water into consideration, a three-dimensional finite element model is established in this paper, which can fully reflect the material characteristics of the rock mass in the dam foundation (shown in Figure 7). And the whole model is divided into 8,856 elements, with a total of 10,847 nodes. The coordinate system of the model takes the transverse direction (pointing to the right bank of the river) as the X axis, the longitudinal direction (pointing downstream) as the Y axis, and the vertical direction (pointing upward) as the Z axis. The Mohr–Coulomb yielding criterion is applied for the analysis of the rock mass in the dam foundation. According to the Code for Seismic Design of Hydraulic Structures of Hydropower Project [11], the standard value of the dynamic elastic modulus of concrete should be increased by 50% compared with that of the static elastic modulus. The material parameters are shown in Table 2.Figure 7
Three-dimensional finite element model of the dam-foundation system.Table 2
Physical and mechanical parameters of concrete and rock mass.
Material
Density
Young’s modulus
Poisson’s ratioν
Mohr–Coulomb parameters
ρ (kg/m3)
E (GPa)
φ
′
c
′ (MPa)
RCC (C15)
2500
23
0.18
—
—
RCC (C20)
2500
25
0.18
—
—
RCC (C25)
2500
26
0.18
—
—
Granite
2650
12
0.21
50.19°
1.50
Sandstone
2625
7
0.24
52.43°
1.70
Of fault
2625
3
0.33
19.29°
0.065The calculated static loads include the weight of the dam, the upstream and downstream hydrostatic pressure, the uplift pressure at the dam foundation, and the sediment pressure. The upstream and downstream water heads are, respectively, 113 m and 33.41 m. The uplift pressure is determined according to the Code for Seismic Design of Hydraulic Structures of Hydropower Project. The floating density of the silt is 6.0 kN/m3, and the internal friction angle is 12°.The basic earthquake intensity at the dam site is VIII degrees, and the seismic category of the water retaining structure is A. The peak value of the seismic horizontal acceleration with a probability exceeding 2% in 100 years is 0.316 g. The vertical peak acceleration takes 2/3 of the horizontal peak acceleration. Taking the site spectrum as the target spectrum to synthesize the horizontal and vertical seismic waves, the calculation time step is 0.01 s, and the ground motion duration is 20 s. The time-history curves of the accelerations are shown in Figure8.Figure 8
The artificial seismic wave. (a) Horizontal direction. (b) Vertical direction.
(a)
(b)The influence of the dynamic water pressure of the reservoir on the seismic response of the dam is considered by applying the Westergaard added mass method. In this paper, the Rayleigh damping method is utilized in order to account for the material damping, and the damping ratio of the dam is 10%.In the case of the foundation model, an appropriate boundary condition for the far end nodes is required to prevent reflection of the outgoing waves back into the system. And the consistent equivalent viscoelastic boundary is adopted in this paper to simulate the ground radiation damping [28–30].


## 3.2. Analysis of Results

Based on the above three-dimensional finite element model, the adaptive discrimination method introduced in Section2 is used for the dynamic stability analysis in this paper and the calculation results with different strength reduction coefficients, k = 2.2, 2.4, 2.6, 2.8, 3.0, 3.2, are selected and analyzed.


### 3.2.1. Analysis for the Dynamic Stability of the Dam Foundation

Since the seismic load changes with time, the displacement at specific points suddenly becoming larger cannot be considered as a basis for judging whether the dam is stable or not, but the final cumulative displacement after the earthquake can be utilized as a criterion for instability discrimination. For this example, two key points, respectively, located at the dam toe and the dam heel, are selected and the displacement-time curves at the two points are drawn (shown in Figures9 and 10).Figure 9
Displacement-time curves of the point at the dam toe.Figure 10
Displacement-time curves of the point at the dam heel.According to the dynamic calculation results with different reduction coefficientsk, the value of the cumulative displacement, y, at the two points under seismic action is obtained. When the reduction coefficient is ki, take the first i values of k and the displacement y to do the quartic polynomial fitting. And then substitute ki into equations (5)–(10) and calculate the mutation factor Δ with different reduction coefficients.The cumulative displacement curves of the point at the dam toe with the specific reduction coefficients are shown in Figure11. When the reduction coefficient equals 2.2, a quartic polynomial fit is performed based on the results with ks that are smaller than 2.2 (Figures 11(a)), and the corresponding mutation factor Δ = 0.082235031 is calculated. Simiarly, the calculation of the mutation factors proceeds when k, respectively, equals to 2.4, 2.6, 2.8, 3.0, and 3.2. It can be seen in Figures 11(a)–11(e) that the goodness of fit is very close to 1, indicating that the accuracy of the quartic polynomial is sufficient and that the cusp catastrophe model can effectively reflect the stability state of the system. In addition, in order to determine the value of the safety factor at Δ = 0, a Bessel function is used to fit the relationship between the dynamic strength reduction factor k and the mutation factor Δ (see Figure 12). It can be seen from the curve in Figure 13 that as the dynamic strength reduction coefficient increases, the value of Δ at the dam toe gradually decreases. And when the value of Δ is less than 0, the system is considered to be unstable; when the curve is nearly close to the horizontal axis, it indicates that the system is on the verge of the critical state. Last but not the least, by writing a macro command, Δ = 0 is substituted into the Bessel relations of cumulative displacement of the point at the dam toe, and the corresponding k is calculated to be 2.58089.Figure 11
(a) The fitting curve of reduction coefficient against displacement of the point at the dam toe (k = 2.2). (b) The fitting curve of reduction coefficient against displacement of the point at the dam toe (k = 2.4). (c) The fitting curve of reduction coefficient against displacement of the point at the dam toe (k = 2.6). (d) The fitting curve of reduction coefficient displacement of the point at the dam toe (k = 2.8). (e) The fitting curve of reduction coefficient against displacement of the point at the dam toe (k = 3.2).
(a)
(b)
(c)
(d)
(e)Figure 12
Relationship between strength reduction.Figure 13
(a) The fitting curve of reduction coefficient against displacement of the point at the dam heel (k = 2.2). (b) The fitting curve of reduction coefficient against displacement of the point at the dam heel (k = 2.4). (c) The fitting curve of reduction coefficient against displacement of the point at the dam heel (k = 2.6). (d) The fitting curve of reduction coefficient against displacement of the point at the dam heel (k = 2.8). (e) The fitting curve of reduction coefficient against displacement of the point at dam heel (k = 3.2).
(a)
(b)
(c)
(d)
(e)Similarly, the curve of cumulative displacement against reduction coefficient of the point at the dam heel is shown in Figure13. It can be seen that the goodness of fit of the quartic polynomial is rather close to 1 as well, indicating that the cusp catastrophe model can well reflect the stability of the dam. And with the bessel function applied, δ = 0 is substituted into the bessel relations of cumulative displacement of the point at the dam heel, and the corresponding k is calculated to be 2.57759 (Figure 14), which is slightly less than that at the dam toe. Besides, it can be drawn from Figures 11–14 that due to the reciprocating effects of the earthquake, the dynamic instability mechanism of the dam is completely different from the static one.Figure 14
Relationship between strength reduction factor and mutation factor.The results of the catastrophe determination of the points, respectively, at the dam toe and dam heel with different reduction coefficients are shown in Table3.Table 3
Catastrophe discrimination results of the points at the dam toe and heel with different reduction coefficients.
Reduction coefficient
Key points
u
v
Δ
Result of catastrophe discrimination
2.2
Dam heel
−0.183078066
−0.0333325
0.005453211
Noncatastrophe
Dam toe
−0.490153037
−0.1431487
0.082235031
Noncatastrophe
2.4
Dam heel
−0.082455007
−0.0122254
0.001793047
Noncatastrophe
Dam toe
−0.339602824
−0.0862048
0.043978566
Noncatastrophe
2.6
Dam heel
−0.093978032
0.00895964
−0.001152577
Catastrophe
Dam toe
−0.160788464
−0.0219337
−0.003638031
Catastrophe
2.8
Dam heel
−0.209970046
−0.0229496
−0.022807694
Catastrophe
Dam toe
−0.214829545
−0.0217247
−0.026916081
Catastrophe
3.0
Dam heel
−0.118377439
−0.0031959
−0.006359630
Catastrophe
Dam toe
−0.146569898
−0.0122332
−0.008554322
Catastrophe
3.2
Dam heel
−0.185233084
−0.0167385
−0.017857578
Catastrophe
Dam toe
−0.19446602
−0.0162687
−0.022270416
Catastrophe


### 3.2.2. Sensitivity Analysis to the Strength Parameters of the Dam Foundation

Due to the constraints of engineering survey technology, it is difficult to accurately obtain the material parameters of the dam foundation. Therefore, it is necessary to perform a sensitivity analysis to the strength parameters of the dam foundation.First, the strength parameters of the rock mass at the dam foundation are increased by 20%, and the dynamic stability analysis of the dam-foundation system is conducted to obtain the values of the mutation factor Δ corresponding to each different strength reduction coefficient (Table4). And then the Bessel function is used to determine the reduction coefficient when Δ equals 0, and the relationship curves between the mutation factor Δ and the strength reduction factor k are plotted (Figures 15 and 16). It can be seen from the figures that after increasing the strength parameters at the dam foundation by 20%, the dynamic stability safety factor at the dam toe is equal to 2.58368, which is slightly larger than that with the parameters not increased (2.58089). And the dynamic stability safety factor at the dam heel is 2.58790, which is slightly greater than that with the parameters not increased (2.57759).Table 4
Catastrophe discrimination results of the points at the dam toe and heel with different reduction coefficients (20% increase in the dam foundation strength).
Reduction coefficient
Key points
u
v
Δ
Result of catastrophe discrimination
2.2
Dam heel
0.036515
0.005371
0.000974
Noncatastrophe
Dam toe
0.072672
0.009195
0.003818
Noncatastrophe
2.4
Dam heel
−0.07127
−0.01014
0.001325
Noncatastrophe
Dam toe
0.05249
0.012883
0.005060
Noncatastrophe
2.6
Dam heel
−0.05448
0.001862
−0.000553
Catastrophe
Dam toe
−0.06151
0.000709
−0.000917
Catastrophe
2.8
Dam heel
−0.13847
−0.00122
−0.010581
Catastrophe
Dam toe
−0.16737
0.001011
−0.018727
Catastrophe
3.0
Dam heel
−0.11747
−0.00111
−0.006450
Catastrophe
Dam toe
−0.11705
−0.00115
−0.006379
Catastrophe
3.2
Dam heel
−0.18426
0.00634
−0.023939
Catastrophe
Dam toe
−0.15204
0.00152
−0.013995
CatastropheFigure 15
Relationship between strength reduction factor and mutation factor at the dam toe (20% increase in the dam foundation strength).Figure 16
Relationship between strength reduction factor and mutation factor at the dam heel (20% increase in the dam foundation strength).Through the above analysis, it can be drawn that when the strength parameters at the dam foundation are increased by 20%, the dynamic stability safety factors become slightly larger, but the change is not obvious.Similar to the above procedures, the strength parameters at the dam foundation are decreased by 20%, and the values of the mutation factor Δ corresponding to each different strength reduction coefficient are obtained (Table5). And then the Bessel function is used to determine the reduction coefficient when Δ equals 0, and the relationship curves between the mutation factor Δ and the strength reduction factor k are plotted (Figures 17 and 18). It can be seen from the figures that after decreasing the strength parameters at the dam foundation by 20%, the dynamic stability safety factor at the dam toe is equal to 2.56634, which is slightly smaller than that with the parameters not increased (2.58089). And the dynamic stability safety factor at the dam heel is 2.55017, which is slightly smaller than that with the parameters not increased (2.57759).Table 5
Catastrophe discrimination results of the points at the dam toe and heel with different reduction coefficients (20% decrease in the dam foundation strength).
Reduction coefficient
Key points
u
v
Δ
Result of catastrophe discrimination
2.2
Dam heel
0.023094
0.005313
0.000811
Noncatastrophe
Dam toe
0.047434
0.011929
0.004269
Noncatastrophe
2.4
Dam heel
0.05625
0.008574
0.002697
Noncatastrophe
Dam toe
0.069669
0.020231
0.012404
Noncatastrophe
2.6
Dam heel
−0.10016
0.009008
−0.001829
Catastrophe
Dam toe
−0.08765
−0.0005
−0.002687
Catastrophe
2.8
Dam heel
−0.15368
−0.00321
−0.014241
Catastrophe
Dam toe
−−0.18609
−0.00098
−0.025751
Catastrophe
3.0
Dam heel
−0.13073
0.00042
−0.008933
Catastrophe
Dam toe
−0.12937
−−0.00341
−0.008347
Catastrophe
3.2
Dam heel
−0.15189
−0.00204
−0.013905
Catastrophe
Dam toe
−0.15459
0.002175
−0.014649
CatastropheFigure 17
Relationship between strength reduction factor and mutation factor at the dam toe (20% decrease in the dam foundation strength).Figure 18
Relationship between strength reduction factor and mutation factor at the dam heel (20% decrease in the dam foundation strength).Through the above analysis, it can be drawn that when the strength parameters at the dam foundation are decreased by 20%, the dynamic stability safety factors become slightly smaller, but the change is not obvious.


### 3.2.3. Determination of the Sliding Area in the Dam Foundation

As for the influence of the weak layer orientation on the stability analysis in this paper, the relevant engineering practice and experimental research results show that there may be two types of instability failure modes for gravity dams on rock foundations: one type is that there are inclined weak layers at the rock mass in the dam foundation. And with the weak layers being subject to external loads, sliding blocks may possibly form at the dam foundation, which will result in deep sliding (Figure19). The other is that there are no obvious weak layers in the dam foundation, making it difficult to assume the possible sliding direction of the dam foundation. And due to the reciprocating action of earthquakes, the cracking and crushing occur at the rock mass in the dam foundation, and the yielding area gradually expands until the final penetration. Finally, the dam body and some part of the foundation suffer instability and failure as a whole, as shown in Figure 20.Figure 19
Failure mode with obvious sliding surface.Figure 20
Failure mode without obvious sliding surface.As for the project in this paper, the rock mass in the dam foundation presents a layered distribution approximately perpendicular to the horizontal plane, with similar material parameters. Besides, there is no obvious weak layer. Therefore, the instability mode of the project belongs to the second type. According to the results through finite element calculation, the rock mass near the upstream dam heel and the downstream dam toe yielded first. And with the increase of the strength reduction factor, the yield area gradually expanded until it was completely penetrated, eventually causing the instability failure of the dam body and some part of the dam foundation as a whole.Due to the complex distribution features of the rock mass in the dam foundation, it is difficult to delineate the possible sliding area by experience. Therefore, the capability of searching for the sliding area for the dynamic stability analysis of the dam is of great concern. With the increase of the strength reduction coefficient, the rock mass in the dam foundation gradually yields, and the plastic zones develop with time until final penetration, which is an important sign of the dam instability. In other words, if the plastic zones in the dam foundation have not been run-through, then even if the displacement at the key point of the dam changes abruptly, the system is not necessarily considered to be unstable. Therefore, when applying the above method, the dynamic changes of the plastic zones with different reduction coefficients need to be considered comprehensively to analyze the dynamic stability of the dam.When the reduction coefficient varies from 1.0 to 2.0, the soft interlayer, which tends to the upstream of the dam, first enters the plastic yielding stage, while the plastic zones gradually and continuously develop at the dam heel, dam toe, and other weak parts. When the reduction coefficient equals 2.4, most of the rock mass in the dam foundation is in the plastic yielding state, with the plastic zones not being entirely run-through. And when the reduction coefficient is 2.6, the rock mass in the dam foundation completely enters the plastic yielding state, and the plastic zones are entirely run-through (shown in Figure21). As described in Section 3.2.1, the displacements at the dam heel and the dam toe are confronted with mutations, with the value of δ less than zero. And according to the results combining the analyses of cumulative displacement mutations and of the dynamic distribution of the plastic zone, it can be concluded that the dam-foundation system is in an instability state. The sliding area searched by the adaptive quantitative method is the plastic zone in the dam foundation in Figure 21(d). This method can solve the problem that the sliding area is difficult to determine by experience in the case of complex geological forms of the dam foundation.Figure 21
Distribution diagram of the plastic zone in the foundation with different reduction coefficients. (a)k = 1.0, (b) k = 2.0, (c) k = 2.4, and (d) k = 2.6.
(a)
(b)
(c)
(d)In addition, because the duration of the earthquake is much shorter than the change time of the seepage field and the earthquake action is reciprocating, the influence of the change of the seepage field is generally not considered during the earthquake action. However, after the earthquake, if only a small part of the dam heel yields and the yield area does not extend to the curtain position, the seepage field will have little effect on the stability of the dam. However, if the yield area at the dam foundation extends to the location of the curtain, possibly leading to the damage of the curtain, the osmotic pressure at the dam foundation will gradually increase, which is detrimental to the stability of the dam and requires danger-reinforcement treatment.


## 3.2.1. Analysis for the Dynamic Stability of the Dam Foundation

Since the seismic load changes with time, the displacement at specific points suddenly becoming larger cannot be considered as a basis for judging whether the dam is stable or not, but the final cumulative displacement after the earthquake can be utilized as a criterion for instability discrimination. For this example, two key points, respectively, located at the dam toe and the dam heel, are selected and the displacement-time curves at the two points are drawn (shown in Figures9 and 10).Figure 9
Displacement-time curves of the point at the dam toe.Figure 10
Displacement-time curves of the point at the dam heel.According to the dynamic calculation results with different reduction coefficientsk, the value of the cumulative displacement, y, at the two points under seismic action is obtained. When the reduction coefficient is ki, take the first i values of k and the displacement y to do the quartic polynomial fitting. And then substitute ki into equations (5)–(10) and calculate the mutation factor Δ with different reduction coefficients.The cumulative displacement curves of the point at the dam toe with the specific reduction coefficients are shown in Figure11. When the reduction coefficient equals 2.2, a quartic polynomial fit is performed based on the results with ks that are smaller than 2.2 (Figures 11(a)), and the corresponding mutation factor Δ = 0.082235031 is calculated. Simiarly, the calculation of the mutation factors proceeds when k, respectively, equals to 2.4, 2.6, 2.8, 3.0, and 3.2. It can be seen in Figures 11(a)–11(e) that the goodness of fit is very close to 1, indicating that the accuracy of the quartic polynomial is sufficient and that the cusp catastrophe model can effectively reflect the stability state of the system. In addition, in order to determine the value of the safety factor at Δ = 0, a Bessel function is used to fit the relationship between the dynamic strength reduction factor k and the mutation factor Δ (see Figure 12). It can be seen from the curve in Figure 13 that as the dynamic strength reduction coefficient increases, the value of Δ at the dam toe gradually decreases. And when the value of Δ is less than 0, the system is considered to be unstable; when the curve is nearly close to the horizontal axis, it indicates that the system is on the verge of the critical state. Last but not the least, by writing a macro command, Δ = 0 is substituted into the Bessel relations of cumulative displacement of the point at the dam toe, and the corresponding k is calculated to be 2.58089.Figure 11
(a) The fitting curve of reduction coefficient against displacement of the point at the dam toe (k = 2.2). (b) The fitting curve of reduction coefficient against displacement of the point at the dam toe (k = 2.4). (c) The fitting curve of reduction coefficient against displacement of the point at the dam toe (k = 2.6). (d) The fitting curve of reduction coefficient displacement of the point at the dam toe (k = 2.8). (e) The fitting curve of reduction coefficient against displacement of the point at the dam toe (k = 3.2).
(a)
(b)
(c)
(d)
(e)Figure 12
Relationship between strength reduction.Figure 13
(a) The fitting curve of reduction coefficient against displacement of the point at the dam heel (k = 2.2). (b) The fitting curve of reduction coefficient against displacement of the point at the dam heel (k = 2.4). (c) The fitting curve of reduction coefficient against displacement of the point at the dam heel (k = 2.6). (d) The fitting curve of reduction coefficient against displacement of the point at the dam heel (k = 2.8). (e) The fitting curve of reduction coefficient against displacement of the point at dam heel (k = 3.2).
(a)
(b)
(c)
(d)
(e)Similarly, the curve of cumulative displacement against reduction coefficient of the point at the dam heel is shown in Figure13. It can be seen that the goodness of fit of the quartic polynomial is rather close to 1 as well, indicating that the cusp catastrophe model can well reflect the stability of the dam. And with the bessel function applied, δ = 0 is substituted into the bessel relations of cumulative displacement of the point at the dam heel, and the corresponding k is calculated to be 2.57759 (Figure 14), which is slightly less than that at the dam toe. Besides, it can be drawn from Figures 11–14 that due to the reciprocating effects of the earthquake, the dynamic instability mechanism of the dam is completely different from the static one.Figure 14
Relationship between strength reduction factor and mutation factor.The results of the catastrophe determination of the points, respectively, at the dam toe and dam heel with different reduction coefficients are shown in Table3.Table 3
Catastrophe discrimination results of the points at the dam toe and heel with different reduction coefficients.
Reduction coefficient
Key points
u
v
Δ
Result of catastrophe discrimination
2.2
Dam heel
−0.183078066
−0.0333325
0.005453211
Noncatastrophe
Dam toe
−0.490153037
−0.1431487
0.082235031
Noncatastrophe
2.4
Dam heel
−0.082455007
−0.0122254
0.001793047
Noncatastrophe
Dam toe
−0.339602824
−0.0862048
0.043978566
Noncatastrophe
2.6
Dam heel
−0.093978032
0.00895964
−0.001152577
Catastrophe
Dam toe
−0.160788464
−0.0219337
−0.003638031
Catastrophe
2.8
Dam heel
−0.209970046
−0.0229496
−0.022807694
Catastrophe
Dam toe
−0.214829545
−0.0217247
−0.026916081
Catastrophe
3.0
Dam heel
−0.118377439
−0.0031959
−0.006359630
Catastrophe
Dam toe
−0.146569898
−0.0122332
−0.008554322
Catastrophe
3.2
Dam heel
−0.185233084
−0.0167385
−0.017857578
Catastrophe
Dam toe
−0.19446602
−0.0162687
−0.022270416
Catastrophe


## 3.2.2. Sensitivity Analysis to the Strength Parameters of the Dam Foundation

Due to the constraints of engineering survey technology, it is difficult to accurately obtain the material parameters of the dam foundation. Therefore, it is necessary to perform a sensitivity analysis to the strength parameters of the dam foundation.First, the strength parameters of the rock mass at the dam foundation are increased by 20%, and the dynamic stability analysis of the dam-foundation system is conducted to obtain the values of the mutation factor Δ corresponding to each different strength reduction coefficient (Table4). And then the Bessel function is used to determine the reduction coefficient when Δ equals 0, and the relationship curves between the mutation factor Δ and the strength reduction factor k are plotted (Figures 15 and 16). It can be seen from the figures that after increasing the strength parameters at the dam foundation by 20%, the dynamic stability safety factor at the dam toe is equal to 2.58368, which is slightly larger than that with the parameters not increased (2.58089). And the dynamic stability safety factor at the dam heel is 2.58790, which is slightly greater than that with the parameters not increased (2.57759).Table 4
Catastrophe discrimination results of the points at the dam toe and heel with different reduction coefficients (20% increase in the dam foundation strength).
Reduction coefficient
Key points
u
v
Δ
Result of catastrophe discrimination
2.2
Dam heel
0.036515
0.005371
0.000974
Noncatastrophe
Dam toe
0.072672
0.009195
0.003818
Noncatastrophe
2.4
Dam heel
−0.07127
−0.01014
0.001325
Noncatastrophe
Dam toe
0.05249
0.012883
0.005060
Noncatastrophe
2.6
Dam heel
−0.05448
0.001862
−0.000553
Catastrophe
Dam toe
−0.06151
0.000709
−0.000917
Catastrophe
2.8
Dam heel
−0.13847
−0.00122
−0.010581
Catastrophe
Dam toe
−0.16737
0.001011
−0.018727
Catastrophe
3.0
Dam heel
−0.11747
−0.00111
−0.006450
Catastrophe
Dam toe
−0.11705
−0.00115
−0.006379
Catastrophe
3.2
Dam heel
−0.18426
0.00634
−0.023939
Catastrophe
Dam toe
−0.15204
0.00152
−0.013995
CatastropheFigure 15
Relationship between strength reduction factor and mutation factor at the dam toe (20% increase in the dam foundation strength).Figure 16
Relationship between strength reduction factor and mutation factor at the dam heel (20% increase in the dam foundation strength).Through the above analysis, it can be drawn that when the strength parameters at the dam foundation are increased by 20%, the dynamic stability safety factors become slightly larger, but the change is not obvious.Similar to the above procedures, the strength parameters at the dam foundation are decreased by 20%, and the values of the mutation factor Δ corresponding to each different strength reduction coefficient are obtained (Table5). And then the Bessel function is used to determine the reduction coefficient when Δ equals 0, and the relationship curves between the mutation factor Δ and the strength reduction factor k are plotted (Figures 17 and 18). It can be seen from the figures that after decreasing the strength parameters at the dam foundation by 20%, the dynamic stability safety factor at the dam toe is equal to 2.56634, which is slightly smaller than that with the parameters not increased (2.58089). And the dynamic stability safety factor at the dam heel is 2.55017, which is slightly smaller than that with the parameters not increased (2.57759).Table 5
Catastrophe discrimination results of the points at the dam toe and heel with different reduction coefficients (20% decrease in the dam foundation strength).
Reduction coefficient
Key points
u
v
Δ
Result of catastrophe discrimination
2.2
Dam heel
0.023094
0.005313
0.000811
Noncatastrophe
Dam toe
0.047434
0.011929
0.004269
Noncatastrophe
2.4
Dam heel
0.05625
0.008574
0.002697
Noncatastrophe
Dam toe
0.069669
0.020231
0.012404
Noncatastrophe
2.6
Dam heel
−0.10016
0.009008
−0.001829
Catastrophe
Dam toe
−0.08765
−0.0005
−0.002687
Catastrophe
2.8
Dam heel
−0.15368
−0.00321
−0.014241
Catastrophe
Dam toe
−−0.18609
−0.00098
−0.025751
Catastrophe
3.0
Dam heel
−0.13073
0.00042
−0.008933
Catastrophe
Dam toe
−0.12937
−−0.00341
−0.008347
Catastrophe
3.2
Dam heel
−0.15189
−0.00204
−0.013905
Catastrophe
Dam toe
−0.15459
0.002175
−0.014649
CatastropheFigure 17
Relationship between strength reduction factor and mutation factor at the dam toe (20% decrease in the dam foundation strength).Figure 18
Relationship between strength reduction factor and mutation factor at the dam heel (20% decrease in the dam foundation strength).Through the above analysis, it can be drawn that when the strength parameters at the dam foundation are decreased by 20%, the dynamic stability safety factors become slightly smaller, but the change is not obvious.


## 3.2.3. Determination of the Sliding Area in the Dam Foundation

As for the influence of the weak layer orientation on the stability analysis in this paper, the relevant engineering practice and experimental research results show that there may be two types of instability failure modes for gravity dams on rock foundations: one type is that there are inclined weak layers at the rock mass in the dam foundation. And with the weak layers being subject to external loads, sliding blocks may possibly form at the dam foundation, which will result in deep sliding (Figure19). The other is that there are no obvious weak layers in the dam foundation, making it difficult to assume the possible sliding direction of the dam foundation. And due to the reciprocating action of earthquakes, the cracking and crushing occur at the rock mass in the dam foundation, and the yielding area gradually expands until the final penetration. Finally, the dam body and some part of the foundation suffer instability and failure as a whole, as shown in Figure 20.Figure 19
Failure mode with obvious sliding surface.Figure 20
Failure mode without obvious sliding surface.As for the project in this paper, the rock mass in the dam foundation presents a layered distribution approximately perpendicular to the horizontal plane, with similar material parameters. Besides, there is no obvious weak layer. Therefore, the instability mode of the project belongs to the second type. According to the results through finite element calculation, the rock mass near the upstream dam heel and the downstream dam toe yielded first. And with the increase of the strength reduction factor, the yield area gradually expanded until it was completely penetrated, eventually causing the instability failure of the dam body and some part of the dam foundation as a whole.Due to the complex distribution features of the rock mass in the dam foundation, it is difficult to delineate the possible sliding area by experience. Therefore, the capability of searching for the sliding area for the dynamic stability analysis of the dam is of great concern. With the increase of the strength reduction coefficient, the rock mass in the dam foundation gradually yields, and the plastic zones develop with time until final penetration, which is an important sign of the dam instability. In other words, if the plastic zones in the dam foundation have not been run-through, then even if the displacement at the key point of the dam changes abruptly, the system is not necessarily considered to be unstable. Therefore, when applying the above method, the dynamic changes of the plastic zones with different reduction coefficients need to be considered comprehensively to analyze the dynamic stability of the dam.When the reduction coefficient varies from 1.0 to 2.0, the soft interlayer, which tends to the upstream of the dam, first enters the plastic yielding stage, while the plastic zones gradually and continuously develop at the dam heel, dam toe, and other weak parts. When the reduction coefficient equals 2.4, most of the rock mass in the dam foundation is in the plastic yielding state, with the plastic zones not being entirely run-through. And when the reduction coefficient is 2.6, the rock mass in the dam foundation completely enters the plastic yielding state, and the plastic zones are entirely run-through (shown in Figure21). As described in Section 3.2.1, the displacements at the dam heel and the dam toe are confronted with mutations, with the value of δ less than zero. And according to the results combining the analyses of cumulative displacement mutations and of the dynamic distribution of the plastic zone, it can be concluded that the dam-foundation system is in an instability state. The sliding area searched by the adaptive quantitative method is the plastic zone in the dam foundation in Figure 21(d). This method can solve the problem that the sliding area is difficult to determine by experience in the case of complex geological forms of the dam foundation.Figure 21
Distribution diagram of the plastic zone in the foundation with different reduction coefficients. (a)k = 1.0, (b) k = 2.0, (c) k = 2.4, and (d) k = 2.6.
(a)
(b)
(c)
(d)In addition, because the duration of the earthquake is much shorter than the change time of the seepage field and the earthquake action is reciprocating, the influence of the change of the seepage field is generally not considered during the earthquake action. However, after the earthquake, if only a small part of the dam heel yields and the yield area does not extend to the curtain position, the seepage field will have little effect on the stability of the dam. However, if the yield area at the dam foundation extends to the location of the curtain, possibly leading to the damage of the curtain, the osmotic pressure at the dam foundation will gradually increase, which is detrimental to the stability of the dam and requires danger-reinforcement treatment.


## 4. Conclusion

(1)
Considering the lack of research on the dynamic instability of dams, a dynamic stability analysis method that can analyze different failure modes is proposed in this paper. This method breaks through the framework of the traditional rigid body limit equilibrium method and deals with the uncertainty of instability discrimination by the inflection points of displacement.
(2)
Through the method proposed in this paper, the mutation of the cumulative displacements at the critical points are taken as the indication of the dam instability, with the catastrophe theory and Bessel function adopted to quantitatively define the displacement catastrophe points, and then the safety factors of the dynamic stability of the dam are obtained. It makes up for the lack of quantitative analysis of the traditional displacement catastrophe criterion to some degree.
(3)
The examples show that the dynamic instability mode of the dam is completely different from the static one. And the research in this paper reveals the inherent law of the dynamic instability of dams and helps provide a new idea and method for the dynamic stability analysis of dams. The method proposed in this paper can also be applied to the dynamic stability analysis of slopes, earth dams, and other projects.



---

*Source: 3702712-2019-12-19.xml*