Abstract
The paper focuses on the failure process and mechanism of the concrete gravity dam considering different nonlinear models under strong earthquakes. By taking a typical monolith of a concrete gravity dam as a case study, a comparative analysis of the failure process and mechanism of the dam considering the plastic damage model and the dynamic contact model, respectively, is performed using the seismic overload method. Moreover, the ultimate seismic capacity of the dam is evaluated for both of the nonlinear models. It is found that the ultimate seismic capacity of the dam is slightly different, but the failure process has significant distinctions in each model. And, the damage model is recommended when the conditions permit.
1. Introduction
Seismic safety of high concrete dams is extremely important. Once a reservoir dam holding hundreds of millions of tons of lake water fails, it will cause unimaginable catastrophic consequences [1]. To ensure the seismic safety of dams, it is necessary to study the seismic response analysis and failure mechanism of high dams in the strong earthquake area further.
The damage, cracking, or even collapse of dams could happen leading to the stress redistribution when the stress exceeds the allowable strength under strong earthquakes. Therefore, the assumption of linear behavior is no longer applicable, and the nonlinear situation of the material should be considered to conform to the actual situation better [2]. The Standard for Seismic Design of Hydraulic Structures [3] also stipulates that, for gravity dams of a seismic fortification of class A with complex structures or complex geological conditions, the nonlinear influence of materials should be considered in the finite element analysis. In view that the concrete material is a kind of quasi-brittle material, the tensile damage occurs before the compressive damage generally. Therefore, the tensile damage of dam concrete is the major concern ignoring the impact of compression damage in the nonlinear seismic analysis [4]. However, for high concrete dams in the strong earthquake area, the upstream water pressure is already huge. Once subjected to strong earthquakes, the static and dynamic comprehensive compressive stress of the dam at the dam toe often exceed the compressive strength of the dam concrete causing compression damage. Thus, it is necessary to study and demonstrate the influence of concrete compression damage on the damage process, mechanism, and ultimate seismic capacity of the dams.
Besides, the tensile strength and shear strength of the RCC gravity dam with the rolled surface are usually lower than that of normal concrete, which may lead to the decrease of stability and safety of the dam [5]. The Code for Seismic Design of Hydraulic Structures [6] also stipulates that the stability of the rolled layer should be checked for RCC gravity dams. Moreover, according to the existing research studies, such as Koyna gravity dam [7], Xinfengjiang dam [8], and Sefid Rud dam [9], the static and dynamic synthesis maximum principal stress of the dam head are relatively large [10]. It can be concluded that the dam head is a weak part of gravity dams, and it is easy to form macrocracks in the upstream and downstream under the strong earthquakes [11]. Although no overall instability occurs under strong earthquakes, penetrating cracks in the upstream and downstream at the head of the dam appear [12]. Therefore, the influence of weak layers and parts of the dams should be considered in the seismic analysis.
At present, a more commonly used model for simulating material nonlinearity is the plastic damage model proposed by Lubliner et al. [13] and improved and developed by Lee and Fenves [14, 15]. The contact model can be used for gravity dams with deep sliding problems and weak parts with a local cracking inside the dam and the dam-foundation interface that may crack and slip during strong earthquakes [16]. Based on the existing research studies, the concrete plastic damage model [17, 18] and the dynamic contact model [19, 20] are used for seismic analysis of a typical nonoverflow monolith of a gravity dam in this paper, respectively. The impact of tensile and compressive damage on the seismic damage process and failure mechanism of the dam under strong earthquakes is analysed and discussed in detail. Moreover, the ultimate seismic capacity of the concrete dam for two models is compared and analyzed, which provides a reference for the engineering application.
2. Plastic Damage Model and Dynamic Contact Model
2.1. Plastic Damage Model
The concrete plastic damage model [21, 22] and both the compression and tensile damage are considered in this paper. The model uses the yield criterion in plastic mechanics to judge whether the concrete is in the plastic state. If it enters the damaged state, the flow rule is used to evaluate the plastic strain (residual deformation), and the stiffness degradation of the material is calculated according to the damage evolution curve.
2.1.1. Yield Function
When the stress is inside the yield surface, the concrete material is in an elastic state. When the stress is on the yield surface, the concrete material begins to enter the plastic state. The yield function can be expressed as follows:withwhere is the effective stress tensor, and are dimensionless material constants, is the effective stress deviator, defined as , is the second-order tensor, is the effective hydrostatic pressure, is the Mises equivalent effective stress, is the maximum principal effective stress, is the effective compressive cohesion stress, is the effective tensile cohesion stress, and and are the tensile and compressive equivalent plastic strain; the typical experimental values of the ratio for concrete are in the range from 1.10 to 1.16; yielding values of are between 0.08 and 0.12, where and .
2.1.2. Flow Rule
The plastic damage model assumes the nonassociated potential flow:where is the plastic strain and is the plastic flow factor.
The flow potential is chosen for the Drucker–Prager hyperbolic function:where is the dilation angle measured in the plane at high confining pressure, is a parameter, referred to as the eccentricity, and is the uniaxial tensile stress at failure. This paper refers to the value of parameters in the concrete gravity dam calculation example given in Reference [23], taking and .
2.1.3. Damage and Stiffness Degradation
When the concrete is subjected to the uniaxial tension and uniaxial compression, the plastic strain can be expressed aswhere and are the tensile damage factor and compression damage factor, is the initial (undamaged) modulus of the material, and are the tensile stress and compressive stress corresponding to the total strain, is the cracking strain, and is the inelastic strain.
The stress-strain relations under uniaxial tension and compression loading are, respectively,
Under the uniaxial cyclic loading conditions, the stiffness degradation mechanisms of the material can be expressed aswhere is the elastic modulus of the material, is the stiffness degradation variable, and are functions of the stress state that are introduced to model stiffness recovery effects associated with stress reversals, and the weight factors and , which are assumed to be material properties, control the recovery of the tensile and compressive stiffness upon load reversal; in this paper, and .
Under the action of multiaxial cyclic loading, the elastic stiffness degradation is assumed to be isotropic, which can be expressed by a single scalar variable . The stress-strain relation of the viscoplastic model is given aswhere is the initial (undamaged) elastic stiffness of the material, is the total strain, and is the plastic part of the strain.
The equivalent plastic strain rates are evaluated according to the expressions:where and are, respectively, the maximum and minimum eigenvalues of the plastic strain rate tensor and is a stress weight factor that is equal to one if all principal stress and are positive and equal to zero if they are negative. The Macauley bracket is defined by . Under uniaxial conditions, since in tension, in compression.
2.2. Dynamic Contact Model
In the dynamic contact model, the contact problem established by adding a supplementary equation constructed by the contact constraints to the dynamic equation is solved in the form of point-to-point contact. The explicit integration method of the central differential method combined with the Newmark constant average acceleration method is adopted, and the detailed procedures are shown in Reference [24, 25]. The displacement expression of th DOF of node can be expressed as follows:where and are the components of and in the j direction at time t, respectively, and are the stiffness matrix and damping matrix of the structure, , , and are the external load, displacement, and velocity column vectors, the subscript represents the component of the vector in the -DOF direction of node , is the concentrated mass on the node , and is the time step.
For the convenience of derivation, formula (13) is simplified into three parts:
To get the normal contact force, taking the contact point pairs and for instance, the normal contact condition including the initial tensile strength should meet the geometric condition of normal mutual inviolability:where is the normal vector of the contact surface at point ; let by formulas (16) and (19), and can be obtained:
When the normal contact between and occurs, the tangential friction force should be calculated; let , and the tangential contact force can be obtained as follows:
When the static friction force exceeds ( is the static friction coefficient), it means that and are in the dynamic friction state:
From the above calculation, , , and can be obtained, and finally, can be obtained.
3. Model Establishment
3.1. Finite Element Model
Taking a typical nonoverflow monolith of a concrete dam as a case study, the dam bottom and crest elevation are 1970 m and 2155 m, respectively. The elevation of the upstream breakpoint is 2040 m, the width of the dam crest is 16 m, and the width of the dam bottom is 165.5 m. The geometric dam section and concrete partition are shown in Figure 1. The finite element model of the dam is shown in Figure 2. The dam foundation system is discretized into 14592 nodes and 14320 elements. The element size for the dam body is about 2 m. Material properties of the dam concrete and foundation rock are described in Table 1 [26].
3.2. Static and Dynamic Loads
The static loads mainly include the self-weight, upstream and downstream hydrostatic pressure, and silting pressure. Westergaard’s added mass without considering the compressibility of the reservoir water is used to consider the hydrodynamic interaction between the dam and the reservoir. The normal upstream water level is 2150 m, and the corresponding downstream water level is 2019.25 m. The elevation of silt is 2023.7 m, and its floating bulk density is 8 with the friction angle of 12°. The horizontal and vertical peak ground accelerations of the design earthquake are both 0.4005 g, and the time histories of normalized acceleration are shown in Figure 3.
(a)
(b)
4. Seismic Damage Analysis considering Tension and Compression Damage
4.1. Plastic Damage Model
Figures 4–7 show the damage evolution curve used in this paper, which is obtained by the concrete material tests. According to the CDP model parameters in Reference [27], the compression damage evolution is converted and adjusted appropriately by formula (6).
(a)
(b)
(a)
(b)
(a)
(b)
(a)
(b)
4.2. Discussion of the Results
The damage evolution process of the concrete gravity dam is discussed in detail through the earthquake overload analysis [28]. The influence of the compression damage on the ultimate seismic capacity and damage mechanism of the gravity dam is mainly studied.
4.2.1. Damage Rupture Process Only considering the Tension Damage Model
Only the tension damage of the dam concrete is considered, and the earthquake overload analysis is performed with increasing overload coefficients of 1.2, 1.4, 1.52, 1.53, and 2.0. The results are shown in Figure 8. The damage occurs at the dam heel first, and the macrocrack depth with a damage coefficient greater than 0.8 is about 18.5 m under the design earthquake. In view that the bedrock material is considered as linear elasticity, the damage at the dam heel is caused by stress concentration and the strong restraint of the foundation. Therefore, the dam heel damage is not discussed in this study. The damage occurs at the downstream face of the dam head for the overload coefficient of 1.2 and gradually expands to the upstream face with the increasing of the overload coefficient. Penetrating cracks (i.e., when a crack extends completely between the upstream and the downstream faces) appear on the upstream and downstream surfaces for the overload coefficient of 1.53. If the penetrating crack of the dam body is taken as the failure criteria, it is suggested that the ultimate seismic capacity of the concrete dam is 1.52 times the design earthquake.
(a)
(b)
(c)
(d)
(e)
(f)
Figure 9 shows the number of damaged elements of the dam head under different overload coefficients, and Figures 10–14 show the time histories of tensile damage variables and relative displacements of point pairs i.e., upper and lower relative nodes of an element under different earthquake overload coefficients. According to these results, the damage begins to appear at the time of 6.47 s for the earthquake overload coefficient of 1.2. Meanwhile, relative displacements of the node pairs change rapidly, and the damage gradually increases. For the earthquake overload coefficient of 1.53, the damage of the dam head slope first appears at the time of 6.43 s. Then, the cracking gradually spreads across the upstream surface, and the relative displacements of the node pairs gradually increase simultaneously with the deformation of the new damage element accumulating to the downstream element. The damage of the upstream surface begins to appear at the time of 16.10 s. The maximum tensile strain of the Gauss point of element 3400 is at the time of 16.11 s, which exceeds the tensile strain corresponding to the ultimate stress. Afterward, the damage expands rapidly with the material entering the softening stage. At the time of 16.16 s, the penetrated cracking through the upstream and downstream is formed. Meanwhile, the relative displacement of the upstream node pairs at the dam head changes rapidly. Moreover, the crack has residual displacement after the earthquake. The damage process will accelerate with the increase of the overload coefficients, and the damage appears and has penetrated the dam head at the time of 3.91 s and 8.18 s, respectively, for the earthquake overload coefficients of 2.0. The relative displacement of the node pairs also has some increase with the expansion of the damage domain.
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
4.2.2. Damage Rupture Process considering Both the Tension and Compression Damage Model
Figure 15 shows the tension damage development process under different overload coefficients considering the tension and compression damage simultaneously. Figure 16 shows the compression damage development process under different overload coefficients.
(a)
(b)
(c)
(d)
(e)
(f)
(a)
(b)
(c)
(d)
(e)
(f)
It can be obtained that the compression damage near the dam toe expands slightly with the increasing of the earthquake overload coefficients (Figure 16). The compression damage mainly occurs at the dam heel for small overload coefficients. As the overload coefficient increases, the compression damage occurs at the upstream and downstream slopes and gradually expands and the damage at the dam toe begins to appear. Therefore, tension damage and compression damage occur simultaneously at the dam heel.
To explain this phenomenon, the results of the overload coefficient of 1.0 are taken as an example of a detailed discussion. The concrete strength grade at the dam heel is , and the initial compressive yield stress is 14.8 MPa in the uniaxial state (Figure 6). Figure 17 describes that the compression damage first appears at 5.82 s for the overload coefficient of 1.0. Table 2 describes the relevant variables at the integral point of dam heel for the overload coefficient of 1.0 at 5.82 s. The maximum compressive stress near the dam heel does not exceed the initial compressive yield stress (Table 1). Meanwhile, the maximum principal stress at the Gauss point of the element at the dam heel is 1.925 MPa, and the minimum principal stress is -0.839 MPa. Hence, the yield surface function F = 3.276 > 0 can be obtained from formula (1),which exceeds the yield surface. According to the flow rule, the plastic strain tensor can be calculated as . Since it is in a biaxial tension and compression state, the stress weight factor is 0.696 by the formula, and the equivalent tension and compression plastic strains are and , respectively, by formulas (10) and (11). The tension and compression damage variables are obtained as 0.63228 and 0.00246, respectively, and the compression damage variable is not 0. Thus, with the increase of the overload coefficient, the tension damage and compression damage appear simultaneously at the upstream and downstream slope.
Figure 18 shows the compressive stress-strain of concrete. Table 3 shows the related variable values of the element 1327 near the dam toe under different overload coefficients. It can be seen that the compression damage has occurred for the overload coefficient of 1.0, but the maximum compressive strain and compressive equivalent plastic strain of the Gauss node are less than the values of the ultimate compressive state. Therefore, the dam concrete is in the hardening stage, and the loading capacity is not decreased resulting in a small damage propagation.
From Figure 19, the damaged area is still slightly even for the earthquake overload coefficient of 3.00. Though the compression damage may appear at the compressive stress concentration region under strong earthquakes, the compression constitutive of the concrete reflects the characteristics of ductile materials with strengthening properties, and the ultimate compressive strain is relatively large. As a result, serious tension damage of the dam has already developed when the compressive strain reaches the softening stage. So, the tension damage of the dam is dominant under strong earthquakes.
Moreover, the tension damage begins to appear at the downstream slope of the dam head for the overload coefficient of 1.2. Penetrating macrocracks through the upstream and downstream appear when the overload coefficient is 1.55. Therefore, it is suggested that 1.54 times the design earthquake can be considered as the ultimate seismic capacity, which is similar to that only considering the tensile damage.
5. Seismic Analysis considering the Dynamic Contact Model
As shown in Figure 20, the joints are preset according to the damage penetrating position of the dam in the damage model. The number of joint pairs from upstream to downstream is 1–16, and its tensile strength is 1.6 MPa, which is consistent with the damage model. Therefore, the dynamic contact model with the tensile and shear strength of joints is adopted for the seismic overload analysis of the concrete gravity dam. The overload coefficients are 1.00, 1.44, 1.47, 1.48, 1.49, and 1.50 times of the design earthquake, respectively.
Results show that the preset contact joints did not open until the overload coefficient reaches 1.49, indicating that no damage and cracking occurs. From Figures 21 and 22, the tensile stress of the downstream preset joints exceeds the tensile strength, and cracks begin to appear for the overload coefficient of 1.49 at the time of 8.82 s. Meanwhile, the displacements of the joint pairs occur rapidly, and the cracking gradually extends to the upstream surface. The upstream preset contact joints start to open at the time of 13.95 s, indicating the appearance of the penetrated cracking through the upstream and downstream. If the penetrating crack of the dam is taken as the failure criteria, the ultimate seismic capacity of the concrete gravity dam is suggested as 1.48 times the design earthquake considering the dynamic contact model.
(a)
(b)
(a)
(b)
6. Comparative Analysis of Two Models
6.1. Comparative Analysis of the Failure Mechanism
According to the results of the plastic damage model, the damage of the downstream surface of the dam head first appears when the overload coefficient is 1.2. With the increase of the earthquake overload coefficient, the damage gradually extends to the upstream surface. The penetrating cracks of the dam appear when the overload coefficients are 1.53 and 1.55 for the damage model only considering the tensile damage and that considering the tension and compression damage, respectively. However, for the contact model, when the overload coefficient is 1.49, the preset contact joint causes tensile cracks and completely breaks after the earthquake.
To sum up, the earthquake overload coefficient of the cracking at the dam head in the damage model is lower than that in the contact model, but the failure process of the contact model is faster than that of the damage model. The main reasons lie in that the damage cracking of concrete in the damage model is mainly judged by the state of the maximum principal stress at the element Gauss point, and there is no restriction on the damage propagation direction. However, for the contact model, the contact status changes among the separation, adhesion, and contact states. The failure of the contact surface is mainly judged by the normal and tangential force, and the influence area of the node is also larger, so the earthquake overload coefficient of the cracking at the dam head of the damage model is lower than that of the contact model. Besides, in the damage model, the dam concrete will enter the softening stage when the maximum principal stress exceeds the tensile strength. In the softening stage, the damage develops and the stiffness gradually degenerates, while the bearing capacity will not be completely lost immediately. Moreover, the elements near the penetrating part of the dam head are also damaged, which plays a role in energy dissipation. However, for the contact model, once the cracks appear in the normal direction, the contact joints cannot bear the normal tensile stress, and the tensile bearing capacity will be completely lost, so the damage and cracking process are faster than those of the damage model.
6.2. Comparative Analysis of the Ultimate Seismic Capacity
If the penetrating crack of the dam is taken as the failure criteria, the ultimate seismic capacities of the damage model only considering tensile damage and considering tensile and compression damage are 1.52 and 1.54 times the design earthquake, respectively. For the contact model, it is suggested as 1.48 times the design earthquake. Therefore, the difference in the ultimate seismic capacities of the concrete gravity dam obtained by the two different models is slight.
Compared with the contact model, the damage model is more consistent with the actual situation without presetting the contact joints. It only needs to consider the damage softening process of the material. However, the damage model requires a refined finite element mesh leading to a large number of computing efforts. The contact model is more suitable for RCC gravity dams with clear weak interfaces. The coarse finite element mesh is sufficient for the analysis, and computing efforts are reduced. Moreover, it can also give a reasonable result to evaluate the ultimate seismic capacity of concrete dams, while, for the dams without clear weak interfaces, the damage model is recommended when the calculation condition permits.
7. Conclusion
In this paper, seismic overload analysis is performed by taking a typical nonoverflow monolith of a concrete gravity dam considering the damage model and contact model. The tension and compression damage development process and failure mechanism of the concrete gravity dam under the strong earthquake including its ultimate seismic capacity are of comparative analysis. The main conclusions are as follows:
For the concrete dam with the material damage model, with the increase of the overload coefficients, the tension and compression damage may occur simultaneously in multiaxial stress conditions, but the damage of the concrete gravity dam is mainly controlled by the tensile damage. A slight difference in the development pattern of the damage is found between the damage model considering only the tension damage and that considering the tension and compression damage.
Compared with the contact model, the earthquake overload coefficient of the cracking of the damage model is earlier, but the damage development process is slower. Based on the failure criteria of the penetrating crack of the dam, the ultimate seismic capacity of the gravity dam with the damage model considering only tensile damage is 1.52 times the design earthquake, and that, with the damage model considering the tension and compression damage, it is 1.54 times the design earthquake. However, the ultimate seismic capacity based on the contact model is 1.48 times the design earthquake. In summary, the ultimate seismic capacities of the gravity dam under the two models are similar, and the contact model is slightly lower.
It is suggested that the contact model should be used for the RCC gravity dam with clear weak interfaces and coarse meshes, while for the dams without clear weak interfaces, the damage model is recommended when the calculation conditions permit.
Data Availability
All data generated or analyzed during this study are included within this article.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
This study was supported by the National Key Research and Development Program (Grant no. 2016YFC0401807) and National Key Research and Development Program (Grant no. 2017YFC0404903) and Research Project of China Three Gorges Corporation (Grant no. XLD/2115). The authors are grateful for these supports.