# Effect of Normal and Shear Interaction Stiffnesses on Three-Dimensional Viscoplastic Creep Behaviour of a CFR Dam

Authors: Memduh Karalar; Murat Çavuşli
Journal: Advances in Civil Engineering
 (2018)

---

## Abstract

Rockfill materials and foundation continuously interact with each other during lifetime of the rockfill dams. This interaction condition alters the viscoplastic behaviour of these dams in time. For this reason, examination of the time-dependent viscoplastic interaction analyses is vital important for monitoring and evaluating of the future and safety of the rockfill dams. In this study, it is observed how the time-dependent displacement and stress behaviour of a concrete-faced rockfill (CFR) dam change by the effect of the normal and shear interaction spring stiffness parameters. Ilısu Dam that is the longest concrete-faced rockfill dam in the world now and has been completed in the year 2017 is selected for the three-dimensional (3D) creep analyses. The 3D finite difference model of this dam is modelled using FLAC3D software that is based on the finite difference method. The concrete slab, rockfill materials, foundation, and reservoir water are separately created for the 3D interaction analyses. A WIPP-creep viscoplastic material model and a burger-creep viscoplastic material model that are special material models for the creep analyses of rockfill dams are used for concrete slab and for rockfill materials and foundation, respectively. Totally 20 different interaction parameters (normal and shear stiffnesses) are separately defined between the rockfill materials and the foundation to represent the interaction condition. According to numerical analyses, the effect of these various interaction parameters on the viscoplastic behaviour of the Ilısu Dam is evaluated for the empty and full reservoir conditions. As a consequence, the most critical normal and shear stiffnesses’ range for creep analyses of the rockfill dams is determined. Afterwards, the long-term viscoplastic interaction behaviour of Ilısu Dam is examined during 35 years considering this important stiffness values. Settlements, horizontal displacements, and principal stresses are evaluated for both reservoir conditions, and these results are compared with each other in detail.


---

## Body

## 1. Introduction

The concrete-faced rockfill dam (CFRD) is a major type of the rockfill dams, and it consists of transition, cushion, main rockfill, and secondary rockfill zones. In recent times, the construction of CFRDs has increased in many countries in the world due to their adaptability to topography and geology, and short construction period, using locally available materials and cost effectiveness. These important water structures have a great importance for the continuation of the human life and vital needs of people. The deformation-stress behaviour of these huge structures should be monitored continuously due to the large deformations and stresses can occur in the dam body by the effect of the hydrostatic pressure in time. Control of these deformations and stresses is one of the most important technical and scientific problems for the CFRDs, and numerical simulation is an effective approach for these problems [1]. During lifetime of the CFRDs, rockfill materials consistently interact with the foundation [2]. Therefore, while the behaviour of rockfill dams is monitored and analysed, this interaction behaviour should not be ignored. Simulation of the interaction condition is obtained defining interface elements between two discrete surfaces while performing numerical simulation. The construction of interface elements is a necessary and efficient technique for the transition between the foundation and rockfill, simplifying the modelling process and increasing the efficiency and accuracy of the computation [3].Dam body and foundation interaction studies have started many years ago. Westergaard is one of the most important pioneers of these studies. He proposed one of the earliest results of the dam-foundation interaction, under external load effects. The importance of estimating the hydrodynamic pressure on rigid dams was examined by him in 1933 [4]. In the following years, many numerical techniques were developed to describe foundation-dam interaction. These techniques can be classified into the three categories. In the first type, interface regards simply as a surface of soil. So, the analytic or numerical solutions can be derived for the relationship between the corresponding foundation-structure displacement and foundation pressure against the structure under specific conditions. In this type of method, the foundation must be considered as an ideal material which is an oversimplification of the nonlinear behaviour of the concrete-faced rockfill dams. The second type is based on the contact mechanics. Many algorithms have been developed to analyse structure-foundation interaction (e.g., Lagrange methods and penalty function methods). This type can represent the discontinuity behaviour of the interface between the rigid structure bodies or continuums. However, creating a model, which can describe the complex behaviour of the interface, is very difficult for this type. The third type is the interface element. It has been widely used in the interaction analysis of the foundation-structure systems. For this type, an element can model the discontinuity between the foundation and neighbouring structure at the interface within a continuum-based numerical method (e.g., finite difference method and finite element method) [5–9]. The interface elements that are defined between the dam body and foundation are very thin and intensively constrained by the structure, and the interface elements usually take into account the three components. These components are one stress component normal to the interface and two shear stress components tangential to the interface [10–13]. Afterwards, an important interaction technique was proposed in 1979 to predict the settlement of concrete-faced rockfill dams during 10 years [14]. This prediction technique was edited in 1984 [15], and this process has provided significant contributions to the literature. Then, investigators were contributed to developments in interaction analysis and deformation behaviours of the CFRDs [16, 17]. An interaction analysis method was proposed for time-dependent interaction analysis of CFRDs. It was revealed that whether the deformation of the CFR dam body is nonlinear or time dependent [18]. There are very few literatures related to time-dependent viscoplastic analysis of CFRDs considering interaction conditions between dam body and foundation. So, it is aimed at filling these important lacks in this study. The effect of various interface element parameters (normal and shear stiffnesses) on the long-term viscoplastic behaviour of Ilısu CFR dam using special material models is examined. According to the analysis results, stiffness value range for creep analyses of the CFR dams is suggested and the deficiencies in the literature have been largely eliminated.


## 2. Background

### 2.1. Ilısu Dam and Mechanical Properties of Rockfill Materials

The Ilısu Dam is the part of the Southeastern Anatolian Project (GAP), and it is currently the largest hydropower project in Turkey. This huge water structure was completed in the year 2017. The project area is located between Siirt, Batman, and Diyarbakır provinces, and it was built in the boundaries of the Ilısu village. Ilısu Dam’s location is shown in detail in Figure1.Figure 1
The location and view of Ilısu Dam.Ilısu Dam that is the longest concrete-faced rockfill dam in the world has 1775 m crest length. Moreover, this structure is one of the largest rockfill dams in the world now. 44 million·m3 filling material is used for constructing of the dam body. 3 diversion tunnels with a diameter of 12 m and length of 1 km and 3 power tunnels with a diameter of 11 m are constructed in the project. Dam’s precipitation area is 35509 km2. The crest width is 8 m, and the dam body height is 130 m. The lake volume is 10.4 billion·m3. Maximum water elevation is 526.82 m, and the reservoir area is 318.5 km2. The dam generates 3.833 GWh power per year in average. The slopes of the upstream and downstream are 1 : 1.4, and the slopes of the rockfill zones and transition zones are 2 : 1.5. The most critical section and depth changes of the dam body are presented in detail in Figure 2.Figure 2
(a) The typical cross section of Ilısu Dam. (b) Change in the dam body depth along crest axis [19].This dam was built as a concrete-faced rockfill dam, and it was constructed using different rockfill materials. While constructing the dam, rockfill materials were compacted by sheepsfoot rollers. This water construction has 3 different filling materials such as basalt (3B), limestone (3A), and bedding zone (2B) and a concrete slab that was constructed on the surface of the dam body to provide impermeability. All materials have different mechanical properties, and mechanical properties of the materials are selected from the laboratory experiments for interaction analyses as given in Table1.Table 1
Material properties of Ilısu Dam body [19].
Characteristics
Specific weight (g/cm3)
Unit weight (g/cm3)
Porosity (%)
Water content (%)
Air content (%)
Material content (%)
2B
2.74
2.23
18.61
4.05
14.56
81.39
3A
2.68
1.99
25.75
7.78
17.97
74.25
3B
3.01
2.25
25.25
1.50
23.75
74.75


### 2.2. Finite Difference Model of Ilısu Dam

Ilısu Dam that is one of the most important rockfill dams in Turkey is selected for the three-dimensional creep analyses in this study. It is modelled using FLAC3D software to examine the time-dependent viscoplastic behaviour of the dam. The details of the 3D finite difference model are shown in Figure3.Figure 3
3D finite difference model of Ilısu Dam.While modelling this water structure, the concrete slab, rockfill materials, and foundation are modelled as the original project of the Ilısu Dam. The finite difference model of Ilısu Dam has 5 different sections and 4 different blocks. Details of sections and blocks are shown in Figure3. After the dam body is modelled, the foundation is modelled in detail. It is extended toward downstream and the valley side as much as dam height. Also, it is extended three times of the dam height at upstream side of the dam. Finally, the height of the foundation is considered as much as the dam height. Total 1304547 finite difference elements are used in the 3D finite difference model. Special material models are taken into account for the numerical analyses. The burger-creep viscoplastic model is characterized with visco-elastoplastic deviatoric behaviour and elastoplastic volumetric behaviour. Its plastic constitutive laws correspond to the Mohr–Coulomb material model. This viscoplastic model was rarely used to simulate the viscoplastic creep analyses of the water structures in the literature [20]. So, rockfill materials and foundation are taken into account in this study. The Drucker–Prager material model is widely used for frictional materials such as concrete. This material model is the most compatible with the WIPP-reference creep law because both models are formulated with the second invariant of the deviatoric stress tensor [21]. Thus, the WIPP-creep viscoplastic model is taken into account for concrete slab. While defining these material models to FLAC3D, the fish codes which are specially defined to FLAC3D software are used. Interface elements are defined between dam body and foundation to provide the interaction. Totally 20 different interaction parameters are used in the time-dependent numerical analyses. Finally, the reservoir water is modelled for the full reservoir water height (130 m), and the effect of the leakage on the behaviour of the dam is taken into account for all time-dependent analyses [22]. The boundary movements of the 3D model must be restricted before the 3D model is analysed. So, the movement of the bottom of the foundation is restricted in three directions (x, y, and z). Moreover, the movement of the side surfaces of the model is allowed only in the vertical direction (z), and it is restricted in the horizontal directions (x, y). The Ilısu Dam has a great number of elements and nodes. So, the creating and meshing of the 3D model took a long time. Many problems are encountered during viscoplastic numerical analyses because the three-dimensional finite difference model of the Ilısu Dam has a great number of elements and nodes. For this reason, the finite difference mesh is changed several times, and a new mesh is created so that the correct result can be achieved. Totally 6 different mesh widths are tried for 3D analyses. These mesh widths are 50 m, 40 m, 30 m, 20 m, 15 m, and 10 m, respectively. It is clearly seen that when mesh range is selected less than 15 meters, the settlements on the Ilısu Dam surface do not change (Figure 4). So, the mesh width is selected 15 m for numerical analyses.Figure 4
Maximum settlements on the dam body for different mesh widths.


## 2.1. Ilısu Dam and Mechanical Properties of Rockfill Materials

The Ilısu Dam is the part of the Southeastern Anatolian Project (GAP), and it is currently the largest hydropower project in Turkey. This huge water structure was completed in the year 2017. The project area is located between Siirt, Batman, and Diyarbakır provinces, and it was built in the boundaries of the Ilısu village. Ilısu Dam’s location is shown in detail in Figure1.Figure 1
The location and view of Ilısu Dam.Ilısu Dam that is the longest concrete-faced rockfill dam in the world has 1775 m crest length. Moreover, this structure is one of the largest rockfill dams in the world now. 44 million·m3 filling material is used for constructing of the dam body. 3 diversion tunnels with a diameter of 12 m and length of 1 km and 3 power tunnels with a diameter of 11 m are constructed in the project. Dam’s precipitation area is 35509 km2. The crest width is 8 m, and the dam body height is 130 m. The lake volume is 10.4 billion·m3. Maximum water elevation is 526.82 m, and the reservoir area is 318.5 km2. The dam generates 3.833 GWh power per year in average. The slopes of the upstream and downstream are 1 : 1.4, and the slopes of the rockfill zones and transition zones are 2 : 1.5. The most critical section and depth changes of the dam body are presented in detail in Figure 2.Figure 2
(a) The typical cross section of Ilısu Dam. (b) Change in the dam body depth along crest axis [19].This dam was built as a concrete-faced rockfill dam, and it was constructed using different rockfill materials. While constructing the dam, rockfill materials were compacted by sheepsfoot rollers. This water construction has 3 different filling materials such as basalt (3B), limestone (3A), and bedding zone (2B) and a concrete slab that was constructed on the surface of the dam body to provide impermeability. All materials have different mechanical properties, and mechanical properties of the materials are selected from the laboratory experiments for interaction analyses as given in Table1.Table 1
Material properties of Ilısu Dam body [19].
Characteristics
Specific weight (g/cm3)
Unit weight (g/cm3)
Porosity (%)
Water content (%)
Air content (%)
Material content (%)
2B
2.74
2.23
18.61
4.05
14.56
81.39
3A
2.68
1.99
25.75
7.78
17.97
74.25
3B
3.01
2.25
25.25
1.50
23.75
74.75


## 2.2. Finite Difference Model of Ilısu Dam

Ilısu Dam that is one of the most important rockfill dams in Turkey is selected for the three-dimensional creep analyses in this study. It is modelled using FLAC3D software to examine the time-dependent viscoplastic behaviour of the dam. The details of the 3D finite difference model are shown in Figure3.Figure 3
3D finite difference model of Ilısu Dam.While modelling this water structure, the concrete slab, rockfill materials, and foundation are modelled as the original project of the Ilısu Dam. The finite difference model of Ilısu Dam has 5 different sections and 4 different blocks. Details of sections and blocks are shown in Figure3. After the dam body is modelled, the foundation is modelled in detail. It is extended toward downstream and the valley side as much as dam height. Also, it is extended three times of the dam height at upstream side of the dam. Finally, the height of the foundation is considered as much as the dam height. Total 1304547 finite difference elements are used in the 3D finite difference model. Special material models are taken into account for the numerical analyses. The burger-creep viscoplastic model is characterized with visco-elastoplastic deviatoric behaviour and elastoplastic volumetric behaviour. Its plastic constitutive laws correspond to the Mohr–Coulomb material model. This viscoplastic model was rarely used to simulate the viscoplastic creep analyses of the water structures in the literature [20]. So, rockfill materials and foundation are taken into account in this study. The Drucker–Prager material model is widely used for frictional materials such as concrete. This material model is the most compatible with the WIPP-reference creep law because both models are formulated with the second invariant of the deviatoric stress tensor [21]. Thus, the WIPP-creep viscoplastic model is taken into account for concrete slab. While defining these material models to FLAC3D, the fish codes which are specially defined to FLAC3D software are used. Interface elements are defined between dam body and foundation to provide the interaction. Totally 20 different interaction parameters are used in the time-dependent numerical analyses. Finally, the reservoir water is modelled for the full reservoir water height (130 m), and the effect of the leakage on the behaviour of the dam is taken into account for all time-dependent analyses [22]. The boundary movements of the 3D model must be restricted before the 3D model is analysed. So, the movement of the bottom of the foundation is restricted in three directions (x, y, and z). Moreover, the movement of the side surfaces of the model is allowed only in the vertical direction (z), and it is restricted in the horizontal directions (x, y). The Ilısu Dam has a great number of elements and nodes. So, the creating and meshing of the 3D model took a long time. Many problems are encountered during viscoplastic numerical analyses because the three-dimensional finite difference model of the Ilısu Dam has a great number of elements and nodes. For this reason, the finite difference mesh is changed several times, and a new mesh is created so that the correct result can be achieved. Totally 6 different mesh widths are tried for 3D analyses. These mesh widths are 50 m, 40 m, 30 m, 20 m, 15 m, and 10 m, respectively. It is clearly seen that when mesh range is selected less than 15 meters, the settlements on the Ilısu Dam surface do not change (Figure 4). So, the mesh width is selected 15 m for numerical analyses.Figure 4
Maximum settlements on the dam body for different mesh widths.


## 3. Mathematical Formulation

### 3.1. Theoretical Background of the Interaction between Discrete Surfaces

In FLAC3D, the interaction condition is represented defining a normal stiffness value and a shear stiffness value between two discrete planes (e.g., dam and foundation) as seen in Figure5.Figure 5
An interface condition between A and B sides.FLAC3D uses a contact logic which is similar in nature to that used in the different element methods, for either side of the interface [23]. As seen in Figure 5, gridpoint N is checked for contact on the segment between grid points M and P. If contact is detected, the normal vector (n) is computed for the contact gridpoint (N). In addition, a length (L) is defined for the contact at N along the interface. This length is equal to half of the distance to the nearest gridpoint to the left of N, irrespective of whether the neighbouring gridpoint is on the same side of the interface or on the opposite side of N. In this way, the entire interface is divided into contiguous segments, each controlled by a gridpoint. During each time step, the velocity is determined as seen at Equation (1). Since the units of velocity are displacement per time step, the calculation of the time step has been scaled to unity to speed convergence. The incremental displacement for any given time step is(1)Δui≡u˙i.The incremental relative displacement vector at the contact point is resolved into the normal and shear directions, and total normal and shear forces are determined by(2)Fnt+Δt=Fnt−knΔunt+1/2ΔtL,Fst+Δt=Fst−ksΔust+1/2ΔtL,where the unit of kn and ks is calculated.Normal and shear stiffnesses (kn and ks) are not easily measured or well-known parameters. Many methods of estimating joint stiffness have been derived in the past. Two important methods are generally used in the numerical analyses. One of them is based on the deformation properties of the rock mass, and second one is adopted from the properties of the joint infilling material. These methods are explained in detail as shown below.


### 3.2. Calculation of Normal and Shear Stiffnesses considering Rockfill Properties

Normal stiffness and shear stiffness values can be calculated from joint structure in the jointed rock mass and information on the deformability and the deformability of the intact rock. If the jointed rock mass is presumed to have the same deformational response as an equivalent elastic continuum, relations can be derived between jointed rock properties and equivalent continuum properties. The following relation applies for uniaxial-loaded rock which contains a single set of uniformly spaced joints oriented normally to the direction of loading:(3)1Em=1Ei+1knL,where Em is the rock mass modulus, Ei is the intact rock modulus, kn is the joint normal stiffness, and L is the mean joint spacing. Equation (3) can be rearranged to obtain the joint normal stiffness as given in the following equation:(4)kn=EiEmLEi−Em.The same expression can be used to derive a relation for the joint shear stiffness as follows:(5)ks=GiGmLGi−Gm,where Gm is the rock mass shear modulus, Gi is the intact rock shear modulus, and ks is the joint shear stiffness. The equivalent continuum assumption, when extended to three orthogonal joint sets, provides the following relations:(6)Ea=1Ei+1Lakna−1a=1,2,3,Gab=1Gi+1Laksa+1Lbksb−1a,b=1,2,3.Several expressions have been derived for two- and three-dimensional characterizations and multiple joint sets [24–27].


### 3.3. Calculation of Normal and Shear considering Joint Infill Properties

Another approach for estimating joint stiffness assumes that a joint has an infill material with known elastic properties. The stiffness of a joint can be evaluated from the thickness and modulus of the infilling material by the following equation:(7)kn=E0h,ks=G0h,where kn is the joint normal stiffness, ks is the joint shear stiffness, E0 is Young’s modulus of the infill material, G0 is the shear modulus of the infill material, and h is the joint thickness or opening.


### 3.4. The Burger-Creep Viscoplastic Material Model

The burger-creep viscoplastic model that is a special model for time-dependent creep analysis is characterized with an elastoplastic volumetric behaviour and a visco-elastoplastic deviatoric behaviour. The viscoplastic strain rate and viscoelastic components are presumed to act in series. The viscoelastic constitutive law corresponds to the Burger model (the Kelvin and Maxwell components), and the plastic constitutive law corresponds to the Mohr–Coulomb model and the burger model. The symbolsSij and eij are used to denote deviatoric stress and strain components [21]:(8)Sij=σij−σ0δij,eij=∈ij−evol3δij,where(9)σ0=σkk3,evol=∈kk.The Kelvin, Maxwell, and plastic contributions are labeled using the superscriptsK, M, and p, respectively. With those conventions, the model deviatoric behaviour may be described by the following relations.Strain rate partitioning:(10)e˙ij=e˙ijK+e˙ijM+e˙ijp.The Kelvin model is expressed as follows:(11)Sij=2ηKe˙ijK+2GKeijK.The Mohr–Coulomb model is expressed as follows:(12)e˙ijp=λ∗∂g∂σij−13e˙volpδij,e˙volp=λ∗∂g∂σ11+∂g∂σ22+∂g∂σ33.The Maxwell model is expressed as follows:(13)e˙ijM=S˙ij2GM+Sij2ηM.In turn, the volumetric behaviour is given by(14)σ˙0=Ke˙vol−e˙volp,where the properties K and G are the bulk and shear moduli and η is the dynamic viscosity (kinematic viscosity times mass density). The Mohr–Coulomb yield envelope is a composite of shear and tensile criteria. The yield criterion is f=0, and in the principal axes formulation, the following formulation is obtained.Shear yielding:(15)f=σ1−σ3Nφ+2CNφ.Tension yielding:(16)f=σt−σ3,where C is the material cohesion, φ is the friction, Nφ=1+sinφ/1−sinφ, σt is the tensile strength, and σ1 and σ3 are the minimum and maximum principal stresses (compression negative). The potential function g has the following form.Shear failure:(17)g=σ1−σ3Nψ.Tension failure:(18)g=−σ3,where ψ is the material dilation and Nφ=1+sinψ/1−sinψ. Finally, λ∗ is a parameter that is nonzero during plastic flow only, which is determined by application of the plastic yield condition f=0 [21].


### 3.5. WIPP-Creep Viscoplastic Material Model

Viscoplasticity can model by combining the viscoelastic WIPP model with the Drucker–Prager plasticity model. The Drucker–Prager model is the most compatible with the WIPP-reference creep law because both models are formulated in terms of the second invariant of the deviatoric stress tensor. The shear yield function for the Drucker–Prager model is(19)fs=τ+qφσ0−kφ,where fs=0 at yield and σ0=σkk/3, and τJ2 τ, where J2 is the second invariant of the deviatoric stress tensor; the parameters qφ and kφ are the material properties:(20)J2=σijdσijd2,where τ may be related to the stress magnitude, σ¯:(21)σ¯=3τ.The plastic potential function in shear,gs, is similar to the yield function, with the substitution of qψ for qφ as a material property that controls dilation:(22)gs=τ+qψσ0.If the yield conditionfs=0 is met, the following flow rules apply:(23)∈˙ijdp=λ∂gs∂σijd,∈˙op=λ∂gs∂σ0,where λ is a multiplier (not a material property) to be determined from the requirement that the final stress tensor must satisfy the yield condition. The superscript p denotes “plastic” and d denotes “deviatoric.”(24)∈˙ijdp=λσijd2τ,∈˙op=λqψ.In elastic/plastic formulation, these equations are solved simultaneously with the conditionfs=0, and the condition is that the sum of elastic and plastic strain rates must equal the applied strain rate. The Drucker–Prager model also contains a tensile yield surface, with a composite decision function used near the intersection of the shear and tensile yield functions. The tensile yield surface is(25)ft=σ0−σt,where σt is the tensile yield strength. The associated plastic potential function is(26)gt=σ0.Using an approach similar to that used for shear yield, the strain rates for tensile yield are:(27)∈˙ijdp=0,∈˙op=λ,where λ is determined from the condition that ft=0. Note that the tensile strength cannot be greater than the value of mean stress at which fs becomes zero.When both creep and plastic flows occur, it is assumed that the associated strain rates act “in series”:(28)∈˙ijd=∈˙ijde+∈˙ijdv+∈˙ijdp,where the terms represent elastic, viscous, and plastic strain rates, respectively. We first treat the case of shear yield fs>0:(29)∈˙ijd=σ˙ijd2G+σijd2σ¯3∈˙+3λ.In contrast to the creep-only model, the volumetric response of the viscoplastic model is not uncoupled from the deviatoric behaviour unlessqψ=0 [21].(30)∈˙kk=3∈˙o=σ˙kk3K+λqψ.


## 3.1. Theoretical Background of the Interaction between Discrete Surfaces

In FLAC3D, the interaction condition is represented defining a normal stiffness value and a shear stiffness value between two discrete planes (e.g., dam and foundation) as seen in Figure5.Figure 5
An interface condition between A and B sides.FLAC3D uses a contact logic which is similar in nature to that used in the different element methods, for either side of the interface [23]. As seen in Figure 5, gridpoint N is checked for contact on the segment between grid points M and P. If contact is detected, the normal vector (n) is computed for the contact gridpoint (N). In addition, a length (L) is defined for the contact at N along the interface. This length is equal to half of the distance to the nearest gridpoint to the left of N, irrespective of whether the neighbouring gridpoint is on the same side of the interface or on the opposite side of N. In this way, the entire interface is divided into contiguous segments, each controlled by a gridpoint. During each time step, the velocity is determined as seen at Equation (1). Since the units of velocity are displacement per time step, the calculation of the time step has been scaled to unity to speed convergence. The incremental displacement for any given time step is(1)Δui≡u˙i.The incremental relative displacement vector at the contact point is resolved into the normal and shear directions, and total normal and shear forces are determined by(2)Fnt+Δt=Fnt−knΔunt+1/2ΔtL,Fst+Δt=Fst−ksΔust+1/2ΔtL,where the unit of kn and ks is calculated.Normal and shear stiffnesses (kn and ks) are not easily measured or well-known parameters. Many methods of estimating joint stiffness have been derived in the past. Two important methods are generally used in the numerical analyses. One of them is based on the deformation properties of the rock mass, and second one is adopted from the properties of the joint infilling material. These methods are explained in detail as shown below.


## 3.2. Calculation of Normal and Shear Stiffnesses considering Rockfill Properties

Normal stiffness and shear stiffness values can be calculated from joint structure in the jointed rock mass and information on the deformability and the deformability of the intact rock. If the jointed rock mass is presumed to have the same deformational response as an equivalent elastic continuum, relations can be derived between jointed rock properties and equivalent continuum properties. The following relation applies for uniaxial-loaded rock which contains a single set of uniformly spaced joints oriented normally to the direction of loading:(3)1Em=1Ei+1knL,where Em is the rock mass modulus, Ei is the intact rock modulus, kn is the joint normal stiffness, and L is the mean joint spacing. Equation (3) can be rearranged to obtain the joint normal stiffness as given in the following equation:(4)kn=EiEmLEi−Em.The same expression can be used to derive a relation for the joint shear stiffness as follows:(5)ks=GiGmLGi−Gm,where Gm is the rock mass shear modulus, Gi is the intact rock shear modulus, and ks is the joint shear stiffness. The equivalent continuum assumption, when extended to three orthogonal joint sets, provides the following relations:(6)Ea=1Ei+1Lakna−1a=1,2,3,Gab=1Gi+1Laksa+1Lbksb−1a,b=1,2,3.Several expressions have been derived for two- and three-dimensional characterizations and multiple joint sets [24–27].


## 3.3. Calculation of Normal and Shear considering Joint Infill Properties

Another approach for estimating joint stiffness assumes that a joint has an infill material with known elastic properties. The stiffness of a joint can be evaluated from the thickness and modulus of the infilling material by the following equation:(7)kn=E0h,ks=G0h,where kn is the joint normal stiffness, ks is the joint shear stiffness, E0 is Young’s modulus of the infill material, G0 is the shear modulus of the infill material, and h is the joint thickness or opening.


## 3.4. The Burger-Creep Viscoplastic Material Model

The burger-creep viscoplastic model that is a special model for time-dependent creep analysis is characterized with an elastoplastic volumetric behaviour and a visco-elastoplastic deviatoric behaviour. The viscoplastic strain rate and viscoelastic components are presumed to act in series. The viscoelastic constitutive law corresponds to the Burger model (the Kelvin and Maxwell components), and the plastic constitutive law corresponds to the Mohr–Coulomb model and the burger model. The symbolsSij and eij are used to denote deviatoric stress and strain components [21]:(8)Sij=σij−σ0δij,eij=∈ij−evol3δij,where(9)σ0=σkk3,evol=∈kk.The Kelvin, Maxwell, and plastic contributions are labeled using the superscriptsK, M, and p, respectively. With those conventions, the model deviatoric behaviour may be described by the following relations.Strain rate partitioning:(10)e˙ij=e˙ijK+e˙ijM+e˙ijp.The Kelvin model is expressed as follows:(11)Sij=2ηKe˙ijK+2GKeijK.The Mohr–Coulomb model is expressed as follows:(12)e˙ijp=λ∗∂g∂σij−13e˙volpδij,e˙volp=λ∗∂g∂σ11+∂g∂σ22+∂g∂σ33.The Maxwell model is expressed as follows:(13)e˙ijM=S˙ij2GM+Sij2ηM.In turn, the volumetric behaviour is given by(14)σ˙0=Ke˙vol−e˙volp,where the properties K and G are the bulk and shear moduli and η is the dynamic viscosity (kinematic viscosity times mass density). The Mohr–Coulomb yield envelope is a composite of shear and tensile criteria. The yield criterion is f=0, and in the principal axes formulation, the following formulation is obtained.Shear yielding:(15)f=σ1−σ3Nφ+2CNφ.Tension yielding:(16)f=σt−σ3,where C is the material cohesion, φ is the friction, Nφ=1+sinφ/1−sinφ, σt is the tensile strength, and σ1 and σ3 are the minimum and maximum principal stresses (compression negative). The potential function g has the following form.Shear failure:(17)g=σ1−σ3Nψ.Tension failure:(18)g=−σ3,where ψ is the material dilation and Nφ=1+sinψ/1−sinψ. Finally, λ∗ is a parameter that is nonzero during plastic flow only, which is determined by application of the plastic yield condition f=0 [21].


## 3.5. WIPP-Creep Viscoplastic Material Model

Viscoplasticity can model by combining the viscoelastic WIPP model with the Drucker–Prager plasticity model. The Drucker–Prager model is the most compatible with the WIPP-reference creep law because both models are formulated in terms of the second invariant of the deviatoric stress tensor. The shear yield function for the Drucker–Prager model is(19)fs=τ+qφσ0−kφ,where fs=0 at yield and σ0=σkk/3, and τJ2 τ, where J2 is the second invariant of the deviatoric stress tensor; the parameters qφ and kφ are the material properties:(20)J2=σijdσijd2,where τ may be related to the stress magnitude, σ¯:(21)σ¯=3τ.The plastic potential function in shear,gs, is similar to the yield function, with the substitution of qψ for qφ as a material property that controls dilation:(22)gs=τ+qψσ0.If the yield conditionfs=0 is met, the following flow rules apply:(23)∈˙ijdp=λ∂gs∂σijd,∈˙op=λ∂gs∂σ0,where λ is a multiplier (not a material property) to be determined from the requirement that the final stress tensor must satisfy the yield condition. The superscript p denotes “plastic” and d denotes “deviatoric.”(24)∈˙ijdp=λσijd2τ,∈˙op=λqψ.In elastic/plastic formulation, these equations are solved simultaneously with the conditionfs=0, and the condition is that the sum of elastic and plastic strain rates must equal the applied strain rate. The Drucker–Prager model also contains a tensile yield surface, with a composite decision function used near the intersection of the shear and tensile yield functions. The tensile yield surface is(25)ft=σ0−σt,where σt is the tensile yield strength. The associated plastic potential function is(26)gt=σ0.Using an approach similar to that used for shear yield, the strain rates for tensile yield are:(27)∈˙ijdp=0,∈˙op=λ,where λ is determined from the condition that ft=0. Note that the tensile strength cannot be greater than the value of mean stress at which fs becomes zero.When both creep and plastic flows occur, it is assumed that the associated strain rates act “in series”:(28)∈˙ijd=∈˙ijde+∈˙ijdv+∈˙ijdp,where the terms represent elastic, viscous, and plastic strain rates, respectively. We first treat the case of shear yield fs>0:(29)∈˙ijd=σ˙ijd2G+σijd2σ¯3∈˙+3λ.In contrast to the creep-only model, the volumetric response of the viscoplastic model is not uncoupled from the deviatoric behaviour unlessqψ=0 [21].(30)∈˙kk=3∈˙o=σ˙kk3K+λqψ.


## 4. Numerical Results

Rockfill and foundation materials have different mechanical properties, and these materials interact during the dam’s lifetime. So, the examination of the dam body-foundation interaction is very important for the evaluation of the safety and future of the rockfill dams. This interaction condition is provided defining the interface elements between the dam body and foundation while modelling these dams. The mechanical parameters of the interface elements are variable for each dam, and the most important parameters for interaction analyses are the normal and shear spring stiffnesses. Effect of these stiffness parameters on the nonlinear behaviour of the structures was examined by very few investigators [28, 29]. They proposed that values of the normal and shear stiffnesses for rock joints typically can range roughly from 10 to 100 MPa/m for joints with soft clay infilling and to over 100 GPa/m for tight joints in granite and basalt. So, mechanical parameters of the interface elements are considered for this range in this study. Mechanical behaviour of the normal and shear spring stiffnesses is explained in Section 2. These stiffness parameters depend on the elasticity modulus and shear modulus of the materials. Thus, the elasticity modulus and shear modulus are assumed as the variables in the numerical analyses. As a consequence, various stiffness values are obtained for each variable (Table 1). Special fish codes are used while defining normal and shear stiffness parameters to FLAC3D software. The effect of these variable parameters on the viscoplastic behaviour of the dam is examined for empty and full reservoir conditions in this section. Totally 20 different normal and shear stiffness parameters are used for both reservoir conditions in the numerical analyses as shown Table 2.Table 2
Normal and shear stiffness parameters for numerical analyses.
Case
Stiffness (kn and ks) (Pa/m)
1
0
2
101
3
102
4
103
5
104
6
105
7
5 × 105
8
106
9
5 × 106
10
107
11
5 × 107
12
108
13
5 × 108
14
109
15
5 × 109
16
1010
17
5 × 1010
18
1011
19
5 × 1011
20
1012Solution algorithm for 3D interaction analyses is shown in Figure6. According to numerical results, the maximum settlements, maximum horizontal displacements, and maximum principal stresses during lifetime of the dam are examined considering 20 different interaction situations. These results are compared graphically for empty and full reservoir conditions. After the most critical stiffness values for the interaction analyses of the rockfill dams are determined, time-dependent viscoplastic analyses of the Ilısu Dam are performed using these critical stiffness values (kn and ks), and it is observed that how the interaction behaviour of the Ilısu Dam will change in the future by the effect of normal and shear stiffnesses.Figure 6
Solution algorithm for three-dimensional numerical analyses.


### 4.1. Vertical Displacements

Monitoring of the dam’s settlement behaviour is one of the most important factors for assessing the future and safety of these water structures, and it provides a warning system for abnormal behaviour of the dams. In this section, time-dependent settlement behaviour of the Ilısu rockfill dam is monitored for empty and full reservoir conditions considering various interaction parameters (normal-shear spring stiffnesses) between the dam body and foundation (Figure7). 4 different significant points are selected from the dam surface to better evaluate the effect of the normal and shear spring stiffnesses on the behaviour of the dam. These critical points are shown in Figure 8.Figure 7
Settlements for selected 4 points on the dam body surface: (a) Point 1, (b) Point 2, (c) Point 3, and (d) Point 4.
(a)
(b)
(c)
(d)Figure 8
Selected points on the dam body surface.According to the numerical analyses, viscoplastic behaviour of the Ilısu Dam changes by the effect of changing normal-shear stiffness parameters for empty and full reservoir conditions. When the full and empty reservoir conditions are compared, more vertical displacements are observed for the full reservoir condition. For the full reservoir condition, the maximum and minimum settlements are observed on Point 3 and Point 1, respectively (Figure7). In other words, the maximum settlement occurs on the approximately middle of the dam surface. But, when the empty reservoir condition of the dam is investigated, maximum settlements occurred at the crest point of the dam. This result clearly indicates the effect of the hydrostatic pressure on the viscoplastic settlement behaviour of the Ilısu Dam. When examining Figure 7, the maximum settlement took place for the smallest stiffness value (knandks=0). In addition, it is clearly seen that when the stiffness parameters (kn and ks) are increased from 0 to 1012 Pa/m for numerical analyses, the vertical displacement values obviously diminish for both reservoir conditions. This result clearly shows that how the stiffness values alter the time-dependent behaviour of the Ilısu Dam. According to Figure 7, the settlements changed continuously until a certain stiffness value, and any significant settlement changes were not observed for larger stiffness values than this value. In other words, when normal and shear stiffness parameters are selected between 108 Pa/m and 1012 Pa/m, the settlement behaviour of the Ilısu Dam obviously does not change. However, if smaller stiffness than 108 Pa/m is chosen for interaction analyses, settlements on the dam body surface increase. This result provides a very important proposal for new dam design. When considering the average settlements of the rockfill dams during its lifetime, it can be understood that the stiffness values can be selected between 106 Pa/m and 108 Pa/m for interaction analyses of rockfill dams. This result provides very important information for modelling these dams.After the most critical stiffness parameters are determined for interaction analyses of the rockfill dams, time-dependent viscoplastic analyses are performed for next 35 years of the Ilısu Dam. In the numerical analyses,kn and ks stiffness values are taken into account as 108 Pa/m. Vertical displacement results are shown for the empty reservoir condition of the Ilısu Dam in Figure 9. According to the numerical results, maximum settlement took place at the crest point of the dam and approximately 37.3 cm settlement is observed on this point. In addition, 25 cm and 10 cm vertical displacements occurred at the middle and bottom of the dam body surface, respectively. The dam body is divided into two parts from the midpoint of the center to better examine the behaviour changes in the dam body, and it is clearly seen that the vertical displacements decreased from the crest to the bottom of the dam. Hydrostatic pressure is an important factor for evaluating of the settlement behaviour [30]. So, effect of the hydrostatic pressure and water flowing on the dam behaviour are examined for the full reservoir condition. First, all displacements and principal stresses, which are obtained from the empty reservoir condition of the dam, are set to zero in order to exclude the stresses and deformations. Then, reservoir water is modelled considering full reservoir height (130 m). According to the numerical results for the full reservoir condition, maximum settlement occurred at approximately middle of the dam body surface, and about 116 cm vertical displacement is observed at this point (Figure 10). Moreover, when the settlements of 4 points on the dam body surface are examined during 35 years, it is seen that the largest settlement is observed at Point 3 and the smallest settlement occurred at Point 1 (Figure 11). Though any settlement did not take place on the surface of the foundation for the empty reservoir condition, approximately 20 cm settlement is observed on the surface of the foundation for the full reservoir condition. When comparing the upstream and downstream sides of the dam, it is obviously seen that more settlements occurred at the upstream side due to the effect of the hydrostatic pressure on the settlement behaviour of the dam.Figure 9
Vertical displacements (m) for the empty reservoir condition of the Ilısu Dam after 35 years.Figure 10
Vertical displacements for the full reservoir condition of the Ilısu Dam after 35 years.Figure 11
Time-dependent vertical displacements for the full reservoir condition during 35 years.


### 4.2. Horizontal Displacements

Horizontal displacements for the rockfill dams generally take place due to the effect of the hydrostatic pressure. These displacements disrupt the stability and functionality of the rockfill dams in time. For this reason, observing of the horizontal displacements for Ilısu Dam is vitally important in order to obtain more information about its stability and future. Numerical results for various interaction conditions between dam body and foundation are shown graphically in this section (Figure12). These graphs are created taking into account the maximum horizontal displacements that may occur during Ilısu Dam’s life. Ilısu Dam is examined for two reservoir conditions to better seen the effect of the hydrostatic pressure. When comparing the empty and full reservoir conditions of the dam, more horizontal displacements are observed for the full reservoir condition as seen in Figure 12. According to the numerical analyses, the interface elements that are defined between the dam body and the foundation clearly altered the horizontal displacement behaviour of the Ilısu Dam. When the horizontal displacements for 4 different points on the dam body surface (Figure 8) are investigated, it is obviously seen that maximum and minimum horizontal displacements took place on Point 2 and Point 1, respectively (Figures 12(a) and 12(b)). In addition, less displacements are observed on Point 4 (crest point of the dam) compared with Point 2 and Point 3 (Figures 12(a) and 12(c)). Moreover, when 20 different interaction conditions are compared with each other, large horizontal displacement changes occurred between Case 1 and Case 12. However, very small displacement changes are obtained between Case 12 and Case 20. These results clearly mean that no matter how much the stiffness values are changed between 108 Pa/m and 1012 Pa/m, large horizontal displacement changes are not observed on the Ilısu Dam body for this stiffness range. But, it is seen from Figure 12 that the creep behaviour of the Ilısu Dam obviously changes for smaller stiffness values than 108 Pa/m. When considering these numerical results, it is understood that the stiffness parameters have a great effect on the horizontal displacement behaviour of the rockfill dams. This information will be a guide for designing and modelling the rockfill dam.Figure 12
Horizontal displacements for selected 4 points on the dam body surface: (a) Point 1, (b) Point 2, (c) Point 3, and (d) Point 4.
(a)
(b)
(c)
(d)After determining the most critical stiffness values for horizontal displacement behaviour of the rockfill dams, these critical stiffness parameters (kn and ks) are defined to FLAC3D software using fish codes. In the numerical analyses, shear (ks) and normal (kn) stiffness parameters are considered as 108 Pa/m, and interaction analyses are performed for empty and full reservoir conditions of the dam. The horizontal displacement behaviour of the Ilısu Dam is investigated for next 35 years of the dam, and numerical results are shown in Figures 13–15. When these results are examined, large displacement changes are observed between the empty and full reservoir conditions. For the empty reservoir condition, very small horizontal displacements occurred on the dam body surface because external loads did not expose to the dam as seen in Figure 13. However, viscoplastic behaviour of the Ilısu Dam clearly changed by effect of the hydrostatic pressure. More displacements are obtained for the full reservoir condition compared with the empty condition, and 83 cm maximum horizontal displacement is observed at the middle of the dam body surface as seen in Figure 14. When 4 critical points (Figure 8) on the dam body surface are examined, maximum and minimum displacement values took place at Point 2 and Point 1, respectively (Figure 15). Moreover, it is understood that large horizontal displacement changes will take place in the dam body surface from the year 2017 to 2040. But, displacement changes will decrease after the year 2040. This result provides very important information about future of the Ilısu Dam.Figure 13
Horizontal displacements for the empty reservoir condition of the Ilısu Dam after 35 years.Figure 14
Horizontal displacements for the full reservoir condition of the Ilısu Dam after 35 years.Figure 15
Time-dependent horizontal displacements for the full reservoir condition during 35 years.


### 4.3. Principal Stresses

In this section, the effect of the interaction parameters on the principal stress behaviour of the Ilısu Dam is examined in detail. The numerical results are presented graphically in Figure16. Graphs are created considering the maximum principal stress values that may take place during dam’s lifetime. Generally, it is clearly understood from the numerical graphs that the principal stress values for the full reservoir condition are larger than those for the empty reservoir condition. This result obviously shows the effect of the hydrostatic pressure on the time-dependent viscoplastic behaviour of the Ilısu Dam. When examining 4 different points that are selected on the dam surface (Figure 8), the maximum and minimum principal stresses are observed at Point 2 and Point 4, respectively (Figures 16(b) and 16(d)). This important conclusion means that when considering the hydrostatic pressure effect, the maximum principal stresses take place at approximately middle of the dam body surface. When comparing 20 different interaction conditions, large stress changes are observed between Case 1 and Case 12 (Figure 16). In other words, if the stiffness parameters are selected between 0 and 108 Pa/m for the interaction analyses of the rockfill dams, great principal stress changes occur on the dam body surface for this stiffness range. However, when these parameters are chosen as greater value than 108 Pa/m, no large stress changes are observed on the dam body surface (Figure 16). According to these results, it is understood that the most critical normal and shear stiffness parameters are 108 Pa/m for principal stress analyses of the rockfill dams. This result provides very important support for modelling and analysing these dams.Figure 16
Principal stresses for selected 4 points on the dam body surface: (a) Point 1, (b) Point 2, (c) Point 3, and (d) Point 4.
(a)
(b)
(c)
(d)After the most critical stiffness range is determined for the principal stress behaviour of the Ilısu Dam, creep analyses are performed considering this critical range. Normal and shear stiffness values are selected as 108 Pa/m in the numerical analyses. Time-dependent results are presented in Figures 17–19. Ilısu Dam is examined for 2 different reservoir conditions to better observe the effect of the hydrostatic pressure on the viscoplastic behaviour of the dam. For the empty reservoir condition, any principal stresses are not observed on the dam body surface, and the stresses increased from the crest to the foundation surface (Figure 17). The 3D dam model is divided into two halves to better see principal stress changes in the dam body. When investigating the split half condition of the Ilısu Dam, approximately 4.61 MPa maximum principal stress is observed at the bottom of the 3D model. In addition, 2.5 MPa principal stress took place at the bottom of the dam body as seen in Figure 17. As soon as the reservoir water interacted with the dam body surface, principal stress behaviour of the Ilısu Dam obviously changed. More principal stresses occurred for the full reservoir condition compared with the empty reservoir condition (Figure 18). For the full reservoir condition, maximum principal stress is 6.78 MPa, and it is observed at the bottom of the foundation. Approximately 1.4 MPa principal stress is obtained on the surface of the foundation, and 2 MPa principal stresses occurred at the middle of the dam body. When examining 4 different points on the dam body surface, the maximum and minimum stresses that may occur on the dam body surface during 35 years are observed at Point 2 and Point 4, respectively (Figure 19). This result indicates the effect of the hydrostatic pressure on the Ilısu Dam principal stress behaviour.Figure 17
Principal stresses for the empty reservoir condition of the Ilısu Dam for next 35 years.Figure 18
Principal stresses for the full reservoir condition of the Ilısu Dam for next 35 years.Figure 19
Time-dependent principal stresses for the full reservoir condition during 35 years.


## 4.1. Vertical Displacements

Monitoring of the dam’s settlement behaviour is one of the most important factors for assessing the future and safety of these water structures, and it provides a warning system for abnormal behaviour of the dams. In this section, time-dependent settlement behaviour of the Ilısu rockfill dam is monitored for empty and full reservoir conditions considering various interaction parameters (normal-shear spring stiffnesses) between the dam body and foundation (Figure7). 4 different significant points are selected from the dam surface to better evaluate the effect of the normal and shear spring stiffnesses on the behaviour of the dam. These critical points are shown in Figure 8.Figure 7
Settlements for selected 4 points on the dam body surface: (a) Point 1, (b) Point 2, (c) Point 3, and (d) Point 4.
(a)
(b)
(c)
(d)Figure 8
Selected points on the dam body surface.According to the numerical analyses, viscoplastic behaviour of the Ilısu Dam changes by the effect of changing normal-shear stiffness parameters for empty and full reservoir conditions. When the full and empty reservoir conditions are compared, more vertical displacements are observed for the full reservoir condition. For the full reservoir condition, the maximum and minimum settlements are observed on Point 3 and Point 1, respectively (Figure7). In other words, the maximum settlement occurs on the approximately middle of the dam surface. But, when the empty reservoir condition of the dam is investigated, maximum settlements occurred at the crest point of the dam. This result clearly indicates the effect of the hydrostatic pressure on the viscoplastic settlement behaviour of the Ilısu Dam. When examining Figure 7, the maximum settlement took place for the smallest stiffness value (knandks=0). In addition, it is clearly seen that when the stiffness parameters (kn and ks) are increased from 0 to 1012 Pa/m for numerical analyses, the vertical displacement values obviously diminish for both reservoir conditions. This result clearly shows that how the stiffness values alter the time-dependent behaviour of the Ilısu Dam. According to Figure 7, the settlements changed continuously until a certain stiffness value, and any significant settlement changes were not observed for larger stiffness values than this value. In other words, when normal and shear stiffness parameters are selected between 108 Pa/m and 1012 Pa/m, the settlement behaviour of the Ilısu Dam obviously does not change. However, if smaller stiffness than 108 Pa/m is chosen for interaction analyses, settlements on the dam body surface increase. This result provides a very important proposal for new dam design. When considering the average settlements of the rockfill dams during its lifetime, it can be understood that the stiffness values can be selected between 106 Pa/m and 108 Pa/m for interaction analyses of rockfill dams. This result provides very important information for modelling these dams.After the most critical stiffness parameters are determined for interaction analyses of the rockfill dams, time-dependent viscoplastic analyses are performed for next 35 years of the Ilısu Dam. In the numerical analyses,kn and ks stiffness values are taken into account as 108 Pa/m. Vertical displacement results are shown for the empty reservoir condition of the Ilısu Dam in Figure 9. According to the numerical results, maximum settlement took place at the crest point of the dam and approximately 37.3 cm settlement is observed on this point. In addition, 25 cm and 10 cm vertical displacements occurred at the middle and bottom of the dam body surface, respectively. The dam body is divided into two parts from the midpoint of the center to better examine the behaviour changes in the dam body, and it is clearly seen that the vertical displacements decreased from the crest to the bottom of the dam. Hydrostatic pressure is an important factor for evaluating of the settlement behaviour [30]. So, effect of the hydrostatic pressure and water flowing on the dam behaviour are examined for the full reservoir condition. First, all displacements and principal stresses, which are obtained from the empty reservoir condition of the dam, are set to zero in order to exclude the stresses and deformations. Then, reservoir water is modelled considering full reservoir height (130 m). According to the numerical results for the full reservoir condition, maximum settlement occurred at approximately middle of the dam body surface, and about 116 cm vertical displacement is observed at this point (Figure 10). Moreover, when the settlements of 4 points on the dam body surface are examined during 35 years, it is seen that the largest settlement is observed at Point 3 and the smallest settlement occurred at Point 1 (Figure 11). Though any settlement did not take place on the surface of the foundation for the empty reservoir condition, approximately 20 cm settlement is observed on the surface of the foundation for the full reservoir condition. When comparing the upstream and downstream sides of the dam, it is obviously seen that more settlements occurred at the upstream side due to the effect of the hydrostatic pressure on the settlement behaviour of the dam.Figure 9
Vertical displacements (m) for the empty reservoir condition of the Ilısu Dam after 35 years.Figure 10
Vertical displacements for the full reservoir condition of the Ilısu Dam after 35 years.Figure 11
Time-dependent vertical displacements for the full reservoir condition during 35 years.


## 4.2. Horizontal Displacements

Horizontal displacements for the rockfill dams generally take place due to the effect of the hydrostatic pressure. These displacements disrupt the stability and functionality of the rockfill dams in time. For this reason, observing of the horizontal displacements for Ilısu Dam is vitally important in order to obtain more information about its stability and future. Numerical results for various interaction conditions between dam body and foundation are shown graphically in this section (Figure12). These graphs are created taking into account the maximum horizontal displacements that may occur during Ilısu Dam’s life. Ilısu Dam is examined for two reservoir conditions to better seen the effect of the hydrostatic pressure. When comparing the empty and full reservoir conditions of the dam, more horizontal displacements are observed for the full reservoir condition as seen in Figure 12. According to the numerical analyses, the interface elements that are defined between the dam body and the foundation clearly altered the horizontal displacement behaviour of the Ilısu Dam. When the horizontal displacements for 4 different points on the dam body surface (Figure 8) are investigated, it is obviously seen that maximum and minimum horizontal displacements took place on Point 2 and Point 1, respectively (Figures 12(a) and 12(b)). In addition, less displacements are observed on Point 4 (crest point of the dam) compared with Point 2 and Point 3 (Figures 12(a) and 12(c)). Moreover, when 20 different interaction conditions are compared with each other, large horizontal displacement changes occurred between Case 1 and Case 12. However, very small displacement changes are obtained between Case 12 and Case 20. These results clearly mean that no matter how much the stiffness values are changed between 108 Pa/m and 1012 Pa/m, large horizontal displacement changes are not observed on the Ilısu Dam body for this stiffness range. But, it is seen from Figure 12 that the creep behaviour of the Ilısu Dam obviously changes for smaller stiffness values than 108 Pa/m. When considering these numerical results, it is understood that the stiffness parameters have a great effect on the horizontal displacement behaviour of the rockfill dams. This information will be a guide for designing and modelling the rockfill dam.Figure 12
Horizontal displacements for selected 4 points on the dam body surface: (a) Point 1, (b) Point 2, (c) Point 3, and (d) Point 4.
(a)
(b)
(c)
(d)After determining the most critical stiffness values for horizontal displacement behaviour of the rockfill dams, these critical stiffness parameters (kn and ks) are defined to FLAC3D software using fish codes. In the numerical analyses, shear (ks) and normal (kn) stiffness parameters are considered as 108 Pa/m, and interaction analyses are performed for empty and full reservoir conditions of the dam. The horizontal displacement behaviour of the Ilısu Dam is investigated for next 35 years of the dam, and numerical results are shown in Figures 13–15. When these results are examined, large displacement changes are observed between the empty and full reservoir conditions. For the empty reservoir condition, very small horizontal displacements occurred on the dam body surface because external loads did not expose to the dam as seen in Figure 13. However, viscoplastic behaviour of the Ilısu Dam clearly changed by effect of the hydrostatic pressure. More displacements are obtained for the full reservoir condition compared with the empty condition, and 83 cm maximum horizontal displacement is observed at the middle of the dam body surface as seen in Figure 14. When 4 critical points (Figure 8) on the dam body surface are examined, maximum and minimum displacement values took place at Point 2 and Point 1, respectively (Figure 15). Moreover, it is understood that large horizontal displacement changes will take place in the dam body surface from the year 2017 to 2040. But, displacement changes will decrease after the year 2040. This result provides very important information about future of the Ilısu Dam.Figure 13
Horizontal displacements for the empty reservoir condition of the Ilısu Dam after 35 years.Figure 14
Horizontal displacements for the full reservoir condition of the Ilısu Dam after 35 years.Figure 15
Time-dependent horizontal displacements for the full reservoir condition during 35 years.


## 4.3. Principal Stresses

In this section, the effect of the interaction parameters on the principal stress behaviour of the Ilısu Dam is examined in detail. The numerical results are presented graphically in Figure16. Graphs are created considering the maximum principal stress values that may take place during dam’s lifetime. Generally, it is clearly understood from the numerical graphs that the principal stress values for the full reservoir condition are larger than those for the empty reservoir condition. This result obviously shows the effect of the hydrostatic pressure on the time-dependent viscoplastic behaviour of the Ilısu Dam. When examining 4 different points that are selected on the dam surface (Figure 8), the maximum and minimum principal stresses are observed at Point 2 and Point 4, respectively (Figures 16(b) and 16(d)). This important conclusion means that when considering the hydrostatic pressure effect, the maximum principal stresses take place at approximately middle of the dam body surface. When comparing 20 different interaction conditions, large stress changes are observed between Case 1 and Case 12 (Figure 16). In other words, if the stiffness parameters are selected between 0 and 108 Pa/m for the interaction analyses of the rockfill dams, great principal stress changes occur on the dam body surface for this stiffness range. However, when these parameters are chosen as greater value than 108 Pa/m, no large stress changes are observed on the dam body surface (Figure 16). According to these results, it is understood that the most critical normal and shear stiffness parameters are 108 Pa/m for principal stress analyses of the rockfill dams. This result provides very important support for modelling and analysing these dams.Figure 16
Principal stresses for selected 4 points on the dam body surface: (a) Point 1, (b) Point 2, (c) Point 3, and (d) Point 4.
(a)
(b)
(c)
(d)After the most critical stiffness range is determined for the principal stress behaviour of the Ilısu Dam, creep analyses are performed considering this critical range. Normal and shear stiffness values are selected as 108 Pa/m in the numerical analyses. Time-dependent results are presented in Figures 17–19. Ilısu Dam is examined for 2 different reservoir conditions to better observe the effect of the hydrostatic pressure on the viscoplastic behaviour of the dam. For the empty reservoir condition, any principal stresses are not observed on the dam body surface, and the stresses increased from the crest to the foundation surface (Figure 17). The 3D dam model is divided into two halves to better see principal stress changes in the dam body. When investigating the split half condition of the Ilısu Dam, approximately 4.61 MPa maximum principal stress is observed at the bottom of the 3D model. In addition, 2.5 MPa principal stress took place at the bottom of the dam body as seen in Figure 17. As soon as the reservoir water interacted with the dam body surface, principal stress behaviour of the Ilısu Dam obviously changed. More principal stresses occurred for the full reservoir condition compared with the empty reservoir condition (Figure 18). For the full reservoir condition, maximum principal stress is 6.78 MPa, and it is observed at the bottom of the foundation. Approximately 1.4 MPa principal stress is obtained on the surface of the foundation, and 2 MPa principal stresses occurred at the middle of the dam body. When examining 4 different points on the dam body surface, the maximum and minimum stresses that may occur on the dam body surface during 35 years are observed at Point 2 and Point 4, respectively (Figure 19). This result indicates the effect of the hydrostatic pressure on the Ilısu Dam principal stress behaviour.Figure 17
Principal stresses for the empty reservoir condition of the Ilısu Dam for next 35 years.Figure 18
Principal stresses for the full reservoir condition of the Ilısu Dam for next 35 years.Figure 19
Time-dependent principal stresses for the full reservoir condition during 35 years.


## 5. Conclusion

In this paper, the effect of the dam body and foundation interaction on the time-dependent viscoplastic behaviour of the Ilısu Dam is examined in detail. The three-dimensional finite difference model of the Ilısu Dam is modelled using special fish codes, and it is created according to the original dam project. The special material models are used for rockfill and foundation materials in the creep analyses. These material models were rarely used for creep analyses of the rockfill dams, previously. Thus, this study is very important in terms of evaluating the effect of the different material models on the viscoplastic behaviour of the CFR dams. 20 different interaction parameters (normal and shear stiffnesses) are used between the dam body and foundation for the interaction analyses of Ilısu CFR dam. Therefore, totally 20 various interaction analyses are performed for the empty and full reservoir conditions of the dam. The effect of these interaction conditions on the viscoplastic behaviour of the Ilısu Dam is evaluated as below:(i)
According to the numerical results, it is clearly understood that the hydrostatic water pressure altered the viscoplastic interaction behaviour of the Ilısu Dam. When comparing the empty and full reservoir conditions, more settlements, horizontal displacements, and principal stresses are observed for the full reservoir condition.(ii)
When examining the vertical displacements on the dam body surface that may occur from the year 2017 to 2052, 0.37 m maximum settlement is observed at the crest of the dam body for the empty reservoir condition. However, as soon as the reservoir water contacts the dam body surface, the settlement behaviour of the Ilısu Dam clearly changes. 1.16 m maximum settlement value is obtained at the approximately middle of the dam body surface for the full reservoir condition. This result obviously indicates the effect of the reservoir water pressure on the creep behaviour of the Ilısu Dam.(iii)
For the empty reservoir condition, very small horizontal displacements are observed on the dam body surface because hydrostatic pressure did not contact the dam body surface. But, horizontal displacements obviously increased on the dam body surface by the effect of the hydrostatic water pressure. 0.83 m maximum horizontal displacement is observed at the approximately middle of the dam body surface for the full reservoir condition.(iv)
According to the stress analyses, any principal stress value is not observed at the dam body surface for the empty reservoir condition of the dam. But, approximately 2 MPa principal stress value is acquired on the dam body for the full reservoir condition.(v)
According to the interaction analyses, it is clearly seen that the most critical shear and normal stiffness values are 108 Pa/m for creep analyses of the CFR dams. When the stiffness parameters (normal and shear stiffnesses) that are defined between dam body and foundation is selected between 0 and 108 Pa/m in the interaction analyses, large changes are observed in the time dependent viscoplastic behaviour of the CFR dams. Moreover, if these parameters are chosen as larger than 108 Pa/m in the creep analyses, small changes are observed in the viscoplastic behaviour of these dams. When considering the average settlement values of the rockfill dams during its lifetime, it can be understood that the normal and shear stiffness values can be selected between 106 Pa/m and 108 Pa/m for interaction analyses of rockfill dams. This conclusion is very important for modelling and analysing the CFR dams.(vi)
It is clearly seen from this study, normal and shear spring stiffness parameters are very important for evaluating the creep behaviour of the rockfill dams. So, these parameters should not be ignored in the time-dependent interaction analyses of CFR dams.



---

*Source: 2491652-2018-09-02.xml*