Reservoir Fluctuation Effects on Seismic Response of High Concrete Arch Dams Considering Material Nonlinearity

Reservoir fluctuation effects on seismic response of concrete arch dams are investigated. Structure nonlinearity is originated from material nonlinearity due to tensile cracking and compression crushing of mass concrete using William-Warnke failure surface in principal stress state and reservoir is assumed compressible. An arch dam is selected as case study and TABAS earthquake record is used to excite the finite element model of dam-reservoir-foundation. It is found that principal stresses on upstream and downstream faces increase significantly with reservoir dewatering. Responses of dam in nonlinear model have special intricacies, but the extension of cracked areas develops meaningfully with decreasing reservoir water.


Introduction
The real seismic response of a high arch dam to earthquake loadings is a very complicated problem and depends on various factors such as interaction of the dam with foundation rock and reservoir water, numerical modelling, selecting appropriate material properties, suitable boundary conditions for system and finally extracting eligible earthquake record for analysing dam-reservoir-foundation system. Therefore, an efficient method is required to properly assess the safety of an arch dam located in the region with significant seismicity. From this point of view, nonlinear seismic analysis of arch dams for earthquake ground motions should be considered based on a detailed analytical model. In addition, this model should be capable of representing nonlinear behaviour of mass concrete in dam body. Interaction between the dam and impounded water is one of key factors affecting the seismic response of arch dams during earthquake ground motion. When dam-reservoir system is subjected to a dynamic effect such as earthquake, hydrodynamic pressures in excess of hydrostatic pressures occur on the dam upstream face due to the vibration of the dam and water in the reservoir. These hydrodynamic pressures and the deformation of the dam interact with each other (Perumalswami and Kar, 1973; [1]). Therefore, the fluctuation of reservoir water level and its effects on the seismic response of arch dams must be considered in nonlinear analysis. The magnitude of hydrodynamic effects on an arch dam during earthquakes is based on the water level in the reservoir. The reservoir water level usually varies based on the operational regime and annual weathering conditions and sometimes reservoir level may be decreased to visit the upstream face or to repair probable damages or even because of drought and they can considerably affect the response of arch dams during earthquakes.
Importance of hydrodynamic effects on behaviour of dams subjected to earthquake ground motions has long been recognized by investigators such as Chopra 1967[2], Chwang 1978[3], Sharan 1985[4], Fok and Chopra 1986[5], Hung and Wang 1987[6], Hung and Chen 1990[7], Tsia and Lee 1990 [8], Dominguez and Maeso 1993[9], Fahjan et al. 2003 [10], Bhattacharyya 2007 and [11,12], Mitra and Sinhamahapatra 2008 [13], Najib and Fei 2009 [14]. They have studied dam-reservoir system with different geometrical features of the dam and the reservoir. In Fok and Chopra studies, the complex frequency functions were used to define dam-reservoir interactions and the concrete material was assumed to be linearly elastic. They found that hydrodynamic pressures play an important role on dynamic response of arch dams. In these studies, the reservoir was generally considered as empty or full. Proulx et al. 2001[15]  At the present study, reservoir level fluctuation effects on nonlinear seismic response of a high arch dam including dam-reservoir-foundation interaction is investigated utilizing Eulerian approach in fluid domain and Lagrangian approach in solid domain. For this purpose the finite element model of dam-reservoir-foundation was developed. In this model dam body is modelled defining concrete materials property and considering smeared crack method (NL-SCM) and dam-reservoir-foundation interaction is taken into account (FSI-SSI). The heights of water in reservoir are considered in four different levels so that in the first 28.9%, in the second 52.1%, in the third 79.9% and in the forth 96.9% of the total height of the simulated dam is filled in order to study responses of the body.

Modeling Material Nonlinearity in Mass Concrete
In theories based on continuum crack models (CCM), the fracture is idealized to propagate as a blunt front smeared over an entire element or a certain band width of the element. After initiation of fracture process, determined by a suitable constitutive model, the pre-crack material stress-strain relation is replaced by an orthotropic relation with material reference axis system aligned with the fracture direction. The tension stiffness across the crack plane is either eliminated suddenly or a gradual stress-release criterion is applied. Thus, only the constitutive relation is updated with propagation of cracks and the finite element mesh is kept unchanged. The advantage of the models lies in its simplicity and cost effectiveness, although the physical nature of crack representation is questionable (Gunn 2001[29]).
One of the most important branches of CCM is nonlinear-smeared crack model (NL-SCM), which considers concrete nonlinear characteristics. These kinds of model must be capable of describing material constitutive behaviour before and after cracking. In the present model it is assumed that concrete material is initially (before cracking) isotropic and linear. The stress-strain matrix is defined by the following equation: where E is isotropic Young's modulus for concrete and ν is Poisson's ratio. The behaviour of concrete elements is assumed linear elastic until it reaches ultimate strength ( Figure 1) and the elasticity modulus of concrete is considered as the average, E, instead of the linear actual, E o . Stress increases linearly along with an increase in strain. In this step, each reloading of elements leads to elastic returning of strain. The tension failure of concrete is characterized by a gradual growth of cracks, which join together and finally disconnect larger parts of the structure. It is a usual assumption that forming cracks is a brittle process and the strength in the tension-loading direction abruptly goes to zero after big cracks or it can be simulated with gradually decreasing strength.
In concrete, cracking occurs when the principal tensile stress in any direction lies outside the failure surface. After cracking, the elastic modulus of the element is set to zero in the direction parallel to the principal tensile stress direction. In this model cracking is permitted in three orthogonal directions at each integration point. When cracking occurs at an integration point, the stress-strain relation is modified by defining a weak plane normal to the crack direction, which is unable to endure any tensile stresses. The presence of crack at an integration point and in special direction, represent through modification of stiffens matrix by exerting shear transfer coefficient in cracked plane. Also crushing occurs when all principal stresses are compressive and lie outside the failure surface; in this condition, the elastic modulus is set to zero in all directions.
Based on the fact that concrete has been cracked in one, two or three orthogonal directions, the stiffness matrix can be represented in following forms: I) Concrete has been cracked in one direction and the crack is open: where E s is the secant modulus of elasticity, υ is Poisson's ratio and β open is open shear transfer coefficient and defined as the factor that represents shear strength reduction across the cracked face.
II) Concrete has been cracked in one direction and the crack is closed: where β close is closed shear transfer coefficient. The shear transfer coefficient, β, represents conditions of the crack face. The value of β ranges from 0.0 to 1.0, with 0.0 representing smooth crack (complete loss of shear transfer) and 1.0 representing a rough crack (no loss of shear transfer).
The value of β used in many studies of reinforced concrete structures, however, varied between 0.05 and 0.25 (Bangash 1989, [30]). Convergence problems were encountered at low loads with β less than 0.2.
III) Concrete has been cracked in two directions and the cracks are open: IV) Concrete has been cracked in two directions and both cracks are closed: VI) Concrete has been cracked in three directions and all cracks are closed. In this situation equation 5 can be written again. It must be noted that all above stress-strain relations are written in coordinate system that is parallel to principal stress directions.

Finite Element Formulation of Reservoir
Hydrodynamic pressure distribution in reservoir is governed by pressure wave equation. Assuming that water is linearly compressible and neglecting viscosity, small-amplitude irrotational motion of water is governed by 3D wave Dams Considering Material Nonlinearity equation given as (Ghaemian and Ghobarah 1998[31]): where P is hydrodynamic pressure and c o is velocity of pressure wave in water. For earthquake excitation, condition at boundaries of the dam-reservoir, the reservoir-foundation, the reservoir-far-end and the reservoir-free-surface as shown in Figure 2 are governed as presented in following sections.

Dam-Reservoir Boundary Condition
At the surface of fluid-structure, there must be no flow across the interface. This is based on the fact that face of concrete dams is impermeable. In the following equation, superscript "s" refers to the structure.
where a s n is normal acceleration of the dam body on the upstream face and n is normal vector on the interface of the dam-reservoir outwards the dam body and finally, ρ is the density of reservoir water.

Reservoir-Foundation Boundary Condition
If there is no energy absorption on reservoir bottom, the same boundary condition represented in equation 8 can be used. Reservoir bottom absorption effect is implemented as: where q is admittance coefficient and relation between q and α (wave reflection coefficient at the reservoir bottom and sides) is expressed as:

Free-Surface Boundary Condition
In high dams, surface waves are negligible and hydrodynamic pressure on the free surface is set to be zero.

Reservoir-Far-End Boundary Condition
For modelling far-end truncated boundary, viscous boundary condition (called as Sommerfeld boundary condition) is utilized to absorb completely the outgoing pressure waves given as:

Finite Element Model of Dam-Reservoir-Foundation System
DEZ is a high double curvature arch dam with a peripheral joint separating the main body from its concrete saddle structure (PULVINO foundation). Its total height is 203m but the height of the dam above concrete plug (the simulated dam) is 194m. Also the height of main body (without PULVINO) is 186m. Thickness of the dam at the crest is 4.5m and its maximum thickness at the base is 21m.

Finite Element Model
The finite element model prepared for dam, foundation rock and reservoir is presented in Figure 3. The model consists of 792 8-node solid elements for modelling concrete dam and saddle and 3770 8-node solid elements for modelling foundation. The 8-node solid elements have three translational degrees of freedom (DOF) at each node. The reservoir is modelled using 3660 8-node fluid elements. Utilized Eulerian fluid elements have three translational DOFs and one pressure DOF in each node. It should be noted that translational DOFs are active only at nodes that are on the interface with solid elements.

Material Properties
Material properties for mass concrete and foundation are described in Table 1. Furthermore, reservoir water density is assumed 1000kg/m 3 , sound velocity is 1440m/s in water and wave reflection coefficient for reservoir around boundary is 0.8 conservatively.

Loading History
Applied loads on the system are dam body self-weight, hydrostatic pressure in various water levels and earthquake load. It should be noted that thermal loads is not considered in the current study because of the main object of the study which is considering fluctuation effects of reservoir water level. The β-Newmark method is utilized to solve the coupled nonlinear problem of dam-reservoir-foundation model.
The system is excited at the foundation boundaries using TABAS earthquake records. Total time duration of earthquake is 48.9s. For specifying significant duration of excitation, Arias intensity on Husid diagram was used and then time interval between 10% and 80% of Arias intensity (totally 70%) was selected (Trifunac and Brady 1975[32]). With this method significant time duration is obtained to be 10s (as shown in Figure 4). Moreover structural damping is taken to be 7% of critical damping.

Numerical Results
Considering the great importance of natural vibration modes, 30 first periods of dam-reservoir-foundation system in four different water levels for monolithic dam body were extracted (as shown in Figure 5). It is found that the amount of period grows with arising reservoir water level which is because of increasing total mass of the system and differences are meaningful for ten first modes. Initial (static) stresses and displacements have strong effects on seismic response of arch dams. Therefore, static analysis of dam-reservoir-foundation rock system under self-weight considering stage construction and impounding is carried out to establish the initial condition for seismic analysis. Subsequently, nonlinear seismic analyses of the system are performed.
Water levels in the reservoir are considered as 56, 101, 155 and 188m that stand to level I, level II, level III and level IV to investigate reservoir level variation effects on nonlinear seismic response of the selected arch dam. Figure  6 shows different reservoir levels.

Study along Crest Nodes
Time-history of crest displacement at the central cantilever resulted from four conducted analyses corresponding to four water reservoir levels are shown in Figure 8. Table 2 represents initial (static) displacement along the crest in stream direction for different reservoir water level. Considering that the positive direction of y axis is toward downstream (DS), it becomes evident that in low levels of the reservoir, the main body of dam tends to move in upstream (US) direction.  Table 3 represents maximum displacement of the crest in US/DS direction during seismic excitation. As can be found, raising the reservoir level leads to increasing the crest displacement in the downstream direction and decreasing it in the upstream. However, this variation is negligible between reservoirs levels I and II.  Table 4 represents absolute maximum velocity and acceleration of the crest centre in the stream direction. It is found that arising reservoir level from I to IV leads to decreasing in velocities and accelerations. However, in high levels, these changes are negligible. The reason can be attributed to free movement of the dam body in low water levels.  Figure 9. Non-concurrent displacement envelope for upstream nodes of central cantilever; NL-SCM Figure 9 represents non-concurrent displacement envelope for upstream nodes of the central cantilever. As shown, displacement in lower half of the central cantilever grows significantly with increasing of reservoir level. Displacements corresponding to levels I and II are approximately the same at upper half of the cantilever. In addition, maximum displacements are the same for levels I, II and III at the crest of central block, while the crest displacement corresponding to level IV of the reservoir water is about 65% more than the others.

Height of central cantilever (m)
level IV level III level II level I Figure 10. Non-concurrent velocity envelope for upstream nodes of central cantilever; NL-SCM Figure 10 represents non-concurrent velocity envelope along the height of the central cantilever. As shown rate of velocity variations in crown cantilever remains same for different levels of the reservoir water. In lower half of the cantilever maximum velocities are approximately the same in all considered cases. Although there is little difference in upper half, rate of growth is similar. A sudden increase can be seen in the height of 140m for all cases but at last at crest of the central block, velocities corresponding to low water levels are more than those obtained for high water levels.  Figure 11 shows non-concurrent acceleration envelope corresponding to nodes considered for the two recent figures. As shown rate of acceleration variations along the crown cantilever are the same for all cases like as shown for velocity. Rate of acceleration increases slowly, up to 150m but a sudden growth occurs in this height for all cases so that accelerations in 194m are about triple of accelerations in 150m. Finally at the crest maximum acceleration is corresponding to minimum considered level of reservoir and vice versa. Table 5 gives hydrostatic, absolute maximum value of hydrodynamic and total pressures at base of the main dam body for different reservoir levels. In addition, Figure 12 shows time-history of total pressure at base of the dam for levels I to IV. As shown absolute maximum of hydrodynamic pressure in level IV is more than the rest and in level I is less than the others. Between these levels hydrodynamic pressure varies according to previous maximum and minimum amounts. On the whole the total pressure grows with increasing reservoir level.   Figures 13 and 14 show non-concurrent envelopes of maximum and minimum principal stresses experienced by dam body nodes. It should be noted that these parameters must be limited to dynamic tensile and compressive strength of mass concrete. As can be seen, the extensions of high tensile stressed regions are more when reservoir water level is decreased. In addition, due to intensive arch action of the dam body in upper parts due to high curvature of horizontal arches, cantilevers act like as fixed support beams and extension of high tensile regions on the downstream face is expected as can be seen in Figure 13. On the other hand, due to arch action of the dam body, raising reservoir water level leads to increasing compressive stresses as shown in Figure 14. However, as expected from arch dams' general behaviour, there is a large margin between experienced compressive stresses and dynamic compressive strength of mass concrete.

Study of Stresses and Crack Profiles within the Dam Body
Cracked areas on the upstream and downstream faces at considered reservoir levels are shown in Figure 15

Discussion on Results
As shown in Figure 16, comparing displacements resulted from NL-SCM with linear elastic (LE) model along the height of central block; rate of growth along with reservoir impounding is the same for both models. Displacements for levels III and IV coincident to each other in both models but by reservoir dewatering this correspondence losses so that for levels I and II, the upper half of the cantilever (dry region) experience larger displacements in NL-SCM in comparison with those obtained from LE. At the crest, percentage of displacement discrimination for levels I and II are about 11% and 4%, respectively. It can be shown that the same conclusion can be obtained comparing velocities and accelerations resulted from LE and NL-SCM models. For discussion on stress distribution and crack profiles within the dam body, it should be reminded that in the current study the William-Warnke failure surface in principal stress state is utilized to identify state of stressed nodes.
Based on obtained results dewatering the reservoir leads to extension of high tensile stress region on downstream face of the dam body, which is due to acting cantilevers like as semi-fixed supports beams. However, raising the reservoir LE / level I NL-SCM / level I level causes higher compressive stresses due to dome like action of the body. This conclusion is confirmed looking at extracted crack profiles from NL-SCM model.

Conclusions
In the present study, effect of reservoir water level fluctuation on nonlinear seismic behaviour of high concrete arch dams is considered. The reservoir is modelled as a compressible medium and the foundation rock is assumed to be mass-less. Mass concrete is assumed to have tensile cracking and compressive crushing ability. In addition, for identifying stress state of each Gaussian point within elements William-Warnke failure surface in stress space is utilized. For case study, DEZ double curvature arch dam in Iran is considered.
Based on extracted results from analysing NL-SCM models reservoir water fluctuation affects intensely on seismic response of high concrete arch dams in such a way that raising the reservoir water level leads to lowering extension of high tensile stressed regions and increasing compressive stresses within the dam body. In fact, based on the results, it can be confirmed that the worst case for high arch dams occurs in low water levels.