Evaluation of the Pseudostatic Analyses of Earth Dams Using FE Simulation and Observed Earthquake-Induced Deformations: Case Studies of Upper San Fernando and Kitayama Dams

Evaluation of the accuracy of the pseudostatic approach is governed by the accuracy with which the simple pseudostatic inertial forces represent the complex dynamic inertial forces that actually exist in an earthquake. In this study, the Upper San Fernando and Kitayama earth dams, which have been designed using the pseudostatic approach and damaged during the 1971 San Fernando and 1995 Kobe earthquakes, were investigated and analyzed. The finite element models of the dams were prepared based on the detailed available data and results of in situ and laboratory material tests. Dynamic analyses were conducted to simulate the earthquake-induced deformations of the dams using the computer program Plaxis code. Then the pseudostatic seismic coefficient used in the design and analyses of the dams were compared with the seismic coefficients obtained from dynamic analyses of the simulated model as well as the other available proposed pseudostatic correlations. Based on the comparisons made, the accuracy and reliability of the pseudostatic seismic coefficients are evaluated and discussed.


Introduction
The seismic stability of earth structures has been analyzed by pseudostatic procedures for many decades in which the effects of an earthquake are represented by constant horizontal and/or vertical accelerations. Stability is expressed in terms of a pseudostatic factor of safety calculated by limit equilibrium procedures. Limit equilibrium analyses consider force and/or moment equilibrium of a mass of soil above a potential failure surface. The first explicit application of the pseudostatic approach to the analysis of seismic slope stability has been attributed to Terzaghi [1]. In their most common form, pseudostatic analyses represent the effects of earthquake shaking by pseudostatic accelerations that produce inertial forces which act through the centroid of the failure mass.
The results of pseudostatic analyses are critically dependent on the value of the seismic coefficient. Selection of an appropriate pseudostatic coefficient (particularly ℎ ) is the most important, and the most difficult, aspect of a pseudostatic analysis. The seismic coefficient controls the pseudostatic force on the failure mass, so its value should be related to some measures of the amplitude of the inertial force induced in the potentially unstable material. If the slope material was rigid, the inertial force induced on a potential slide would be equal to the product of the actual horizontal acceleration and the mass of the unstable material. This inertial force would reach its maximum value when the horizontal acceleration reached its maximum value. In recognition of the fact that actual slopes are not rigid and that the peak acceleration exists for only a very short time, the pseudostatic coefficients used in practice generally correspond to acceleration values well below the maximum value. Terzaghi [1] originally suggested the use of ℎ = 0.1 for sever earthquakes (Rossi-Forel IX), ℎ = 0.2 for violent and destructive earthquakes (Rossi-Forel X), and ℎ = 0.5 for catastrophic earthquakes. Seed [2] listed pseudostatic design criteria for 14 dams in 10 seismically active countries and 12 required minimum factors 2 The Scientific World Journal of safety of 1.0 to 1.5 with pseudostatic coefficients of 0.10 to 0.12. Marcuson [3] suggested that appropriate pseudostatic coefficients for dams should correspond to one-third to onehalf of the maximum acceleration, including amplification or deamplification effects, to which the dam is subjected. Using shear beams models, Seed and Martin [4] and Dakoulas and Gazetas [5] showed that the inertial force on a potentially unstable slope in an earth dam depends on the response of the dam and that the average seismic coefficient for a deep failure surface is substantially smaller than that of a failure surface that does not extend far below the crest. Seed [2] also indicated that deformations of earth dams constructed of ductile soils with crest accelerations less than 0.75 g would be acceptably small for pseudostatic factors of safety of at least 1.15 with ℎ = 0.1 ( = 6.5) to ℎ = 0.15 ( = 8.25). This criteria would allow the use of pseudostatic accelerations as small as 13 to 20 percent of the peak crest acceleration. Hynes-Griffin and Franklin [6] applied the Newmark sliding block analysis to over 350 accelerograms and concluded that earth dams with pseudostatic factors of safety greater than 1.0 using ℎ = 0.5 max / would not develop dangerously large deformations.
As can be seen from above discussions, there are no hard and fast rules for selection of a pseudostatic coefficient for design. However, it seems that the pseudostatic coefficient should be based on the actual anticipated level of acceleration in the failure mass and that it should correspond to some fractions of the anticipated peak acceleration, although engineering judgment is required for all cases.
Representation of the complex, transient, dynamic effects of earthquake shaking by a single constant unidirectional pseudostatic acceleration is obviously quite crude. Detailed analyses of historical and recent earthquake-induced landslides have illustrated significant shortcomings of the pseudostatic approach. Results of pseudostatic analyses of some earth dams (e.g., Upper San Fernando dam, Lower San Fernando dam, Sheffield dam, and Tailing dam) show that pseudostatic analyses produced factor of safety well above 1.0 for a number of dams that later failed during earthquakes. Romo and Seed [7] collected many of the destructed dams since 1900 to 1980 which had been designed using pseudostatic method. These are summarized and listed in Table 1. These cases illustrate the inability of the pseudostatic method to reliably evaluate the stability of slopes susceptible to weakening instability. Nevertheless, the pseudostatic approach can provide at least a crude index of relative, if not absolute, stability.
Despite the above-mentioned limitations, the pseudostatic approach has a number of attractive features. The The Scientific World Journal 3 analysis is relatively simple and straightforward. Indeed, its similarity to the static limit equilibrium analyses routinely conducted by geotechnical engineers makes its computations easy to understand and perform. It produces a scalar index of stability (the factor of safety) that is analogous to that produced by static stability analyses. It must always be recognized, however, that the accuracy of the pseudostatic approach is governed by the accuracy with which the simple pseudostatic inertial forces represent the complex dynamic inertial forces that actually exist in an earthquake. Difficulty in the assignment of appropriate pseudostatic coefficients and in interpretation of pseudostatic factors of safety, coupled with the development of more realistic methods of analysis, has reduced the use of the pseudostatic approach for seismic slope stability analyses. Methods based on evaluation of permanent slope deformation are being used increasingly for seismic slope stability analysis.

Upper San Fernando Dam Geometry and Performance.
The Upper San Fernando dam located northwest of Los Angeles and north of the Lower San Fernando dam was built in 1922 using semihydraulic fill technique [8]. The dam was about 24.4 meters high and was constructed upon 15.5 m of alluvial deposits overlying bedrock. The 1971 San Fernando earthquake had a moment magnitude of 6.7 and an epicenter about 13 km from the dam site. The peak horizontal acceleration at the dam site was estimated to be around 0.6 g. Several longitudinal cracks were observed along almost the full length of the dam on the upstream slope slightly below the preearthquake reservoir level. The crest of the dam settled 0.76 m down and moved 1.5 m downstream. The maximum amount of horizontal displacements was about 2 m. Sand boils below the toe and increased water levels in three standpipe piezometers suggested that soil liquefaction had occurred. Water overflowed from two of piezometers. The reservoir level at the time of earthquake was at elevation 369.51 m, 1.83 m below the crest of the dam.

Material Properties of Upper San Fernando Dam.
As shown in Figure 1, the dam is divided into 9 different layers, each representing a different soil material zone. These soil units of dam section have been classified using (1) suggested by Sawada and Takahashi [9]. The equation indicates the effect of height variation on the shear wave velocity and shear modulus of the layers materials located at different levels of the dam section. Consider where is the shear wave velocity and is the dam height below the crest. The relationship between shear wave velocity, shear, and Young's module ( and ) is as follows: where ] is the Poisson ratio. The layers' soil parameters associated with the nine soil zones obtained using the above equations and detailed results given by Seed et al. [8] are listed in Table 2.  is an earth dam located about 33 kilometers northeast of the epicenter. It was minor damage that did not affect the structural safety and the water storage functions of the dam. But this was the first time that an embankment dam designed based on design standards and filled by the engineered rolling compaction was damaged in this way in Japan. The sliding failure zone with a depth of 1.5 to 2 m was confirmed and the length in the dam axis direction of the sliding was about 100 m ( Figure 2). The Kobe University motion is the earthquake motion observed during the Kobe earthquake. The observation station is located on a weathered rock with an S-wave velocity of 340 m/s and 24 km far from the epicenter. The Kobe University motion has a peak acceleration of 270.4 gal and duration of 20 sec. During the earthquake, the reservoir water level was at the top of the sliding failure block. The level differences on the top of the sliding block were from 1 to 1.5 m.  Figure 3, the dam is divided into 24 different layers, each representing a different soil material zone. These soil units of dam section have been also classified using (1). The layers' material properties associated with the 24 soil zones determined using the mentioned equations and detailed results given by Sakamoto et al. [10] are listed in Table 3.

Study Methodology
Most engineers consider the seismic coefficient as a means of designating the magnitude of a static force which is equivalent in effects (i.e., produces the same deformations of the earth dam) to the actual dynamic inertia forces induced by the earthquake. But how would the seismic coefficient denoting this equivalent static force be determined? It would seem that The Scientific World Journal 5 Figure 4: Equivalent static forces acting at layers gravity center. the determination of an appropriate value would necessarily involve two steps: (1) determination and specification of deformations and degree of instability of dam induced by the earthquake; (2) evaluation of equivalent static force with the capability to make the same displacements or instabilities.
It would appear that any attempt to select a final value of such a seismic coefficient without going through step (1) and without a large backlog of experience to guide the selection could have little reliable basis.
In order to determine exact results for stage (1), it will be preferable to utilize dynamic analyses based on finite element method, and hence the Plaxis software seems to be an appropriate choice. High accuracy of dynamic analysis puts it at high point of view. The results obtained from twodimensional dynamic analyses of dams under corresponding earthquake, such as horizontal and vertical displacements, almost justify the observed displacements. Then an equivalent static force is determined for each layer and seismic coefficient is obtained for those layers. In order to reach this aim, the static forces were activated to each layer's gravity center (as shown in Figure 4) and displacements and dam deformations were gained. The importance of this study shines in evaluating the varying seismic coefficient for dams and that is relevant to differentiation of each layer's coefficient of the dam. Assuming a constant seismic coefficient would be applicable for rigid structures and using this current method for earth dams which have not rigid-body response is not rational. Destruction of Lower San Fernando dam and Oshima Tailing dam confirms the invalidity of pseudostatic analysis with constant coefficient, since both of them had been designed using pseudostatic analysis having seismic coefficients of 0.15 and 0.2, respectively. In this study, the equivalent seismic coefficients for different soil zones of Upper San Fernando and Kitayama dams have been determined using two-dimensional dynamic analyses and large backlog, and then the results have been compared with the design seismic coefficient (0.15) of the dams.

Dynamic Analyses
Most of the problems encountered in the area of geotechnical engineering such as retaining walls, tunnels, earth dams, and embankments are studied using two-dimensional dynamic analyses based on the finite element method (FEM) which is one of the available powerful numerical methods. The fundamental stages required to create a FE model include selecting an appropriate element, dividing the model into elements and nods, extending equations of each element and determining element's stiffness matrix, combining element's matrix, and creating a single matrix for model. Elements movement equation is given by in which [ ] is the whole mass matrix, [ ] is the whole damping matrix, [ ] is the model nods axial movement, and { ( )} is the axial force of model points.
One of the current methods used to solve the movement equation is the Newmark step-by-step method. Newmark [11] provided this method for dynamic analysis of earthquake loading. In this method, displacement and velocity are determined using the following equations: where Δ is time pace and and are controlling parameters for numerical integration accuracy, according to the implicit the Newmark scheme [12]. In order to obtain a stable solution, these parameters have to satisfy the following condition: In the classical Lagrange method in [13], = 0.5 leads the calculations to rational results. Despite Newmark's damping method, taking advantage of = 0.6 and = 0.3025 values, in this study, average acceleration method is being used to solve movement equations, as well as Newmark's method. Special boundary conditions have to be defined in order to avoid the spurious reflections of the waves on the model boundaries. These boundaries are based on the Lysmer-Kohlmeyer model. According to this model, the normal and shear stress components absorbed by a damper are determined as follows: where is the mass density, is the shear wave velocity, are the longitudinal wave velocity,̇anḋ= velocity of particle motion in the direction of and , respectively, and 1 and 2 are relaxation coefficients used to improve the wave absorption on the absorbent boundaries. 1 corrects the dissipation in the direction normal to the boundary and 2 in the tangential direction. The research and experience findings recommend to choose 1 = 1 and 2 = 0.25 for the best results [14].

Dam Simulation
The process begins with specifying the clusters (Upper San Fernando dam involved 9 clusters and Kitayama dam involved 24 clusters) and defining the properties relevant to   each cluster. Figure 5 shows the generated mesh sections and Figure 6 exhibits upstream water level and phreatic line of the dams.
The numerical calculations using Plaxis software involve 3 phases. First phase is dam plastic analysis conducted for the time when the construction is over. Second phase includes dam plastic analyses under own body load and finally the last step consists of dynamic analysis under earthquake loading. The third phase loading is applied in the form of a file (accelerogram) input to the program. The whole deformations and horizontal and vertical displacements of the dams are obtained in the output of the program, but considering the importance of horizontal displacements, and for the sake of space saving, only the whole deformations and horizontal displacement diagrams are illustrated in Figures 7 and 8.

Stress-Time Analyses
Before starting calculation step, stress points are chosen on cross section of the dam. These points are located in the direction with three points at each level as shown in Owing to the generation of numerous curves and again for the sake of space saving, only some stress-time curves related to points at various levels of the sections are provided. Figures 11 and 12, respectively, show these envisaged curves for Upper San Fernando and Kitayama dams.

Determination of Equivalent Force
For each curve, a maximum value is deliberated for a period and is considered due to its conservative value. Though the maximum value of each curve is multiplied by 0.7, The Scientific World Journal 7     Figure 13. The equivalent force for each part of layer across stress point is determined using the product of layer height to approximate stress value. The results are illustrated in Figure 14.
Then each layer's equivalent force can be determined by employing the following equation: in which is force of the part (upstream, middle, and downstream), is the effective length of the part, and is The Scientific World Journal the layer's length. Finally seismic coefficient is calculated by dividing layer force to its weight. The results are summarized and listed in Table 4.

Conclusion
The   respectively. It seems that the rational value for seismic coefficient exists between these evaluated values.
(2) Kitayama Dam (i) The maximum value of seismic coefficient offered by JCOLD is 0.25 which is not applicable in the case of this dam too. (ii) Although the minimum acceleration recorded in the seismograph has a value about 0.221 g but considering the rigid-body response method ( = 0.221), which has a maximum seismic coefficient of 0.47, is not reliable.
(iii) In this case, the maximum and minimum values determined using Ambraseys' method are 0.27 and 0.23, respectively. As mentioned above, the rational value for seismic coefficient exists between these computed values.
Finally this study shows that considering inconstant seismic coefficient in earth dam design is more realistic and rational than considering a constant seismic coefficient. Furthermore the procedure employed in this study can be utilized for evaluation of design seismic coefficient of constructed earth dams designed using pseudostatic analyses.