Numerical Simulation of Hydraulic Fracturing in Earth and Rockfill Dam Using Extended Finite Element Method

Geotechnical Engineering Department, Nanjing Hydraulic Research Institute, Nanjing 210024, China Key Laboratory of Failure Mechanism and Safety Control Techniques of Earth-Rock Dam of the Ministry of Water Resources, Nanjing Hydraulic Research Institute, Nanjing 210029, China Key Laboratory of Ministry of Education for Geomechanics and Embankment Engineering, Hohai University, Nanjing 210098, China


Introduction
In recent years, with the needs of hydropower construction, hydropower energy development is increasing rapidly in China.erefore, many high earth and rockfill dams have been designed or under construction, such as the Changhe (240 m in height) in Sichuan province, Shiziping (136 m in height) in Hebei province, and Nuozhadu (261.5 m in height) in Yunnan province.In addition, a number of super high (up to 300 m in height) earth and rockfill dam have been built, such as Lianghekou (295 m in height) in Shanxi province and Shuangjiangkou (322 m in height) in Guizhou province.With the rapid development of earth and rockfill dam, many accidents have occurred.Although no devastating disaster is caused by the accidents, safety operation of the dam is seriously affected.Consequently, many researchers have gradually paid attention to the safety design of high earth and rockfill dams.
Hydraulic fracturing is one of the most important factors affecting the safety of earth and rockfill dam.A series of studies have been carried out to simulate hydraulic fracturing.Ng and Small [1] presented a method of predicting hydraulic fracturing by using special joint elements that allow fluid flow and fracture to be modeled.Secchi and Schrefler [2] proposed a cohesive fracture model to simulate crack propagation of a concrete dam, which needs mesh updating of the crack tip.Barani and Khoei [3] presented a 3D cohesive crack propagation method in saturated porous media.By using double-nodded zero-thickness cohesive interface elements, the fracture behavior was simulated.Wang et al. [4] thought that the mechanism of hydraulic fracturing can be explained by fracture mechanics, and the crack propagation may follow any of mode I, mode II, and mixed mode I-II.Liu et al. [5] presented a stabilized extended finite element formulation to simulate the hydraulic fracturing process in an elastoplastic medium.e fracture propagation process was governed by a cohesive fracture model.Li et al. [6] used the smeared cracking theory to simulate the process of hydraulic fracturing in a certain earth and rock ll dams.Tang et al. [7][8][9] gave di erent mechanical parameters to di erent material units by using a certain statistical distribution.us, the crack propagation can be simulated.In general, the existing simulation methods require speci c means, like grid re nement and estimate the position and direction of the fracture or set cohesive element, which result in the large calculating quantity and complex simulation process.
In 1999, Belytschko and Black [10] and Moës et al. [11] presented the extended nite element method (XFEM).e unique advantages of simulating crack propagation in homogenous material such as rock and concrete materials have been demonstrated [12].By using the XFEM, the crack is allowed to pass through the element; thus, cracks with complex shapes can be calculated in regular mesh without remeshing.Because of these advantages, the XFEM has been considered as a more accurate and e cient method in many areas, such as dynamic crack propagation, shear zone evolution, and multiphase ow.Yet, no researchers have used the XFEM to simulate crack problem in earth and rock ll dam engineering.In theory, it is possible to simulate hydraulic fracturing in earth and rock ll dam by using the XFEM.
In this paper, the XFEM is used to simulate the hydraulic fracturing in an actual high earth and rock ll dam.e possibility of hydraulic fracturing occurrence is analyzed, and the critical crack length is obtained when hydraulic fracturing occurs.
en, the crack propagation path and length is analyzed by inserting di erent length initial cracks at di erent elevation.Furthermore, based on the results, the possibility of crack penetration upstream and downstream of the core is analyzed.Related research achievements are of great signi cance for the prevention and control of hydraulic fracturing in actual earth and rock ll dam construction.

The Mechanism of Hydraulic Fracturing
A number of researchers did some meaningful works toward the mechanism of hydraulic fracturing in the core of earth and rock ll dam [13].e "water wedging" action has been approved by many engineers, which means that the core is not completely homogeneous and there exists di erent planes of weakness during the construction of the dam.us, the "permeable weak surfaces/bodies" occur at the upstream surface of the core as can be seen by ABC area in Figure 1. e permeability of these weak areas is higher than that of the soil around these weakens which allows water to enter into these weak areas rapidly as the reservoir water level rises.Consequently, there exists a directional water pressure on the AC surface and a downward water pressure on the BC surface as shown in Figure 1.If the water pressure is large enough and the core wall has poor anticrack ability, cracks may occur at the point C. Furthermore, the entry of high-pressure water will cause crack propagation, which will result in hydraulic fracturing.

Simulation Method of XFEM
3.1.Displacement Mode.Based on the partition of the unity method, the XFEM displacement formula can be derived by introducing an additional function which re ects the local characteristic of the crack on the basis of the conventional FEM displacement mode.
e displacement formula can be expressed as where u I is the nodal displacement vector, α I is the nodal enriched degree of freedom vector for crack penetration elements, H(x) is the jump function which represents the gap between the crack surfaces, b i I is the nodal enriched degree of freedom vector for crack tip elements, and ψ i (x) is the crack tip functions which represent the crack tip singularity.
e jump function H(x) for a general crack is de ned as e crack tip enrichment function for an isotropic elastic material is 2 Advances in Civil Engineering e Application of Water Pressure.e hydrostatic pressure in the crack is applied on the crack surface as a facial pressure.As the crack expands and propagates continuously, the upstream pressure water penetrates into the crack tip slowly.
erefore, the hydrostatic pressure needs to be updated with the propagation of the crack in real time.In this paper, to simplify the problem, it is assumed that the water pressure is evenly distributed on the surface of the crack and the application of water pressure is carried out for each crack segment.According to the hydrostatic water pressure at different elevations, the water pressure is applied on the nodes with additional degrees of freedom (DOFs) in the crack segments.

Critical Water Pressure Inducing Hydraulic Fracturing.
Investigators carried out a great deal of works on hydraulic fracturing criterion.An elastic-plastic solution to fracture initiation pressure is derived based on the Mohr-Coulomb shear failure criterion and the theory of cavity expansion, which can be given by [14] the following equation: where η � (r a /r c ) 2 and M � 3η + (2/η) − 3η sin φ + 7 sin φ − 1.Actually, critical water pressure p f is the minimum value of u a .Based on the theory of expansion of a hollow cylinder in finite length, the radius of elastic-plastic zone boundary is obtained as follows [15]: where Y � (2c cos φ)/(1 − sin φ), α � (1 + sin ϕ)/(1 − sin ϕ), c is the cohesion, φ is the internal friction angle, r a and r b are internal and external radii of the hollow cylinder, respectively, r c is the radius of the plastic zone, p 0 is the radial stress of the external boundary, and p is the water pressure applying on the internal surface of cavity.

Damage Evolution Criterion.
e damage evolution is determined by the definition of damage variable D. e expression of the normal and tangential stress components affected by the damage can be defined as where D is the damage variable, which represents the average damage value between the two cracks and varies between 0 and 1, t n , t s , and t t are the normal stress component and two tangential stress components, respectively, T n is the normal stress component of the unit, and T s and T t are the first and second tangential stress components of the unit, respectively.

Simulation of Hydraulic Fracturing in Earth and Rockfill Dam
4.1.Mesh, Constitutive Model, and Parameters.To obtain the length of crack propagation accurately, the local mesh refinement method is conducted at different elevation of the initial crack.e average length of the element is about 0.5 m after refined.Figure 2 shows the schematic diagram of the mesh with an initial crack at 1/5 height of the core wall.e bottom of the model is constrained in three directions, while the sides of the model are constrained in the horizontal direction.e vertical direction is free to allow the settlement of the model.To simulate the construction of the dam, the load applied on the model is divided into ten stages.
Quadrilateral plane strain elements are used in the mesh, and a small number of triangular elements are used for transition.In addition, the reduced integration method is used in the XFEM enrichment area, and the convergence is controlled by the hourglass.Differently, the full integration method is used in other areas, which is similar to the conventional finite element method.
e Duncan-Chang E-v nonlinear constitutive model is used in this paper.e constitutive parameters are shown in Table 1.
e fracture energy is 4.8 N/m and the tensile strength is 47 kPa of the core material obtained from uniaxial tensile test.

Research Framework and Scheme.
Firstly, the possibility of hydraulic fracturing occurrence without initial crack in the core is analyzed by comparing the water pressure on the upstream with the calculated values p f using ( 4). e criteria can be implemented in the XFEM platform automatically en, apply "permeable weak surfaces" (assumed as initial cracks) of di erent lengths at di erent elevations in the core.Similarly, ( 4) is used to determine the occurrence of hydraulic fracturing.If the calculated critical water pressure is less than the water pressure at the crack tip, then the new crack is determined to appear.e propagation of the initial crack is determined by whether the accumulated energy of the crack propagation is more than the fracture energy of the core material.us, the critical crack length under the full reservoir water level can be obtained by establishing the relationship between the initial crack length and the reservoir water level when the hydraulic fracturing occurs.Finally, according to the initial cracks of di erent lengths at di erent elevation, the crack propagation path and length are studied by increasing the critical water level to the full reservoir water level gradually.e calculation scheme is shown in Table 2, in which the reservoir water level is applied in two stages.If the water level is less than the maximum reservoir water level (697 m), then the water pressure on the crack surface is equal to the pressure on the surface of the core at the same elevation.If the water level is greater than the reservoir water level, the water pressure on the surface of the core wall is the maximum reservoir water level while the water level on the crack surface is applied according to the actual water level.e main purpose of this method is to research the critical length of initial crack under full water level.

e Possibility of Hydraulic Fracturing Occurrence.
Figure 3 shows the relationship between the initial fracturing pressure p f (calculated from formula (4)) and the hydrostatic pressure at the elevation where the initial crack locates.It can be seen that the initial fracturing pressure calculated at each elevation is greater than the hydrostatic pressure at its corresponding elevation.e result

Advances in Civil Engineering
indicates that hydraulic fracturing will not occur in the earth and rockfill dam in this paper without permeable weak surfaces.
As the reservoir water level rises, the pressure water enters into the crack; thus, "water wedge effect" appears [16].
e water pressure on the crack surface increases with the  Advances in Civil Engineering increasing water level.Consequently, hydraulic fracturing is likely to occur.e larger the length of the initial crack and the water pressure is, the stronger the water wedge e ect is.In other words, there exists a critical length of the initial crack for an earth and rock ll dam under the normal water level.
To obtain the critical length of the initial crack, the relationship between the length of initial crack and the water level required for the propagation of the initial crack is established in Figure 4.
It can be seen that under the condition of full reservoir water level, the critical length of the initial crack required for hydraulic fracturing is 5.3 m, 5.9 m, and 7.9 m at the elevation of 507 m, 577 m, and 649 m, respectively.e results illustrate that the minimum length of the initial crack required for hydraulic fracturing is 5.3 m.Since the water pressure on the bottom is the largest, the possibility of hydraulic fracturing occurrence at the bottom of the core wall is greater than that of the upper part.

Analysis of the Crack Propagation.
In the dam engineering, it is meaningful to analyze the path and length of the crack propagation when hydraulic fracturing happens under extreme conditions.Figures 5-7 show the nephogram of crack propagation with di erent lengths of initial crack at di erent elevation.XFEMSTATUS is a variable to describe crack propagation and changes between 0 and 1 and 1 represents that the crack is completely open.Considering the ratio of the crack to the whole dam is small, local areas of nephograms are given in order to show the propagation of the cracks clearly.
It is obvious that, under the constant driving of water pressure, the crack propagates to the inner part of core continuously.e larger the length of initial crack is, the longer the crack propagation will be.e direction of crack propagation changes from the horizontal direction to the upstream slowly.e reason is that the water pressure on the lower part of the core is larger than that of the upper part, which leads the core has a tendency to de ect.
With the increase of initial crack elevation, the length of crack propagation tends to decrease.For the reason that, the water   Advances in Civil Engineering pressure is decreasing as the elevation grows, when the initial crack is located at 4/5 height of the dam (elevation of 649 m), the propagation length is merely 5 m with a 12 m length initial crack.e relation curves between the length of crack propagation and reservoir water level at 507 m elevation are shown in Figure 8. e reservoir water level is increased to the full reservoir level (697 m) gradually.Obviously, the crack propagation length and the reservoir water level have a linear relationship.e longer the length of initial crack, the greater the degree of crack propagation is.When the initial crack length is 10 m, the crack propagation reaches a maximum value, 14.7 m, which accounts for 15% of the width of the core wall.erefore, the crack propagation is not enough to penetrate the core wall.
Figures 9 and 10 show the relation curves of the length of crack propagation versus the reservoir water level at the initial crack elevation of 577 m and 649 m. e reservoir water level required for the crack propagation is relatively close to the full reservoir water level, and the change of water pressure is relatively small.erefore, the crack propagation length is longer than that in the lower part of the core wall.However, it is noteworthy that the width of the core wall decreases with the increasing elevation, which is merely 28.6 m at the elevation of 649 m. e critical length of crack is 7.9 m, which accounts for 27.6% of the core width.
As a consequent, although the occurrence probability of hydraulic fracturing at the bottom of the core wall is larger, the possibility of hydraulic fracturing to penetrate the upstream and downstream is extremely high when the upper part of the core wall reaches the critical crack length.

Conclusions
On the basis of the research on the possibility of hydraulic fracturing occurrence in the earth and rock ll dam, the hydraulic fracturing behavior is simulated by using the extended nite element method.
(1) For the calculation model in this paper, hydraulic fracturing will not occur without the permeable weak surface (initial crack), and the critical initial crack length required for hydraulic fracturing is 5.3 m. (2) e propagation length of hydraulic fracturing decreases with the increase of elevation, and the average propagation length decreases from 9.4 m to 3.4 m. (3) After the water level gradually increased to full reservoir level, there is no crack penetration through the core wall of all calculation schemes.e direction of crack propagation has a certain angle with the horizontal plane toward the downstream.(4) Considering the up-narrow and down-wide types of the core wall, although the occurrence probability of hydraulic fracturing at the bottom of the core wall is larger, the possibility of hydraulic fracturing to penetrate the core is extremely high when the upper part of the core wall reaches the critical crack length.

Figure 1 :
Figure 1: Sketch diagram of hydraulic fracturing in the earth and rock ll dam.

Figure 3 :Figure 4 :
Figure 3: Relation curves of the initial fracturing pressure versus elevation.

Figure 5 :Figure 6 :
Figure 5: Crack propagation nephogram at 507 m elevation.(a) Initial crack with a length of 6 m.(b) Initial crack with a length of 8 m.(c) Initial crack with a length of 10 m.

Figure 7 :Figure 8 :
Figure 7: Crack propagation nephogram at 649 m elevation.(a) Initial crack with a length of 8 m.(b) Initial crack with a length of 10 m.(c) Initial crack with a length of 12 m.

Table 1 :
Constitutive parameters of the materials in the earth and rock ll dam.

Table 2 :
Calculation schemes of hydraulic fracturing.