Finite element analysis for identifying locations of cracking and hydraulic fracturing in homogeneous earthen dams

Abstract

This paper reports a finite element study to identify the locations of crack initiation in a homogeneous earthen dam at its post-construction and reservoir operation stages. The steady state and transient analyses, including reservoir rise-up and drawdown conditions, are simulated to identify the favorable conditions and locations of crack generation. The behavioral response of the dam is represented in terms of the developed total stress, horizontal deformations of the faces, strain accumulation in the dam, and the differential settlement of the dam base. The locations of post-construction cracks are identified based on the negative minor principal stresses developed on the dam faces. For both single-lift and multiple-lift modeling techniques, the upstream face of the dam is found to be the most favorable location for the crack generation. For transient reservoir operations (rise-up and drawdown scenarios), it is identified that hydraulic fracturing may occur on either faces of the dam at specific heights, governed by the minor principal stresses becoming lesser than the developed pore-water pressure. Depending on whether the reservoir drawdown occurs before or after the attainment of steady-state phreatic level within the body of the dam, the pore-water pressure distribution within the dam are found notably different. This results in hydraulic fracturing occurring at different faces of the dam and at different heights. It is important to have a thorough understanding of the tentative location of the cracks developed in homogeneous earthen dams so that proper mitigation measure can be adopted as per requirement.

Introduction

The detrimental effects of cracks on earthen dams and embankments was first observed by Casagrande [1], wherein the importance of dam safety, even with the occurrence of cracks, was highlighted. The different causes of crack formation in earthen dams was summarized by Sherard [2]. As per the existing literature, the cracks in earthen dams can be broadly classified as ‘longitudinal’ cracks and ‘transverse’ cracks (that open parallel or transverse to the dam axis, respectively). They can further be classified as interior or exterior cracks, which may form in vertical, horizontal or intermediate planes under different circumstances. The primary causes of crack formation in earthen dams are namely: (a) differential settlement resulting from comprising elements with different deformability characteristics, (b) hydro-mechanical forces causing redistribution of stresses in the dam during rapid filling and emptying of the reservoir, thereby promoting the formation of loosened or softened zones in the body of an earthen dam, (c) presence of soils in the body and core of earthen dams that has a potential to undergo piping, (d) foundation bed comprising compressible soils, and (e) marked changes in the topography of the footing and side abutments of the dam. There are instances where cracking is initiated by seismic forces as well. Most often, the cracks are formed by a combination of several of the aforementioned causes. The occurrence of cracks in earthen dams and embankments can lead to progressive failure, thereby resulting in catastrophic consequences and hazard to life and property. Thus, a detailed study to reliably identify possible locations of crack generation is the need of the hour, which would aid in developing preventive and healing measures to combat the impending diasasters. Such understanding can help to adopt suitable mitigation techniques and extend the performance life of an earthen dam. Finite element method has already proven its potential to analyze cracking in dams, as has been used by different researchers in the past [3,4,5,6]. The presence and existence of cracks in earthen dams is not considered in the conventional designs that mostly realizes the system as an intact fill. In reality, the formation of cracks during the construction of dams is inevitable. The existence of cracks increases the susceptibility of the soil slopes to erosion, shear strength loss, water seepage, and consequent failures. Millions of dollars are spent every year in the repair and maintenance of these slopes all around the world [7]. The issue of cracking in dams has attracted the attention of many researchers in the past. With the aid of numerical analysis, the settlement behaviour of concrete-faced rockfill Shuibuya dam was evaluated by Zhou et al. [8] during construction, initial impounding, and 2 years after its initial operation. Roosta and Alizadeh [9] assessed the nonlinear behavior of rockfill material through numerical analyses and laboratory tests. In addition, the effect of impounding process on the overall stability of a high dam was addressed by Luo et al. [10]. Nobari et al. [11] studied the embankment dam behaviour for prediction of arching and cracking potential.

The existing literature hints that it is very important to comprehend the mechanisms of cracking, hydraulic fracturing and the subsequent interstitial flow of water so that unfortunate failures can be avoided. In real practice involving the design of earthen dams, cracking and its effects are seldom considered. Lack of proper information regarding the crack geometry, its formation and its propagation through the body of the dam makes it even more difficult to deal with this problem. Thus, in this study, an attempt has been made to identify the locations of crack initiation during the post-construction and operational stages of a homogeneous earthen dam. In this regard, a study on the cracking and hydraulic fracturing in homogeneous earthen dam is attempted with the aid of finite element analysis, following the propositions of Ohne and Narita [12]. Based on saturated–unsaturated material models for the earthen material, various analyses were conducted, which included the steady-state analyses for the post-construction cracking and transient analyses for incorporating the reservoir drawdown and rise-up conditions to identify hydraulic fracturing. Based on the analyses, the article summarizes the conditions that would possibly lead to cracking and hydraulic fracturing of homogeneous earthen dams, and their possible locations. This understanding is necessary to incorporate in the design of such dams for controlling and preventing the failure of dams.

Cracking and hydraulic fracturing in earthen dams

The possibility of cracking and hydraulic fracturing of a homogeneous earthen dam follows two different scenarios during the reservoir impounding and drawdown conditions [12]. The first scenario relates the cause of cracking in the earthen dam to the differential settlement at the base of the earthen dam after its construction. The crack, once initiated, results in erosion during the reservoir operations. When these differential settlements are dominant at the base of the earthen dam, it results in the development of tensile strains on the surface, or in the interior, of the earthen dam. The tension cracks jack open as the minor principal stress (σ₃, in terms of total stress) tends to decrease locally with the development of tensile strains. The criterion for the possibility of hydraulic fracturing in this case is given by the following condition

$$\sigma_{3} < - p_{t}$$

(1)

where, p_t refers to the tensile strength in terms of total stress. The corresponding stress state is shown in Fig. 1a. It is observed that the decrease in σ₃ results in the enlargement of initial stress circle (represented by Circle I) towards the left, until it finally touches the failure envelop (represented by Circle II) to open up tension cracks. This type of cracking is mostly associated with the post-construction settlement of the earthen dam before reservoir filling takes place.

Fig. 1

Mechanism of hydraulic fracturing in earth dams explained through Mohr’s circles based on a differential settlement b development of seepage-induced pore-water pressure [12]

The second scenario is the case where the developed pore-water pressure governs the situation. As the reservoir filling takes place, the pore-water pressure increases towards the centre of the homogeneous dam. This phenomenon results in the decrease of the minor principal stress ($\sigma_{3}^{^{\prime}}$, in terms of effective stress) to the extent that it becomes equal to the effective tensile strength ($p_{t}^{^{\prime}}$) to jack open the latent cracks. The stress state for this scenario is shown in Fig. 1b. In this case, the initial stress circle (I) shifts towards the left without any change in the diameter, and eventually touches the failure envelope as shown by the circle (II). The criterion for the possibility of hydraulic fracturing in this case is given by the following condition

$$\sigma_{3}^{^{\prime}} < - p_{t}^{^{\prime}}$$

(2)

Hence, although having a similar manifestation, the term ‘cracking’ is used to denote the conditions under total stress (as in post-construction cracking of earthen dam before reservoir operation); whereas, the term ‘hydraulic fracturing’ is used to denote the conditions under the influence of pore-water pressure generation and effective stress (as in hydraulic fracturing of earthen dam during transient reservoir operations).

Finite element modeling

Analyses types used in the present study

Finite element (FE) analysis is conducted for the study using the commercially available software GeoStudio, and its modules Seep/W and Sigma/W. Seep/W is capable of handling the analysis of groundwater seepage and problems related to dissipation of excess pore-water pressure within a porous media such as rock and soil. Sigma/W is a powerful tool to address stress–strain and deformation analysis, with special provision to conduct ‘Coupled stress-PWP analysis’. Different FE models are simulated in Seep/W to conduct the initial steady-state seepage analysis, while establishing the existing pore-water pressures and total head conditions in the earthen dam. The initial stress conditions are generated using ‘In-situ analysis’ of Sigma/W. The construction of earthen dam (represented by placing the fill material in single lift or sequential multiple lifts) is simulated using the ‘Load/Deformation analysis’ in Sigma/W. Eventually, the ‘Coupled Stress/PWP analysis’ is used to simulate the ‘rapid rise-up’ and ‘rapid drawdown’ conditions of the reservoir. In both types of reservoir operation, the pore-water pressure generated in the ‘steady-state analysis’ of Seep/W and stresses generated in the ‘In-situ analysis’ of Sigma/W are used to define the initial conditions for the transient reservoir conditions.

Discretization and meshing

Discretization or meshing is one of the most fundamental aspects of finite element modeling, which is used to classify any given region, surface, or line, into smaller parts to address its global response based on the localized responses. Based on the available default ‘fully automatic’ scheme of meshing, there is no need to draw individual ‘finite elements’. A default mesh would be generated once the geometry regions are specified. However, the generated default mesh can be altered, globally or locally, as per necessity. Each of the finite elements comprises ‘nodes’ that are used to describe the distribution of the primary unknowns within the element, and thus aid in developing the finite element equations. In a finite element formulation, it is necessary to adopt a model describing the linear or nonlinear distribution of the primary variable within the element (e.g., total head). For a linear distribution of the primary unknown, the nodes are required only at the corners of the element. With three or more nodes defined along an edge of an element, a quadratic or any other higher order equation can be utilized for describing the distribution of the primary unknown. The compatibility between different elements are established through their common nodes, thereby considering a common distribution of the primary unknown along the common edge of the elements. A special integer-based algorithm is also employed in Geostudio to check the compatibility between regions. The robustness of the meshing algorithm influences the number of divisions along the edge of a region that may be actually higher than what was specified as ‘default’. There are different finite element mesh patterns available as default in Geostudio, namely (a) Quads and Triangles (b) Triangles only (c) Rectangular grid of Quads, and (d) Triangular grid of Quads/Triangles. In general, the global element size (e.g., 0.5) is used for the main regions, while refined element size (expressed as ‘ratio of the default size’) is used along the intersections of different regions with different material constituents. The total number of nodes and elements used in the finite element model depends on the number of mesh elements used and the adopted refinements. In the present study, the model domain is discretized into an unstructured mesh, comprising first order quadrilateral and triangular elements. Based on the geometry and adopted refinements, the domain comprises 1636 nodes and 1531 elements. The global element size used for the basic model is 1 m, while the local refinements yielded the mesh elements having size 0.3 m to the least. Figure 2 shows a typical model of a homogeneous earthen dam resting on its foundation and highlighting the adopted meshing and discretization.

Fig. 2

A typical model of the homogeneous earthen dam resting on its foundation and highlighting the adopted meshing and discretization

Boundary conditions

The boundary conditions are the driving force in any numerical analysis, and the solutions to the numerical problems are considered as a direct response to the applied boundary conditions. The boundary conditions may themselves be predefined or transient that evolve with the iterative solution. Due to their extreme importance, it is essential to have a clear understanding of the physical significance of the various types of boundary condition used in a model. In GeoStudio, all the boundary conditions must be applied directly on geometry items such as region faces, region lines, free lines or free points.

In Seep/W (used for seepage analysis), the total hydraulic head at each node is the primary unknown that is computed relative to the total hydraulic head (H) and/or the flow quantities (Q) predefined or specified at some nodes in the domain. The specified magnitudes of H or Q serve as the boundary conditions. In a steady-state analysis, at least one node in the entire mesh must have a specified H. Seep/W facilitates applying different boundary conditions such as total head (H), total flux (Q), unit flux (q), unit gradient (i) and pressure head (P). In the present study, either the constant total head (H) boundary condition or total head boundary varying with time (H,t) is used to specify the steady state or fluctuating reservoir levels, respectively. Figure 3 exhibits a typical application of the ‘head’ boundary condition for a seepage analysis through a homogeneous earthen dam. Along the upstream end, the full supply level of the reservoir is 24 m and thus a total head boundary of H = 24 m is applied. Along the downstream end of the dam, the water table is at the ground surface. The ground elevation at the bottom of the reservoir being 10 m, a head boundary condition of H = 10 m is applied at the downstream end. Along the downstream face of the dam, a potential seepage face is also flagged as a boundary condition through a zero-flux (Q = 0), which indicates that no additional flux is going to be added or removed at these nodes. This boundary condition ensures that there would be no accumulation or ponding of water on the downstream face, in case the phreatic surface or any flow lines intersect with the slope.

Fig. 3

A typical model of the homogeneous earthen dam with different boundary conditions as applied for Seep/W analysis

In Sigma/W module, multiple boundary condition types are implemented to support different types of load-deformation modelling scenarios. Fundamentally, there are only two types of boundary conditions that can be applied to a stress-deformation analysis; namely, ‘stress’ and ‘force or displacement’ type boundary conditions. In the ‘stress type’ category, the hydrostatic pressure boundary condition (either maintained constant or fluid elevation varying with time) is used to simulate steady or fluctuating reservoir condition, respectively. In all stress-deformation models, it is critical to put a ‘bound’ the problem, which means that some parts of the geometry must be defined as ‘zero displacement’ boundary conditions. Commonly, the left and right far boundaries of a numerical model are restricted from displacing in the horizontal (Cartesian x-) direction, while the bottom edge is restricted from displacing in both horizontal and vertical (Cartesian x–y) directions. Accordingly, the suitable zero-displacement boundary conditions are specified on the boundaries of the domain. The portions of the domain boundaries that should remain free to exhibit displacement are not ‘bound’ in any direction. Further, in order to generate the water-pressure due to the reservoir, a fluid pressure boundary is defined on the upstream face of the dam by specifying a water surface elevation. Sigma/W facilitates both static and transient fluid pressure boundaries to represent steady state and fluctuating reservoir levels respectively. Figure 4 shows a typical representation of the boundary condition adopted for stress-deformation analysis of an earthen dam in Sigma/W.

Fig. 4

A typical of the homogeneous earthen dam with different boundary conditions as applied for Sigma/W analysis

FE models and materials used in the analysis

The flow and deformation responses of homogenous earthen dams are studied to identify the locations of possible initiation of hydraulic fracturing or cracking. The entire set of analyses is divided under two broad categories, Case I and Case II, for representing and simulating cracking during the construction stage and hydraulic fracturing during operational stages of earthen dams. Case I represents those scenarios of crack initiation that are associated with the construction of earthen dam. In general, the field construction of earthen dam is associated with compaction of sequential multiple layers of earth fill. However, in many instances, the numerical models consider the earthen dams to be constructed of single lift, i.e. the dams are considered ‘wish in-place’ during the numerical analysis. In such a case, a full-sized numerical model of the dam is analyzed, which does not represent the field construction scenario. Hence, in this study, analysis pertaining to numerical modeling of Case I considers both single and multiple lifts to represent an earthen dam construction, and assess the influence of commonly adopted single-lift technique on the post-construction cracking of earthen dams. For Case II, the transient reservoir operations (namely the rapid rise-up and rapid drawdown conditions) were simulated to assess their influence on the initiation of hydraulic fracturing of a homogeneous earthen dam. The numerical models used in the present study, along with its relevant dimensions, are shown in Fig. 5.

Fig. 5

Numerical models of earthen dam constructed in a simplified single lift b real-field multiple lifts

The materials of the earthen dam and the foundation are modeled using the ‘elastic–plastic’ material model in Sigma/W to perform the stress-deformation analysis. In Seep/W analysis, the material for earthen dam is modeled using ‘saturated/ unsaturated’ material model, while the foundation is modeled using ‘saturated only’ material model. The volumetric water content function and the hydraulic conductivity function that are used as input for the materials chosen for earthen dam and foundation are shown in Fig. 6. These curves are estimated following the procedures proposed by van Genuchten [13] and are required to solve the seepage equation in the ‘coupled stress-PWP analysis’.

Fig. 6

Hydraulic characteristics of materials for earthen dam and foundation a Volumetric water content function b Hydraulic conductivity function

van Genuchten [13] proposed the following closed form equation to describe the hydraulic conductivity of a soil as a function of matric suction:

$$k_{w} = k_{s} \frac{{\left[ {1 - \left( {a\Psi^{{\left( {n - 1} \right)}} } \right)\left( {1 + \left( {a\Psi^{n} } \right)^{ - m} } \right)} \right]^{2} }}{{\left( {\left( {1 + a\Psi } \right)^{n} } \right)^{m/2} }}$$

(3)

where, k_s is the saturated hydraulic conductivity, Ψ is the required suction range, a and m are the curve fitting parameters, and $n = {1 \mathord{\left/ {\vphantom {1 {\left( {1 - m} \right)}}} \right. \kern-\nulldelimiterspace} {\left( {1 - m} \right)}}$. From the above equations, the hydraulic conductivity function of a soil can be estimated once the saturated conductivity and the two curve fitting parameters, a and m, are known. The curve fitting parameters can be estimated graphically based on the volumetric water content function of the soil. According to van Genuchten [13], the best point to evaluate the curve fitting parameters is the halfway point between the residual and saturated water content of the volumetric water content function. The slope of the function is calculated as:

$$s_{p} = \frac{1}{{\theta_{s} - \theta_{r} }}\left| {\frac{{d\theta_{p} }}{{d\left( {\log \Psi_{p} } \right)}}} \right|$$

(4)

where, θ_s and θ_r is the saturated and residual volumetric water contents respectively, θ_p is the volumetric water content at the halfway point of the volumetric water content function curve, and Ψ_p is the matric suction at the same point.

After the determination of s_p, the following formulae are used to estimate m and a:

$$m = \left\{ {\begin{array}{*{20}c} {1 - e^{{ - 0.8s_{p} }} \quad \,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\quad \forall \,\,\,0 < s_{p} < 1} \\ {1 - \frac{0.5755}{{s_{p} }} + \frac{0.1}{{s_{p}^{2} }} + \frac{0.025}{{s_{p}^{3} }}\quad \,\,\,\,\,\forall \,\,\,s_{p} > 1} \\ \end{array} } \right.$$

(5)

$$a = \frac{1}{\Psi }\left( {2^{1/m} - 1} \right)^{{\left( {1 - n} \right)}}$$

(6)

The various other soil properties used in the modeling are obtained from extensive literature survey [14,15,16,17,18]. The material properties of the soils used in the models are listed in Table 1.

Table 1 Material properties used in the finite element simulation of hydraulic fracturing of homogeneous earthen dam

Interpretations and discussions

Analysis for case I: post-construction cracking of homogeneous earthen dam

The simulation of earthen dam construction is carried out using the ‘Load/Deformation analysis’ in Sigma/W by placing the embankment fill with either single lift or sequential multiple lifts. For the analysis pertaining to multiple lifts, the earthen dam of height 15 m was divided into five equal lifts, the height of each lift being 3 m. The foundation is simulated using ‘elastic–plastic material model’ with the ‘effective drained parameters’, as it remains submerged by the water table during any stage of analysis. The earthen dam is also simulated with the ‘elastic–plastic material model’ while considering ‘total stress parameters’, as this stage of analysis considers only the constructed earthen dam before the commencement of the reservoir filling operation. Hence, the water table remains at the base of the earthen dam during the analysis. The outcomes are presented in terms of the variation and distribution of minor principal stresses, deformations and strain concentrations at various locations in the faces, body of the earthen dam and the foundation soil. It is worth mentioning that GeoStudio 2012 describes the minor principal stress as ‘minimum total stress (σ_min)’, and the same is expressed as

$$\sigma_{\min } = \frac{{\left( {\sigma_{x} + \sigma_{y} } \right)}}{2} - \sqrt {\left[ {\frac{{\left( {\sigma_{x} + \sigma_{y} } \right)}}{2}} \right]^{2} + \tau_{xy}^{2} }$$

(7)

where, σ_x and σ_y are the horizontal and vertical stresses, while τ_xy is the shear stress on any soil element.

For an earthen dam simulated by both single and sequential multiple lifts, Fig. 7 exhibits the variation of minimum total stress along the upstream and downstream face of the dam. As mentioned earlier, conventionally, single-lift modeling is a common technique of finite element simulation of earthen dams. Through this study, an attempt is made to check the verity of the adopted conventional technique, and identify whether there is any marked difference in the response of the dam, in terms of stresses and deformations, if the sequential multiple lift mode of construction is adopted in numerical modeling. Figure 7 clearly indicates that the maximum negative stresses occur along the mid-height of the upstream and downstream faces of the earthen dam. Although for multiple lift simulation, the variation of minimum total stresses along the faces of the dam are more jagged, yet the overall variation of stresses along the faces are largely comparable both in magnitude and in trend. For the earthen dam simulated through single and sequential multiple lift modeling, the maximum negative stress in the upstream face of the dam are recorded to be 26.31 kPa and 19.29 kPa, respectively, at heights of 17.02 m and 17.28 m, respectively. Similarly, along the downstream face, the maximum negative stresses are obtained as 24.61 kPa and 20.5 kPa, respectively, for single and multiple lift modeling, occurring at a height of 17.05 m and 17.28 m, respectively. Tensile stresses obtained from single lift modeling are noted to be marginally higher than the multiple lift modeling. Hence, it can be stated that the sequential multiple lift modeling of the construction of the earthen dam hardly produces any significant difference in the stress characteristics of the dam, and hence, the single lift modeling is an acceptable model for representing an earthen dam. In both the types of modeling, the maximum negative stresses appear at a height of approximately 17–17.3 m, which is approximately 0.5 times the height of the earthen dam. This makes the mid-height of the earthen dam to be the most vulnerable for post-construction crack initiation along the upstream and downstream faces. It is worth mentioning that in this study, the tensile strength of the constituent material is assumed negligible, and hence, generation of negative minor principal stress would be considered as an indicator for the initiation of the cracking in the earthen dam.

Fig. 7

Variation of minimum total stress along the dam faces just after the completion of construction a Upstream face, b Downstream face

During the construction or just after the construction of the homogeneous earthen dam, the effect of reservoir does not come into play and the soil within the dam body is under compression due to its self-weight, the possibility of crack initiation lies along the side slopes of the dam. The horizontal deformation of the earthen dam explains the occurrence of the tensile stresses at its mid-height. Figure 8 displays the horizontal deformation contours in the earthen dam modeled using single lift and sequential multiple lifts. As can be observed, prior to reservoir impounding, owing the assumed symmetric compaction control in the field and symmetric distribution of gravity load, symmetric contours of the horizontal displacements are obtained. The contours clearly reflect that the maximum horizontal displacement occurs at mid-height of the earthen dam, along both the upstream and downstream faces, irrespective of the way the construction of earthen dam being modeled. For the earthen dam modeled using single lift, a maximum horizontal displacement of 50 mm is observed along the upstream and downstream faces, while the same is observed to be approximately 8 mm for the dam modeled using sequential multiple lifts. The horizontal deformation response of the dam is obtained to be slightly higher for modeling with multiple lifts. It is expected that for multiple-lift modeling, the lower layers would show larger overall deformation as multiple lifts results in higher sequential stress transfer to the lower fill layers. For the same purpose, in practice, the compaction of earthen fills is always carried out in thin layers (mostly 3 m as adopted in the present study), so that the compaction stresses are transferred properly to the lower layers and uniform distribution of compaction stresses are achieved. Hence, based on the overall results, it can be aptly stated that the single lift mode of modeling the earthen dam is acceptable enough to represent the stress characteristics within the dam, and hence, further analyses reported in this dissertation follows the stated form of modeling of earthen dam using single-lift technique.

Fig. 8

Variation of horizontal displacements along the faces of the homogeneous earthen dam a Single lift, b Sequential multiple lift

Another important factor that plays a very important role in the determination of crack initiation in a homogeneous earthen dam is the post-construction settlement of the earthen dam. The differential settlement of the dam/foundation may cause detrimental cracking in the dam. Figure 9 displays the settlement along the base of the earthen dam due to the single lift and multiple lift simulations. From the figures, it can be observed that for an earthen dam simulated by single lift technique, the differential settlements at its base are noticeably higher than that obtained from simulation using sequential multiple lift technique. The maximum settlement for the latter case is 10 mm, whereas for the former case, the settlement is observed as 70 mm.

Fig. 9

Settlement profile along the base of the earthen dam due to adoption of single-lift and multiple-lift techniques for simulating the construction scenario

Figure 10a, b exhibit the variation of strains and their concentrations within the earthen dam and its foundation for two different types of adopted modeling scenarios. For the single lift technique, the maximum strain concentration is observed within the foundation, while for the multiple lift scenario, the strain is observed to be concentrated within the body of earthen dam. This is an expected observation as single lift technique considers the dam body to be made of a single entity, and transfers all the bearing stresses generated in the foundation soil. However, for multiple lift technique, each of the sequential lift considers the previous lifts to be the bearing stratum. Hence, for the first lift, the foundation soil acts as the bearing stratum that is subjected to a much lesser load pertaining to a 3 m height of the first fill. However, for the successive lifts, all the previous lifts and foundation are subjected to sequential stress from 3 m lifts. Hence, the lesser stress is transferred to the foundation, and the intermediate fill layers bear more stresses. Consequently, the strain concentration also remains restricted within the body of the dam. Hence, although the stress characteristics along the face of the dam remains approximately similar for the single-lift and multiple-lift modeling, the horizontal deformation and the strain concentrations are observed to be more realistic pertaining to the multiple-lift modeling. The multiple lift modeling, at the same time, is also more realistic of the field construction technique. Further, the results reveal that the multiple-lift technique will exhibit more deformation and strain concentrations, thereby indicating the dam to be more susceptible to cracking due to post-constructional stresses and deformations.

Fig. 10

Maximum strain contours within the earthen dam and its foundation as obtained from a Single-lift technique of simulating construction of earthen dam, b Multiple-lift technique of simulating construction of earthen dam

Hence, it is suggested that the multiple-lift technique would be more appropriate to analyze the post-constructional response of the earthen dam, while the single-lift technique is acceptable for further studies subjected to hydraulic fracturing of earthen dam due to reservoir operations, since the initiation of hydraulic fracturing is more dependent on initial stresses than on the strains and deformations.

Analysis for case II: hydraulic fracturing of homogeneous earthen dam during reservoir rise-up condition

A reservoir rise-up (RU) condition is simulated using the ‘Coupled Stress/PWP analysis’ in Sigma/W module of GeoStudio. In this analysis, the simulation of dam construction is first carried out using the ‘Load/Deformation analysis’, in which the dam construction was simulated in a single lift. At the end of construction, the first impounding of the reservoir is carried out at a rate of 2.8 m/day. Figure 11 shows the rate of reservoir rise-up, in which the reservoir of height 24 m was filled in 5 days, leaving a freeboard of 1 m. ‘Stress and deformation’ analysis is carried out to find out the possibility of crack initiation and its location during the first reservoir filling. The reservoir filling is completed in 5 days, and the analysis is carried out until 1 year after the completion of filling to find out the response of the homogeneous earthen dam during and after the filling process.

Fig. 11

Reservoir rise-up (RU) condition and rate adopted in the present analysis

Figure 12 highlights the variation of the minimum total stress along the upstream and downstream slopes of the earthen dam at different times during and after the reservoir filling operation. It is observed that along the upstream face, as shown in Fig. 12a, the sudden introduction of the water load from the reservoir filling creates a stress inequilibrium and generates a negative stress at the exposed portion of the upstream face. This is particularly observed in the 2nd day when the reservoir filling has reached a height of 5.6 m from the base of the earthen dam. However, with the passage of time, as the reservoir filling progresses towards its completion, on the 5th day, the minimum total stresses on the upstream face increases to the extent that the entire upstream face experiences high positive magnitude of minimum total stress. Further, if the reservoir level is maintained in its completely filled condition for a long time, say 146 days, the saturation front gradually moves towards the downstream of the dam. Based on the low hydraulic conductivity of the material of the earthen dam, the gradual dissipation of excess pore-water pressure towards the body of the dam takes time to release the stress, due to which a decrease in the minimum total stress is observed along the upstream face. The process continues until a steady-state seepage flow is attained, and a marginal change in the pore-water pressure (PWP) profile is observed for further days to come (i.e. 365 days from reservoir filling). In case of the downstream face of the dam, the reservoir filing and steady-state condition has the least effect. The minimum total stresses at some portions of the downstream face will increase only when the phreatic surface intersects with the face. Hence, until then, the minimum total stress will change marginally with the reservoir filling operations, as shown in Fig. 12b. Figure 13 elucidates the variation in the PWP contours, and the progression of the phreatic surface during and after the reservoir filling operation, as is being explained in this section.

Fig. 12

Variation of minimum total stress in the earthen dam during and after the reservoir rise-up a upstream face, b downstream face

Fig. 13

PWP contours and the migration of phreatic surface during and after the reservoir filling operations a 2nd day, b 5th day, c 146th day, d 1 year

Thus, from the stress analysis, it can be observed that during the reservoir filling process (i.e. in the initial days of impounding of reservoir until it is filled), there is a high possibility of hydraulic crack initiation in the upstream face as the stress conditions abruptly changes due to the sudden addition of water load. To examine further the possibility of hydraulic fracturing during the rise-up condition, the pore-water pressure is examined along both the slope faces of the earthen dam. As already coined by Ohne and Narita [12], one criterion for the occurrence or initiation of hydraulic fracture during the reservoir rise-up condition is related to the reduction of effective stress below the effective tensile strength of the earthen dam material, which results in the initiation of cracks. In this study, the identification of the location of hydraulic fracture is carried out by comparing the minimum total stress (i.e. the minor principal stress in terms of total stress) with the pore-water pressure along the dam slopes, as shown in Fig. 14. In this technique, if the PWP exceeds the minor principal stress, the effective minor principal stress becomes negative, thereby inducing a possibility of development of hydraulic fracturing if the magnitude becomes lesser than the effective tensile strength of the constituent material. It is worth mentioning that in this study, the tensile strength of the constituent material is assumed negligible, and hence, generation of negative effective minor principal stress would be considered as an indicator for the initiation of the hydraulic fracturing. From Fig. 14a, it is clearly recognized that in the upstream face of the dam, the pore-water pressure exceeds the minimum total stress values, thereby rendering the effective minor principal stress to be negative and inducing the possibility of hydraulic fracturing along the upstream face.

Fig. 14

Comparison of total minor principal stress with the developed PWP to identify the possible location of hydraulic fracturing along a upstream face, b downstream face

The observation is consistent for both the stages of partial and complete reservoir filling operation. The tentative location of hydraulic fracturing is observed approximately at 0.2–0.25 times the height of the earthen dam, measured from its base. During the reservoir filling operation, it is observed from Fig. 14b that the possibility of hydraulic fracturing does not exist in the downstream face. This is attributed to the fact that the phreatic surface takes long time to reach the downstream face, which occurs only at the steady-state condition. During the reservoir rise-up condition, the pore-water pressure along the downstream face always remain lesser than the minor principal stress, and hence, there is no possibility of negative stresses and consequent initiation of cracking.

Figure 15 highlights the horizontal deformation contours as observed with the progress of reservoir rise-up. The earthen dam exhibits horizontal deformation along with the filling of the reservoir due to the introduction of the water load on the upstream face. Such horizontal displacements might lead to the development of cracks on the crest of the dam as well as facilitate favorable conditions for cracks on the upstream face. If the deformations become considerably large and beyond the tolerance limits of the upstream face, it would result in opening of the latent cracks that might have developed during the construction of the dam. Such a mechanism can lead to significant strain accumulation on the upstream face as shown in Fig. 16. Similar observation is also reported in the past by Khalid et al. [19]. The vertical displacement along the base of the dam is shown in Fig. 17. Due to the reservoir filling and consequent additional load in the upstream side, additional settlement is observed along the upstream side of the dam leading to marginal heaving towards the downstream side. Based on the combined horizontal and vertical deformations of the dam, it is understood that the reservoir filling renders a marginal anticlockwise rotation of the dam.

Fig. 15

Horizontal displacement contours within the earthen dam and foundation during the reservoir rise-up condition at the a 2nd day, b 5th day

Fig. 16

Accumulated maximum strain within the earthen dam and foundation during the reservoir rise-up condition at the a 2nd day, b 5th day

Fig. 17

Vertical displacement profile along the base of the dam during reservoir filling

Analysis for case II: hydraulic fracturing of homogeneous earthen dam during reservoir drawdown condition

The drawdown condition is a classical scenario in slope stability and it is a very common problem encountered by the riverbanks, dams, embankments and levees, when they are subjected to rapidly declining water levels in their upstream side. For earthen dams, the sudden removal of the external hydrostatic pressures may lead to the destabilization of the upstream face and, thus, the drawdown condition is rendered to be critical for the stability of upstream face. There are a number of case histories associated with reservoir drawdown and instability of the upstream face of the earthen dam, as discussed by different researchers [20, 21]. However, it is not customary that every drawdown scenario will lead to failure of the upstream face; rather, the failure is strongly governed by the rate of drawdown along with the shear strength and permeability characteristics of the material comprising the earthen dam, especially near the upstream face. The rapid drawdown, if not failure, might induce tension cracking in the upstream of the earthen dam, which may instigate failure of the earthen dam in successive cycles of reservoir filling and drawdown during its performance period. Hence, in this section, only the influence of reservoir drawdown condition on the possibility of crack initiation is assessed and reported. The reservoir drawdown condition is simulated for two scenarios. The first scenario pertains to the condition when the drawdown is initiated after the reservoir has been in a steady state for a considerable time so that the phreatic surface within the body of the earthen dam attained a steady state. The second scenario pertains to the condition when the drawdown is initiated before the phreatic surface within the body of the dam has attained a steady state (i.e. the transient migration of the phreatic surface is yet to be complete) under the stead-state reservoir level in the upstream of the dam.

Reservoir drawdown after attaining the steady-state phreatic surface

In order to simulate the drawdown scenario from an initial state wherein the phreatic surface within the body of earthen dam is assumed to have attained a steady state, the reservoir is assumed at its highest level for the in-situ steady-state seepage analysis. Thereafter, the reservoir drawdown (DD) is simulated to the study the flow, stress–strain and deformation characteristics of the earthen dam. The rate of reservoir drawdown is considered 2.8 m/day, for which the variation of the reservoir water level with time is shown in Fig. 18.

Fig. 18

Reservoir drawdown (DD) condition and rate adopted in the present analysis where reservoir drawdown is conducted after attaining a steady-state phreatic surface in the earthen dam

Figure 19 depicts the variation of minimum total stress along the upstream and downstream face of the dam. It is observed that along the upstream face, tensile stresses are developed very near to the heel region of the earthen dam. The reduction of minimum total stress initiated as soon as the reservoir drawdown commenced, which can be noted by comparing the minimum total stress profile in the 2nd and 5th day from Fig. 19a. The reduction of minimum total stresses continued until the completion of drawdown at the 5th day. On the other hand, in the downstream face of the dam, the tensile stresses were observed almost along the entire slope face during and after the completion of drawdown process.

Fig. 19

Variation of minimum total stress in the earthen dam pertaining to reservoir drawdown operation initiated after the steady-state phreatic surface is attained a upstream face, b downstream face

Further investigation for the identification of possible hydraulic fracturing and its tentative location during the reservoir drawdown was conducted. Figure 20 presents the comparison of the PWP and the minimum total stress developed along the upstream and downstream faces during partial and complete drawdown of the reservoir. From Fig. 20a, it can be clearly stated that for all conditions, the possibility of hydraulic fracturing along the upstream face is extremely marginal as the PWP profile is observed to be lesser than the minimum principal stresses. If at all hydraulic fracturing occurs, the tentative location would be near the heel, as the PWP is observed to marginally higher than the minimum total stress at this concentrated region. However, as the drawdown attains completion, along the downstream face of the earthen dam (Fig. 20b), the pore-water pressure exceeds the minimum total stress, thereby generating a favorable condition for hydraulic fracturing in the downstream face. The tentative location for possible hydraulic fracturing is identified approximately at 0.5–0.75 times the height earthen dam, measured from its base.

Fig. 20

Comparison of total minor principal stress with the developed PWP to identify the possible location of hydraulic fracturing pertaining to reservoir drawdown operation initiated after the steady-state phreatic surface is attained a upstream face, b downstream face

Figure 21 displays the horizontal deformation contours at the end of reservoir drawdown. It is clearly highlighted that owing to the removal of the stabilizing hydraulic stress from the upstream side, the earthen dam depicts a deformation towards the upstream face. The fact is also supported by the significant accumulation of strain along the upstream face, as shown in Fig. 22. Such significant deformation and strain accumulation towards the upstream side is the triggering phenomenon for upstream face failure of earthen dams due to raid drawdown of the reservoir, although the actual occurrence of instability depends on several factors as stated earlier. Along with the horizontal movements and deformations, if the earthen dam undergoes differential settlement, the outcome will become further favorable for crack initiation on the crest of the dam.

Fig. 21

Horizontal displacement contours within the earthen dam upon completion of reservoir drawdown on the 5th day

Fig. 22

Maximum strain contours within the earthen dam upon completion of reservoir drawdown on the 5th day

Figure 23 presents the vertical displacement profile along the base of the dam during the reservoir drawdown process. The plot reflects the occurrence of differential settlement at the base of the dam, with minor heaving on the upstream end. This behavior is attributed to the sudden removal of the water load, and generation of stress inequilibrium towards the upstream end. The removal of reservoir load would result in a negative settlement, or heaving, to reach a state of stress equilibrium at the end of drawdown. As the drawdown does not have a significant influence on the accumulation of stresses and strains towards the downstream side (as can be seen from Figs. 21 and 22), the differential settlement due to reservoir drawdown is not observed towards the downstream side of the earthen dam.

Fig. 23

Vertical displacement profile along the base of the dam during reservoir drawdown after the steady-state phreatic surface is achieived with the earthen dam

Reservoir drawdown before attaining the steady-state phreatic surface

This case represents a reservoir drawdown scenario much representative of the real-field condition. In this case, the reservoir is assumed empty in the initial stage, i.e. the earthen dam is in its post-construction condition. Further, the reservoir is filled up in a certain number of days. The first filling of the reservoir is carried out in 5 days. Thereafter, the reservoir water level is maintained constant for 2 months. During this period, the transient migration of the phreatic surface through the earthen dam is not complete, and it could not intersect with the downstream face. The adopted steady-state reservoir condition is followed by a drawdown (DD) simulation to study the flow and stress-deformation characteristics of the earthen dam. The rate of reservoir drawdown is considered 2.8 m/day, in which the variation of fluid elevation with time is shown in Fig. 24.

Fig. 24

Reservoir drawdown (DD) condition and rate adopted in the present analysis where reservoir drawdown is conducted before attaining a steady-state phreatic surface in the earthen dam

Figure 25 displays the comparative of the minimum total stress and the PWP profiles along the upstream and downstream faces of the dam, when the drawdown is implemented before the steady-state phreatic surface is attained. It is observed that there is a substantial possibility of hydraulic fracturing along the upstream face, due to a significant release of PWP that accumulated during the reservoir rise-up followed by a constant reservoir level for a certain duration. The tentative location of possible hydraulic fracture along the upstream face is identified approximately within a height of 9 m from the base of the earthen dam, which approximates to 0.6 times the height of the dam measured from the base (Fig. 25a). However, along the downstream face of the dam, Fig. 25b reveals that favorable conditions leading to hydraulic fracturing is not generated. This is attributed to the fact that as the migration of phreatic surface within the earthen dam did not attain a steady state, it did not reach the downstream face within the time for which the reservoir was maintained at constant level. Hence, sufficient pore-water pressure was not generated, and it remained sufficiently lower than the minimum total stress, thereby nullifying the possibility of hydraulic fracturing. However, the minimum total stress is found to reach negative magnitudes at 0.2–0.6 times the height of the dam measured from its base. This scenario can lead to some tension cracks along the downstream face owing to the reservoir drawdown.

Fig. 25

Comparison of total minor principal stress with the developed PWP to identify the possible location of hydraulic fracturing pertaining to reservoir drawdown operation initiated before the steady-state phreatic surface is attained a upstream face, b downstream face

Figure 26 presents the horizontal deformation contours within the earthen dam at the 68th day when the dam is subjected to partial reservoir drawdown. It is observed that owing to the removal of hydraulic load from the upstream face, the earthen dam depicts a movement towards the upstream face. The removal of the hydraulic load also resulted in strain accumulation along the upstream face as shown in Fig. 27. The settlement behaviour of the dam, during drawdown phase, is represented through vertical displacement contours as shown in Fig. 28. It is clearly observed that differential settlement occurred along the base of the dam, accompanied by marginal heaving on the upstream end. The reason for such a response is already explained in the previous section for Fig. 23.

Fig. 26

Horizontal displacement contours within the earthen dam upon partial drawdown of the reservoir on the 68th day

Fig. 27

Maximum strain contours within the earthen dam upon partial drawdown of the reservoir on the 68th day

Fig. 28

Vertical displacement contours within the earthen dam upon partial drawdown of the reservoir on the 68th day

Conclusions

The paper presents the response of homogeneous earthen dam subjected to cracking and hydraulic fracturing. Finite element based seepage analysis and stress analysis are conducted to decipher the responses of the homogeneous earthen dam under various conditions. Based on the existing criteria provided by Ohne and Narita [12], the tentative locations of possible cracking and hydraulic fracturing of homogeneous earthen dam are identified. The influence of single-lift and multiple-lift modeling technique on the response of homogeneous earthen dams on the post-construction cracking are assessed. Based on the negative minor principal stresses generated on the faces of the dam, the tentative locations of post-construction cracking are assessed. The influence of transient reservoir operations, i.e. reservoir rise-up and reservoir drawdown, on the tentative generation of hydraulic fracturing and their possible locations on the dam faces are investigated. Based on the investigations and numerical simulations conducted for the present study, the following conclusions are drawn:

• Numerical modeling of homogeneous earthen dams can be carried out by following the conventional single-lift technique or more detailed realistic multiple-lift technique. It is concluded that both types of modeling technique results in similar stress distribution along the face of the dam. However, sequential multiple-lift technique of modeling uses lifts of lesser thickness that reduces the unrealistic generation of stress-based cracking along the dam faces. In such cases, if cracking occurs, the most probable location is at the mid-height of the slope face.

• During the transient reservoir rise-up condition, the tentative location of hydraulic fracturing along the upstream face of a homogeneous earthen dam is approximately at 0.2–0.25 times the height of the earthen dam, measured from its base. During the reservoir filling operation, the possibility of hydraulic fracturing does not exist in the downstream face.

• The possibility of hydraulic fracturing along the upstream face of a homogeneous earthen dam is extremely marginal during drawdown condition generated after a steady-state phreatic surface is attained. However, the possibility of hydraulic fracturing exists along the downstream face, and the tentative location for possible hydraulic fracturing is approximately at 0.5–0.75 times the height of the earthen dam, measured from its base.

• During the reservoir drawdown condition generated before a steady-state phreatic surface is attained, the tentative location of possible hydraulic fracturing along the upstream face of a homogeneous earthen dam is approximately within 0.6 times the height of the dam, measured from its base. Along the downstream face of the dam, no favorable conditions leading to hydraulic fracturing is generated. However, the minimum total stress attained negative magnitudes at 0.2–0.6 times the height of the dam, measured from its base. The negative magnitudes indicate the possibility of tension crack formation along the downstream face owing to the reservoir drawdown.