Investigation of seismic performance of concrete gravity dams using probabilistic analysis

The effect of the modulus of elasticity of concrete on seismic behavior of Koyna gravity dam in India is studied in the paper using probabilistic analysis. Numerical model based on the finite element method is used to analyse the base-case scenario involving the dam-reservoir-foundation interaction. The results show that the modulus of elasticity significantly affects seismic behaviour of concrete gravity dams. The results of the analysis are presented as bilinear curves.


Introduction
Gravity dams can be found all over the world and some of them have been built in earthquake-prone areas.The need to build gravity dams in high seismicity regions is expressed even today due to an increasing demand for both water supply and flood protection.On the other hand, the accuracy of risk evaluation associated with existing dams, as well as the success of design of future dams, are highly dependent on proper understanding of their behaviour with regard to earthquake action.Because the dam failure consequences can be disastrous, the seismic design of such structures has been widely recognized as being of particular significance.As a result, further study of seismic behaviour of dams remains a topical issue for engineers [1][2][3].Some of the important aspects that may affect the response of gravity dams subjected to earthquakes have been recognized through the dam-reservoir-foundation interaction.The damwater interaction must be taken into account since the dam undergoes deformation, which influences the motion of water in the reservoir.The interaction between the dam and water stored behind the dam leads to an increase in the period of vibrations at the dam.The reason for this is that the dam can not move without spatial variability of water tangent.The water that moves along the dam increases the total mass moved due to earthquake.The added mass increases the natural vibrations period of the dam, and affects inertial forces created due to earthquake action [3][4].The interface of the dam with the impounded water is an important boundary where the static and hydrodynamic forces are applied to the dam structure.These forces provide a significant contribution to seismic response analysis and design of dams.The maximum hydrodynamic pressure on the dam when subjected to a moderately strong earthquake ground motion may reach the magnitude of hydrostatic pressure [4][5][6].The seismic performance of dams with regard to interaction effects is dependent on some uncertain parameters.Specific concrete properties used in the design of concrete gravity dams include unit weight, compressive, tensile, and shear strengths, modulus of elasticity, creep, Poisson's ratio, coefficient of thermal expansion, thermal conductivity, specific heat, and diffusivity.Out of these properties, the modulus of elasticity is considered to be the most important for the design of dams.In fact, the modulus of elasticity is the main factor affecting concrete strength.Concrete strength should satisfy the early load and construction requirements, and the stress criteria.Thus, selecting an optimum value of the modulus of elasticity is very important in the design of concrete gravity dams.The probabilistic and sensitivity analysis is used in the paper to identify the modulus of elasticity of concrete as a particular parameter that has a significant effect on seismic behaviour of gravity dams.The Monte Carlo simulation with Latin hypercube sampling (LHS) is applied as the probabilistic analysis technique for evaluating intensity of parameter influence on seismic behaviour of concrete gravity dams.Monte Carlo simulations are typically characterized by a large number of unknown parameters, many of which are difficult to obtain experimentally.Some of the parameters that are used in seismic analysis of concrete dams are the concrete modulus of elasticity, concrete density, and modulus of elasticity of the bed rock or foundation.The parameter sensitivity analysis can be used to quantify the effect of unknown parameters.Sensitivity analysis, as applied in risk assessment, is dependent on the variability and uncertainty of factors contributing to risk.In short, sensitivity analysis identifies what is "driving" the risk estimates.It is used in both point estimate and probabilistic approaches to identify and rank important sources of variability, as well as important sources of uncertainty.The quantitative information obtained by sensitivity analysis is important for properly dealing with complexity of the analysis, and for communicating important results.The sensitivity analysis focuses on a set of graphical and statistical techniques that can be used to determine which variables present in risk model contribute most to the variation of risk estimates.This variation in risk could represent variability or uncertainty, depending on the type of risk model and characterization of input variables [7].The Monte Carlo Simulation method is the most common method for probabilistic analysis.One simulation loop represents one system component that is subjected to a particular set of loads and boundary conditions.The method is always applicable regardless of the physical effect modelled in a finite element analysis.Assuming the deterministic model is correct and a very large number of simulation loops are performed, Monte Carlo techniques always provide correct probabilistic results.Monte Carlo simulations can be employed using either the Direct Sampling method or the Latin Hypercube Sampling method.Monte Carlo simulation can be improved using the LHS as described by McKay et al. [8].The Latin Hypercube Sampling may be viewed as a stratified sampling scheme designed to ensure that the upper or lower ends of the distributions used in the analysis are well represented.Latin hypercube sampling is generally recommended over Direct Sampling method when the model is complex or when time and resource constraints are an issue [8].The nature of LHS does not determine the sample size needed to achieve a certain confidence level.There is no specified value for sample size N to achieve a certain confidence level in LHS [9].By sampling N times from the parameter distributions, this procedure creates a population of N possible instances of the structure, each of which needs to be analysed.The use of relatively high N that is substantially larger than the number of parameters will always result in reasonably accurate estimates for practical purposes.The optimal N to use is a function of the number of random variables and their influence on the response [10].In this research, Monte Carlo with LHS is performed and the value of N equal to 30 is chosen to allow sufficient accuracy in the estimates, and the sensitivity of the dam performance to the main parameter is investigated.For this purpose, the Koyna dam in India is considered as a case study and a complete analysis is conducted using the finite element method in the Investigation of seismic performance of concrete gravity dams using probabilistic analysis time domain.The effects of concrete modulus of elasticity on seismic performance of gravity dam are also discussed using probabilistic and sensitivity analysis.

Governing equations
The structural and hydrodynamic aspects of the problem involving the base-case model are formulated separately in this section.The dam and foundation are considered as an elastic solid with linear and plane stress behaviour.Water is taken as a compressible, inviscid fluid, and the dam as an elastic solid [11][12][13].It must be mentioned that the dam and foundation are different in material characteristics such as the modulus of elasticity, density, and Poisson factor.So, the dam and foundation are the structural section of system, and they are modelled together but with relative material characteristics.

Dam and foundation model
The governing equation for the dam and foundation model is the motion equation.But, in order to completely describe the fluidstructure interaction problem, the fluid pressure load acting at the interface is now added to structural equation.The reason is that the dam and foundation move with spatial variability of water tangent.The moved water increases the total mass and inertia due to an earthquake action [11][12][13].So, the equations of the system subjected to the ground motion including the effects of reservoir are written as (1) where: M -the mass matrix C -the structural damping matrix K -the structural stiffness matrix.K can be obtained from material and strain-displacement matrix, and C is proportional to the mass and stiffness matrix according to the Rayleigh method [14][15].u is the vector of displacement relative to ground, is the vector of velocity, and is the vector of ground acceleration.The fluid pressure load F Pr at the dam-reservoir and foundation-reservoir interface is induced because of interaction, and it is obtained by integrating the hydrodynamic pressure over the dam and foundation wetting of the unit length surface.The hydrodynamic pressure is induced at the solid-fluid interface because of the following boundary condition [16][17][18]: (2) where: n -denotes the inward normal direction to interface a n -the normal component of acceleration ρ w -he mass density of water.

Reservoir model
In acoustical fluid-structure interaction problems, the structural dynamics equation needs to be considered along with the Euler equations of fluid momentum and the flow continuity equation.Assuming that the water in the reservoir is inviscid, compressible and irrotational, and that its motion is of small amplitude, the fluid momentum and continuity equations are simplified to get the acoustic wave equation as follows [16][17][18]: where: C -the speed of sound in fluid medium P -the hydrodynamic pressure.
Since the viscous dissipation has been neglected, equation ( 3) is referred to as the frictionless wave equation for propagation of sound in fluids.The discretized structural equation and the frictionless wave equation have to be considered simultaneously in the fluid-structure interaction problems.The fluid pressure acting on the structure at the fluid-structure interface will be considered to form the coupling stiffness matrix.The Sommerfeld boundary condition is implemented for a truncated boundary at the tail of reservoir.This boundary condition is assuming a damper in the end of reservoir with the following equation [17][18]: (4) where: C -the sound velocity in water n -denotes the outward normal direction to the far end of reservoir.

Finite element formulation
Governing equations of the fluid-structure system can be expressed in the matrix form using the finite element method [14].In order to completely describe the fluid-structure interaction problem, the fluid pressure load acting at the interface is now added to structural equation.The finite element approximating shape functions for the spatial variation of displacement components can be expressed as [15]: where: {N} i {N'} -the element shape functions for pressure and displacements {P e } -the nodal pressure vector {u e } -the nodal displacement component vector.

Majid Pasbani Khiavi
From equation ( 5) and equation ( 6), the first and second time derivative of the variables and the virtual change in the pressure can be written as follows: δP = {N} T {δP e } (10)

Finite element model of dam and foundation
The discretized structural dynamics equation can be formulated by means of structural elements.The structural equation is rewritten here as follows [15]: where: -

Finite element model of reservoir
The matrix operator applied to the fluid element shape function {N} is defined by [15]: So, the finite element statement of the wave equation is given by: (17) In which {n} is the normal vector at dam-reservoir and foundation-reservoir interface boundary.Terms which do not vary over the element are taken out of the integration sign.{δP} is an arbitrarily introduced virtual change in nodal pressure, and it can be factored out in equation ( 17).Since {δP} is not equal to zero, equation ( 17) becomes: Equation ( 18) can be written in matrix notation to get the discretized wave equation: (19) where: ρ 0 -Water density -Fluid mass matrix -Fluid stiffness matrix -Mass matrix related to fluidstructure interaction The finite element discretization and numerical time integration procedures developed in previous section have been implemented into a finite element model.In the model, the hydrodynamic pressure, solid displacement, and stresses are nodal unknown variables.In this paper, the Newmark method is applied to solve discretized dynamic equations.The step by step solution based on the Newmark integration method is defined as follows: -Form stiffness matrix K, mass matrix M, and damping matrix C -Select time step ∆t and parameters α and δ and calculate integration constants [15]: -Form effective stiffness matrix: -Triangularize: -For each time step calculate effective loads at time t + ∆t: Investigation of seismic performance of concrete gravity dams using probabilistic analysis -Solve for displacement and pressure at time t + ∆t: or -Calculate acceleration and velocity at time t + ∆t: where: ∆t -the time step α i δ -parameters that can be determined to obtain integration accuracy and stability.
When α= 0.25 and δ = 0.5 the Newmark method is unconditionally stable [15].So, these values are selected in the paper as integration constants.

Case study
To demonstrate effectiveness of the analytic procedure presented in this paper and the effect of the modulus of elasticity on seismic performance of concrete gravity dams, the response of Koyna Dam to the horizontal and vertical component of El Centro earthquake is presented.
The Koyna Dam is one of the largest dams built on the Koyna river in Maharashtra State in western India.The dam has withstood many earthquakes in recent past.This large engineering structure is located in one of the earthquakeprone zones of Maharashtra and is, therefore, an ideal site for deformation monitoring studies.The structure is subjected to the effect of crustal movements, seasonal changes in water load in the reservoir, and to the effect of seismicity in the region.Consequently, the Koyna Dam has been selected as case study for the seismic and risk analysis by many researchers.
The 1940 El Centro earthquake (or the 1940 Imperial Valley earthquake) occurred at Pacific Standard Time on May 18 in the Imperial Valley in Southern California near the international border of the United States and Mexico.It was the first major earthquake to be recorded by a strong-motion seismograph located next to a fault rupture.The earthquake was characterized as a typical moderate-sized destructive event with a complex energy release signature.It was the strongest recorded earthquake to hit the Imperial Valley, and it caused widespread damage to irrigation systems and resulted in nine fatalities.The El Centro earthquake is usually used to evaluate seismic performance and safety of structures during earthquake action.Figure 1 and Figure 2 show ten seconds of the horizontal and vertical components of the El Centro site records, as selected for the purpose of seismic analysis.The values of integration parameters according to the Newmark method were taken as α= 0.25 and δ = 0.5 with the time step (∆t) equal to 0.02 second.The modulus of elasticity, density and Poisson`s ratio of concrete were taken as 20.7 GPa, 2643 kg/m 3 and 0.2 for the dam structure, while for foundation, the modulus of elasticity and Poisson s ratio were taken as 16.86 GPa and 0.2.For Majid Pasbani Khiavi reservoir, the water density and bulk modulus of water were taken as 1000 kg/m 3 and 2.1 GPa [16].
The stiffness and mass proportional damping (Rayleigh damping) is used in the analysis.The velocity of pressure wave in water was taken as 1438.66m/s.The lognormal distribution, which is a basic and commonly used distribution, was used to describe scatter of the input data.The lognormal distribution is very suitable for phenomena that arise from the multiplication of a large number of error effects.
It is also correct to use the lognormal distribution for a random variable resulting from multiplication of two or more random effects [7,19,20].
The simulations performed here comprised 30 simulations per random seed number.These simulations consisted of a basecase simulation with the modulus of elasticity (E C ) amounting to 20.7 GPa.The maximum of horizontal displacement at dam crest, 1 st principal stress at heel, and 3 rd principal stress at dam toe, were assumed to be critical responses during earthquake and were selected as output parameters affected by E C .Their sensitivity was investigated.

Model analysis
First, the base-case model was analysed.Then the probabilistic analysis was done to show the effect of variation of the modulus of elasticity as related to seismic performance of the dam-reservoir-foundation system.Figures 5 to 8 show the results obtained from time history analysis of the base-case model with the mean value of modulus of elasticity.In Figure 5, positive and negative values of horizontal displacement denote dam crest movement in the downstream and upstream directions along the river, respectively.Investigation of seismic performance of concrete gravity dams using probabilistic analysis

Probabilistic analysis
In this section, the probabilistic analysis is illustrated using cumulative distribution functions of output parameters.The cumulative distribution functions are alternatively referred to in literature as the distribution function, cumulative frequency function, or the cumulative probability function.The cumulative distribution function expresses the probability that a random variable will assume a value lower than or equal to some value.For continuous random variables, the cumulative distribution function is obtained from the probability density function by integration, or by summation in the case of discrete random variables   The value of the cumulative distribution function is the probability that the values will stay below a specific value.For example, it can result from curves: -There is about 50 % probability that the horizontal maximum displacement of dam crest will remain below 7.5 cm during earthquake.
-There is about 50 % probability that the maximum hydrodynamic pressure at heel of dam will remain below 340 KPa during earthquake.
-There is about 50 % probability that the 1 st principal stress at the heel of dam will remain below 4.5 MPa.
-There is about 50 % probability that the 3rd principal stress at the toe of dam will remain below 3.6 MPa.
If cumulative distribution is considered, the comparison between stresses induced in dam body with allowable stress of concrete can be possible.For example, if the allowable tensile stress is equal to 3.56 MPA, it can be concluded that there is about 15 % probability that the tensile principal stress at the heel of dam will remain below the allowable value.

Conclusion
The modulus of elasticity is of highest significance in the design of dams because it varies with concrete strength.While concrete strengths should satisfy the early load and construction requirements and the stress criteria for safety design, it is very important in the design of concrete gravity dams to select an optimum value of the modulus of elasticity.In this paper, a probabilistic and sensitivity analysis of a dam-reservoir-foundation system was performed using Monte Carlo simulation to show the effect of the modulus of elasticity on seismic performance of concrete gravity dams.The method was used to identify the modulus of elasticity (E C ) as a particular parameter that has a significant effect on the responses.Maximum values of horizontal displacement of dam crest, hydrodynamic pressure at heel, 1 st principal stress at heel, and 3 rd principal stress at dam toe, were assumed as critical responses during earthquake, and selected as output parameters that are strongly dependant on the variation of E C .Because of the relation between the modulus of elasticity and concrete strength, conclusions can be made about uncertainty in the design of gravity dams using cumulative distribution function of results.For example, it is possible to evaluate the stresses induced in dam body with the allowable stress of concrete for safety design.Finally, it must be mentioned that the model is applicable in probabilistic analysis of other parameters as a means to show sensitivity of responses.Its use is particularly recommended for the realistic analysis of large dams.
Solid mass matrix -Solid stiffness matrix -Strain-displacement matrix -Matrix operator applied to element shape functions {F e } -Force of seismic loading -Matrix related to fluidstructure interaction ρ -Concrete density [C] -Material matrix of plane stress element.

Figure 7 .Figure 8 .
Figure 7. Time history of 1 st principal stress at dam heel

Figure 9 .Figure 10 .Figure 11 .
Figure 9. Sensitivity of horizontal displacement curve to concrete modulus of elasticity [7].The cumulative distribution function also visualizes what the reliability or failure probability would be if it were necessary to change admissible limits of the design.Figures 12 to 14 show the cumulative distribution function of selected responses of the model.

Figure 12 .
Figure 12.Cumulative distribution function of horizontal displacement of dam crest

Figure 13 .Figure 14 .
Figure 13.Cumulative distribution function of hydrodynamic pressure at heel of dam