Abstract
The concrete face of a rockfill dam is a long and thin slab structure, which is highly susceptible to fracture when subjected to the settlement of the dam. The study of the generation and propagation of cracks in the concrete face of rockfill dams is of great significance to dam construction and face crack prevention. In this study, the initiation and propagation of cracks in the concrete face of a rockfill dam are investigated using an extended finite element method (XFEM) and ABAQUS software for the Gongboxia concrete-face rockfill dam. A numerical model for this dam is established using a finite element method, and the face stress and deformation distributions are obtained. Based on the results, a numerical model is built to find the location where cracks are initiated in the face. The displacement of the entire model is treated as the equivalent displacement for the numerical model. XFEM is utilized throughout the modeling process to obtain the stress concentration, crack initiation, and crack propagation in the concrete face, and an analysis of crack initiation and propagation is conducted. Finally, the effects of the thickness of concrete covers and reinforcement layers on the stress intensity of crack tips are also discussed using the established numerical model, and techniques for controlling the fracturing of the concrete face have been proposed in this paper.
1. Introduction
Monitored data from projects show that several concrete-face rockfill dams display varying degrees of cracks during construction and operation, most of which are transverse cracks in the middle or lower parts of the face [1]. Shrinkage stress, temperature stress, and structural stress caused by dam deformation are the main reasons for face cracks [2–4], among which structural stress is the most common. A deflection of the face, which is a part of the antiseepage body of a concrete-face rockfill dam, due to a large deformation can lead to face cracking. This leads to severe leakage when the cracks penetrate and connect with one another. Fine particulates flow along with water during seepage, posing a great risk to the stability and safety of the dam. Thus, the study of crack propagation in the concrete face of rockfill dams is crucial for the prevention of face cracking.
In general, cracks on the concrete face of rockfill dams can be classified into two major groups: (i) microcracks with a random distribution, which impact the compression and tensile strengths of concrete, and (ii) macrocracks, which result in anisotropy of the mechanical properties of concrete. The process of cracking involves the generation of microcracks, followed by their expansion and perforation leading to the formation of macro cracks and subsequent failure. Many researchers studied the process of propagation of cracks and the methods for simulating concrete-face cracking. Jin and Arson [5] used conventional finite element methods to simulate the process of concrete-face cracking based on the crack model and obtained stress distributions during the cracking process. However, for the simulation of discontinuities in crack propagation, separating the crack model requires continuous changes in finite element meshes to model the developing interface discontinuities. Aghajanzadeh and Mirzabozorg [6–8] applied the smeared crack model to simulate the cracking and distribution of cracks in the concrete face of a rockfill dam, but the smeared crack model does not avoid stress locking caused by the lack of a discontinuous displacement model in the approximate displacement field. Jirásek [9] introduced an enhanced assumed strain mode to locate discontinuous interfaces, thereby avoiding the problems arising from utilizing preset discontinuous interfaces during the separation of the crack model and simulating a relatively accurate distribution of cracks on the concrete face. However, the embedded crack model is based on an element level and has the deficiency of not describing the appropriate discontinuous field between elements in a harmonious manner.
The extended finite element method (XFEM) has been widely applied for the simulation of cracks in the concrete face of rockfill dams to overcome the deficiencies of crack propagation modeling using finite element methods. XFEM [10] uses a discontinuous function as the basis of a shape function to describe the discontinuous displacement field in the computational domain. Belytschko et al. used a near-field displacement solution for crack tips to describe crack propagation and achieved relatively good results for crack propagation in concrete. Dolbow et al. [11] introduced discontinuous functions and near-tip crack functions to study crack propagation in the concrete face of rockfill dams and obtained the distribution of the crack tip intensity factor and the elements that affect its change. Singh and Bhardwaj [12] combined XFEM and geometric analysis to predict the fatigue life and crack growth of concrete faces, studied the stress intensity factor of interfacial cracks, and partially explained the mechanism for the initiation of cracks on concrete faces. Fang and Jin [13] proposed a virtual node method for coupling XFEM with the commercial software ABAQUS and successfully simulated continuous-discontinuous processes in the propagation of face cracks. The previous studies are restricted to the detection of one single flaw or crack in the dam as a simple geometry. The level set front propagation method [14], which is implemented in the extended finite element method (XFEM) [15], allows the accurate localization of the flaws. Chao et al. [16] investigated the effects of various dynamical test loads on the crack identification, and fracture is modeled by the XFEM. Zheng et al. [17] used a cohesive crack model (CCM-based XFEM) to simulate and monitor the cracking process in concrete arch dams. Alalade et al. [18, 19] proposed an effective and efficient method, which employs the inverse analysis of a dynamic coupled hydromechanical problem. The numerical model is based on the extended finite element model. The proposed method is able to identify cracks which may be detrimental to structural performance and reliability. Despite the merit of XFEM being able to model and identify crack propagation in concrete, not much has been done so far with regard to the simulation of this damage/cracks for concrete facing rockfill dam storage period using XFEM. The understanding of crack propagation and evolution law of face slab of rockfill dam is still limited.
In this study, the factors and principles governing crack propagation in the concrete faces of rockfill dams are analyzed using XFEM and ABAQUS. The Gongboxia rockfill dam in Qinghai Province, China, is analyzed in this study, for which an overall computational model yields the deformation of the concrete face, which is then treated as an equivalent displacement load acting on the cracks on the concrete face of the numerical model. Based on the cracks on the concrete face in the numerical model, the propagation and evolution of cracks are studied, including stress concentration, crack initiation, and crack propagation. The effects of concrete cover thickness and reinforcement layers on the stress intensity of crack tips are also discussed with the established numerical model.
2. Methods and Principles
2.1. Displacement Model
XFEM not only inherits the advantages of conventional finite element methods in modeling crack propagation but also includes enhancements for addressing fracture problems. This method does not require the meshing and mesh size to be related to the surrounding stress field, and remeshing is not necessary for the simulation of crack propagation. XFEM is based on conventional finite element methods and includes enhanced functions for cracks and crack tips to describe the discontinuity of cracks [20, 21]. Figure 1 shows the propagation of an arbitrary crack within the elements of an XFEM model (elements through which the crack propagates have their nodes annotated by hollow squares, and the elements where the crack tips are located are indicated by hollow circles). Since different functions are introduced in the model [22], the functions for describing an approximation of the displacement fields of standard, penetration, and crack tips are also different, and the displacement in XFEM can be written aswhere are standard nodal degrees of freedom, are standard shape functions, is the number of nodal points for all finite elements, which are represented by the hollow squares and hollow circles in Figure 1, are shape functions of enhanced parts, refers to enhanced shape functions, are enhanced degrees of freedom (DOF), and is the number of nodal points for discontinuous elements.
For the crack penetration unit, the displacement field is given bywhere is the number of nodes in the crack penetration unit, which are shown as hollow squares in Figure 1, and is the Heaviside function given by
2.2. The Description of Cracks in the Level Set Method
Since XFEM allows a discontinuous surface to cross an element, the element mesh is independent of the fracture surface [22]. It is necessary to provide a geometrical description of the discontinuous surface. In XFEM, the level set method is commonly adopted to describe the interface, wherein interfaces are associated with spatial and time domains that are used to describe the discontinuity which develops over time. Another dimension is added to the level set method, which is a one-dimensional interface in R2 [19] expressed aswhere is a level set function.
The level set function is frequently defined by the signed distance function as follows:
The right side of the equation is positive when x lies above the crack and negative elsewhere.
The crack propagation can be described by a function :
When is known, represents the propagation velocity along the unit normal vector of the moving interface, which depicts the propagation and development of the crack.
2.3. Governing Equation
Figure 2 shows the equilibrium of cracks under external loads and boundary conditions, where is the region containing the fracture, crack can propagate in any direction, the force acts on the boundary , and a displacement is generated at the boundary .
The strong form of the equilibrium equation is
The boundary conditions are as follows: Displacement boundary : Force boundary : Inner boundary :
u is the displacement vector, are the known displacements, is stress tensor, n is the normal vector perpendicular to the boundary, is the force of external load on crack body, and the stress is perpendicular to the normal direction in inner boundary, and then the dot product is zero. is the force of crack body.
2.4. Stress and Displacement Fields of the Crack Tip
An open crack is considered to have the most damaging crack propagation pattern and is usually used to simplify a compound crack for simulation. The deformation of the elastic area of a crack tip is an important factor that controls crack propagation, while the mesh density in the model heavily influences the integral point value of the crack tip. When the mesh size is gradually decreased, stress approaches infinity; hence, conventional strength theory, where the strength criterion is built based on stress, is not applicable. In XFEM, a stress intensity factor (K) is introduced to describe the stress field of a crack tip. This factor depicts the singularity and propagation of the crack tip. The theoretical value of K for an open crack in a face with an initial crack under a uniform tensile stress σ is given bywhere K is the stress field intensity factor for the crack, with a unit of ; σ is the tensile stress of the face in all directions; α is a constant related to the face thickness and initial depth of the crack; is the initial crack depth.
A cracked concrete face, subjected to a unidirectional tensile force with L = 5.7 m, B = 2 m, and T = 0.5 m, is chosen for numerical analysis to verify whether the values of K for the crack tip under various tensile stresses σ obtained from XFEM satisfy the theoretical equation describing the relationship between σ and K (Figure 3(a)). A crack with a length of 0.1 m is first assumed to exist at the center of the face. For concrete with elastic modulus E = 28000 MPa and Poisson’s ratio μ = 0.2, a numerical calculation model for concrete cracking is established as shown in Figure 3(b). When , the coefficient α in the theoretical equation for the stress intensity factor [23] becomes
(a)
(b)
Figure 4 shows the values of K versus stress and length-thickness ratios (i.e., the ratios of crack length to plate thickness) for different tensile stresses. It can be seen that K increases with crack length under the same tensile stress, and the trends of both numerical and theoretical values are similar. For the same length-thickness ratio, K of the crack tip increases proportionally with tensile stress and is in accordance with the theoretical value of K for an open crack tip. The above results prove the feasibility of the ABAQUS and XFEM method in obtaining the stress intensity factor of a crack tip. Figure 5 displays different values of K of a crack tip under a tensile stress of 3.0 MPa using XFEM, a closed contour integral, and the theoretical equation. The results show that, compared to the closed contour integral method, the XFEM is closer to the theoretical method, and the difference is less than 10% for various initial crack lengths, thereby illustrating the advantage of XFEM in calculating the value of K for a crack tip.
3. Numerical Analysis of Cracks in the Concrete Face of the Gongboxia Rockfill Dam
3.1. Face Deformation Analysis
3.1.1. Numerical Model of Face Deformation
The stress and deformation of the concrete face should be obtained before investigating the initiation and propagation of cracks in the concrete face of a rockfill dam. Simulating concrete face cracking using the overall model of a concrete-face rockfill dam warrants the consideration of the contact of the face with the cushion, wave wall, toe board, and the hydrostatic water pressure. This usually leads to a large volume of calculation that is very unlikely to converge. To overcome these problems, this study obtains the face displacement after impoundment based on the overall dam model and then considers the face displacement as an external load applied to the cracking of the concrete face model.
The Gongboxia concrete-face rockfill dam (CFRD) has a crest elevation of 2010.00 m, a width of 10 m, a maximum height of 132.2 m, an upstream dam slope of 1 : 1.4, and a downstream dam slope of 1 : 1.55. A cushion and transition layer with a thickness of 0.3 m is placed between the concrete face and rockfill, and the concrete face has a thickness of 0.3 m on the top and 0.68 m on the bottom. The rockfill consists of a main rockfill area and a downstream rockfill area. The cross-sectional view of the dam is shown in Figure 6.
The numerical calculation model is built with ABAQUS and shown in Figure 7. The center of the dam bottom is chosen as the origin of the coordinate frame, the positive X-axis is from upstream to downstream, the positive Z-axis is from the bottom to the top of the dam, and the positive Y-axis is from the right to the left side of the dam. The foundation extends 150 m along the upstream, downstream, and thickness directions, and its width is 10 m. There are a total of 5222 eight-node hexahedral elements and 7764 nodes in the model. The thickness is divided into ten layers to accurately depict the changes in face displacement.
The loading process of the dam, pouring of face, and impoundment is carried out per the actual construction and impoundment. The bottom of the model is fixed, the left and right banks are treated to be simply supported in the X direction, the upstream and downstream are also simply supported in the Y direction, and the other boundaries are free.
The Duncan-Chang [24] model is used to model the stress-strain relationship of the material of the Gongboxia concrete-face rockfill dam, with parameters acquired from laboratory tests as presented in Table 1; the Goodman unit [25] without thickness is used to represent the contact effect between the face and cushion, with parameters as listed in Table 2; and a linear elastic model is used for the concrete face, toe board, wave wall, and bedrock [26]. The main parameters for concrete are elastic modulus E = 28 GPα, Poisson’s ratio μ = 0.2, and density . The bedrock has an elastic modulus E = 10 GPα, Poisson’s ratio μ = 0.25, and density .
3.1.2. Face Deformation Results
The distribution of the maximum displacement of the Gongboxia concrete-face rockfill dam during the impoundment period is shown in Figure 8. The dam undergoes a smaller displacement on the upstream side after impoundment with a maximum value of 0.0022 m that occurs at approximately 1/3 height of the dam from the bottom, compared with a maximum displacement of 0.4384 m on the downstream side that occurs at approximately 1/2 height of the dam from the bottom. The maximum displacement of the dam during impoundment is 0.8664 m at approximately 1/2 height of the dam from the bottom near the downstream side. The maximum displacement accounts for 0.655% of the total dam height and the distribution of the maximum displacement is in accordance with the general displacement distribution of the face rockfill dam.
The calculated and monitored face deformations after impoundment are shown in Figure 9. The face moves towards the downstream direction due to the water pressure and settlement of the rockfill. The calculated deflections at elevations of 1950 m and 1980 m are 8.80 cm and 10.95 cm, respectively, whereas the respective monitored values are 10.08 cm and 11.70 cm, indicating that the respective calculated values are 85.45% and 92.15% of the monitored values. The numerical analysis results show that the face deformation pattern is convex in the upper-middle part and concave in the lower part, which is in good agreement with monitored data. Rheology and humidification have not been considered in the analysis, resulting in small face deformations. Since the output of the inclinometer is greatly affected by temperature, using a fitting integral to treat the monitored data leads to a large cumulative error. For these reasons, the results for the deformation of the concrete face using the numerical model match the results from the monitored data. From the calculations, downstream of the location where the maximum deflection of the face occurs has a large tensile stress, which implies that it is the critical location for the initiation of cracking. The downstream side of the face, which is subjected to large bending moments, is more susceptible to cracking due to large deflections.
3.2. Structural Crack Propagation in the Concrete Face and Its Evolution
3.2.1. Crack Numerical Model
The location where the maximum deformation of the Gongboxia concrete face occurs is selected for the analysis of cracking (Figure 9), and a new numerical model for simulating the cracking of the concrete face is shown in Figure 10. The contact point of the face bottom and cushion is chosen as the origin, the positive X-axis is along the diagonal from the face to the dam crest, the positive Z-axis is from the bottom to the top of the dam, and the positive Y-axis is along the direction of the thickness towards the water. To consider the structural crack at the location where the maximum deformation caused by displacement occurs, the diagonal length of the model is selected as 5.7 m, and the thickness of the model is the actual thickness of the face (0.5 m). The axial width of the face is 2.0 m. The model consists of 45600 eight-node hexahedral elements and 51865 nodes. The concrete face is made of C25 concrete with tensile and compression strengths of 1.3 MPa and 12.5 MPa, respectively. A linear elastic model is used to model the concrete face, and the elastic modulus, Poisson’s ratio, and density are 28 GPa, 0.2, and 2450 kg/m3, respectively.
(a)
(b)
Based on the calculated face displacement, an equivalent displacement is considered to act on the model. The upper and lower displacement boundaries are fixed in the Z direction, the displacements near the dam crest and bottom are 0.459 m and 0.464 m, respectively, and all the other boundaries are free. The displacement in the Z direction after impoundment is represented by two equations: Z = 0.464 + 0.003X (0 < X < 2.85 m) and Z = 0.484 – 0.004X (2.85 m < X < 5.70 m). The displacement is then imposed on the concrete face as an external displacement to simulate the stress state of the concrete face after impoundment. Based on XFEM, the crack displacement is set as the face cracking pattern, and the maximum stress loss is set as 1.3 MPa, considering the tensile strength of concrete.
3.2.2. Crack Propagation Analysis
The maximum principal stress nephograms and crack distribution curves during the cracking of the concrete face for four typical loading phases are shown in Figure 11. Stress concentration occurs at 1/3 height of the concrete face model upstream, as the displacement increases. Since the face on the upstream side bears a compressive stress, it does not demonstrate an obvious compression deformation because of the high compression strength of the concrete face. However, the concrete face on the downstream side bears a tensile stress due to the displacement load. When the tensile stress is lower than the concrete tensile strength (i.e., ), the face does not show cracks. When the tensile stress reaches the concrete tensile strength (i.e., ), the face shows the first horizontal crack on the outer surface. With an increase in the displacement load (i.e., ), the crack tip stress concentrates and extends along the direction of the face thickness. At the location of the downstream side of the face where tensile stress concentration occurs, uniform horizontal cracks with the same spacing gradually emerge. With a continuous increase in the displacement load, the number of cracks stops increasing but the cracks gradually increase in depth and width. The data shows that the depths of the cracks perpendicular to Gongboxia concrete face are in the range of 11 to 25 cm and tend to increase and propagate with time. The cracks are wide on the surface but narrow at the bottom, and no through cracks are found. Figure 12 displays the monitored depths and distributions of cracks at the location with the maximum displacement of the Gongboxia concrete face. The numerical simulation results are in good agreement with the monitored data, implying that the numerical model in this study can well simulate the distribution and propagation of cracks in the concrete face.
(a)
(b)
(c)
(d)
Concrete face slab plays an important role in dam seepage control; when cracks occur in the concrete face slab, the rigidity, shear strength, tensile strength, and bending strength of the concrete face slab will decrease, which will lead to stress redistribution and further destruction. At the same time, the water entering concrete panel leads to the corrosion of steel bar, and the increase of corrosion substances volume results in a great tensile force in the concrete around the reinforcement, which leads to the formation of longitudinal fracture. As the cracks propagate until the panel is penetrated, once the reservoir water passes through the concrete face slab and leak into the rockfill, the rockfill will be taken away by water and cause panel void, thus threatening the overall safety of the dam. In view of this, we can reduce the differential settlement of the dam by optimizing the design of the dam structure size, slowing down cracking due to panel void, and optimizing mix proportion of concrete from the aspect of the concrete foundation materials, increase the resistance to cracking, and reduce the crack propagation speed.
3.3. Factors Affecting Structural Crack Propagation in a Concrete Face
This study further investigates the factors affecting structural crack propagation in a concrete face. Based on the numerical simulation results of the crack model, the deepest crack is chosen as the preset crack of the concrete face. An analysis of the effects of the reinforcement layout and concrete cover thickness on the crack propagation and stress intensity factor of the crack tip has been conducted. Two-way distribution reinforcements of one, two, and three layers are successively selected for investigating the effects of reinforcements on crack tip stress intensity; the effect of thickness of a concrete cover with a two-layered distribution of reinforcements on the crack tip stress intensity is also studied. A separated finite element model is utilized to analyze reinforced concrete. The concrete and reinforcement are modeled using C3D8R and T3D2 elements, respectively. The reinforcement and concrete are assumed to have rigid connections, which means that there is no relative movement between them. The reinforcement layout and numerical calculation model for cracking of the reinforced concrete face are shown in Figure 13. The size, load, and properties of concrete for the model are the same as those in Section 3.2. The elastic modulus, yield strength, and diameter of the reinforcement are 210 GPa, 0.21 GPa, and 8 mm, respectively. The size of the mesh should not be too large for crack tip stress analysis and is selected as 0.1 mm × 0.1 mm.
(a)
(b)
Figures 14–17 show the distributions of the maximum principal stresses when no reinforcement and different layers of two-way reinforcement are applied. The results of the numerical calculations indicate that as the displacement increases, the unbalanced stress in the cracked element converts into unbalanced nodal stress acting on the structure. The numerical calculation models the redistribution of stress in this study. Under the displacement load, cracks continue to develop along the original direction, and no cracks develop with a certain area surrounding the cracks. The above-mentioned crack propagation is in good agreement with the real crack propagation of reinforced concrete. As observed, the maximum principal stresses of the crack tip are 1.946 MPa, 1.945 MPa, and 1.945 MPa when one, two, and three layers of reinforcement are applied, respectively. However, the maximum principal stress is 1.947 MPa when no reinforcement is applied. This indicates that the application of reinforcement in two ways can improve the stress state of the crack tip, and more layers result in more improvement. The distribution of reinforcement has a small impact on the crack stress when it is far from the crack, and the layer of distribution reinforcement only needs to satisfy the overall structural stress requirement, given the cost of steel.
Figure 18 shows the relationship between concrete cover thickness and the maximum principal stress for both the crack tip and distribution reinforcement. The two-way distribution reinforcement bears most of the stress in the concrete face. When the reinforcement is close to the crack tip (with a small concrete cover thickness), for an equivalent displacement, the reinforcement deforms due to large tensile stresses, thereby controlling the crack propagation in the concrete to a certain extent. When the reinforcement is far away from the crack tip, the stress caused by the deformation in the reinforcement gradually reduces and the stress field induced by concrete cracking recovers at the crack tip. The cohesion effect of reinforcement on cracks alleviates stress concentrations on crack tips. In general, the cover thickness of the concrete face has a substantial effect on the reinforcement stress but only a small effect on the crack tip stress intensity. When a double-layered two-way reinforcement is added to the concrete face, the stress intensity of the face crack tip is reduced and the anticracking capacity of the face increases when a reasonable cover thickness is chosen. Considering the interaction between the reinforcement and concrete, a separated crack model treats concrete as a different element. The calculation results show that the separated crack model can simulate cracking of a reinforced concrete face with a satisfactory outcome but its efficiency is low and requires substantial computing resources, making its application questionable.
4. Conclusions
Based on the verification of XFEM in obtaining the stress intensity factor of crack tips being in good agreement with the theoretical value, this study focuses on the Gongboxia concrete-face rockfill dam project and combines XFEM and ABAQUS to investigate face cracking and crack propagation in concrete-face rockfill dams. The results show that, for a displacement load, the first horizontal crack occurs on the downstream side of the face when the tensile stress in concrete face reaches the tensile strength; as the displacement load increases, the crack propagates along the direction of the width of the face, and uniform horizontal cracks with the same spacing emerge gradually at the location of tensile stress concentration. The stresses on the concrete face are concentrated in a two-way distribution of reinforcement, which improves the stress state of the crack tip, and more layers of reinforcement result in more improvement. The thickness of the concrete cover has a great impact on the stress state of the reinforcement but a small impact on the crack tip stress. Hence, a reasonable concrete cover thickness can reduce the stress intensity of face crack tips and heighten the anticracking capacity of the face.
Data Availability
The data used to support the findings of this study are available from the corresponding author upon request.
Conflicts of Interest
The authors declare no conflicts of interest.
Acknowledgments
The research described in this paper was funded by the National Natural Science Foundation for Excellent Young Scientists of China (Grant no. 51722907) and the National Natural Science Foundation of China (Grant no. 51979224).