1. Introduction

Site response analysis is usually the first step of any seismic soil-structure interaction study and the problem refers to predict soil response at selected locations of the profile such at the free surface or foundation depth, resulting from the prescribed seismic excitation. On mathematics point of view, the problem is wave propagation in a continuous medium. Modeling of nonlinear soil behavior in 3D wave propagation exactly is extremely difficult and for this reason, the main response in the soil deposit can properly approximated with 1D or 2D vertical propagation of shear waves (Vs). There are large varieties of computer programs available for solving the 1D wave equation.

Nonlinear site response has discussed by Seed and Idriss (1970). However, linear theory has generally been used to estimate site response (Borcherdt, 1970; Aguirre et al., 1994). With more data being recorded in the near source areas of recent large earthquakes, it has become clear that nonlinear site response can have a major effect (Idriss, 1991; Aki, 1993; Iwasaki and Tai, 1996; Field et al., 1997, 1998; Bonilla et al., (2005). Several studies have estimated at what level of shaking one should expect a noticeable nonlinear siteresponse. Idriss (1991) indicates that if PGA is more than 0.4g at soft soil sites, the effect will be noticeable. Midorikawa (1993) and Beresnev (2002) also estimate nonlinear site response if PGA exceeds 200gal, or peak ground velocity (PGV) exceeds 15cm/s.

2. Geological Setting of Semnan Province

The portion of Alpine belt from Iran in the west to Burma in the east, seismically, is one of the most active intercontinental regions of the world. Semnan region is part of the Iranian plateau that is subjected to many tectonic activities, including active folding, faulting and volcanic eruptions. Semnan province is placed in southern flank of Alborz Mountains and north of Kavir desert of Iran, so its geology setting belongs to Alborz and central Iran structural zones. Semnan fault (passes from north of Semnan) and Attari fault (30 km east of Semnan) are two major structural features with general NE-SW trend which pass through the northern part of Semnan province. Geological data specially the comparison of geological characteristics of northern part (Alborz) and southern part (Central Iran) have demonstrated that in fact there are not contrasting differences between northern and southern parts, in other words northern parts of this province is in fact the folded margin of central Iran. However, block tectonics has caused creation of sedimentary basins in tectonic realms of these fault zones so that the lithology and thickness of lithostratigraphic units of the same age have some differences in both sides of these fault zones.

Figure 1. quick glance of geological features of Semnan province

The historical context of earthquakes in this region has been compiled by Ambraseys and Melville (2005). First-motion data indicate that these earthquakes are almost all generated by reverse slip and subsidiary strike-slip movements (Dehghani and Makris, 1983). The strike-slip component of deformation probably has a complex origin. Earthquakes are characteristically shallow, of large magnitude but are discontinuous, with long recurrence periods (Berberian 1977). Figure1 shows the general geological aspects of the Semnan province.

3. Mathematical Modeling of Nonlinear Soil Behavior

In site response analysis, soil properties including shear modulus and cyclic soil behavior are required. Shear modulus is estimated using field tests such as seismic down hole or cross hole tests. Cyclic soil behavior is characterized using laboratory tests such as resonant column, cyclic triaxial or simple shear tests. The maximum shear modulus is defined as G_max and corresponds to the initial shear modulus. Insitu measurement of V_s using geophysical methods is the best method for measuring the G_max (Rolling et al. 1998). Geophysical methods are based on the fact that the velocity of propagation of a wave in an elastic body is a function of the modulus of elasticity, Poisson ratio and density of material (Hvorslev 1949).

By consideration, a uniform soil layer lying on an elastic layer of rock that extends to infinite depth and the subscripts s and r refer to soil and rock, the horizontal displacement due to vertically propagation harmonic S wave in each material can be written as:

(1)

(2)

u: displacement, ω: circular frequency of the harmonic wave, k^*: complex wave number

No shear stress can exist at the ground surface (zs=0), so

(3)

Where is the complex shear modulus of the soil. Schnabel et al. (1972) explained that within a given layer (layer j); the horizontal displacements for two motions (A and B) may be given as:

(4)

Thus, at the boundary between layer j and j+1, compatibility of displacements requires that:

(5)

Continuity of shear stresses requires that:

(6)

The effective shear strain of equivalent linear analysis is computed as:

(7)

(8)

γ_max: maximum shear strain in the layer, R_γ: strain reduction factor, M: magnitude of earthquake

The motion at any layer can be easily computed from the motion at any other layer (e.g. input motion imposed at the bottom of the soil column) using the transfer function that relates displacement amplitude at layer i to that the layer j:

(9)

The nonlinear hyperbolic model used in this paper was developed by Konder and Zelasko (1963) to model the stress-strain soil behavior of soils subjected to constant rate of loading. The hyperbolic equation is defined as:

(10)

τ: shear stress, γ: shear strain, G_mo: initial shear modulus, τ_mo: shear strength, : reference shear strain

The reference shear strain is strain at which failure would occur if soil were to behave elastically. It has been considered a material constant by Hardin and Drnevich (1972). The reference strain can also be represented as function of initial tangent modulus and undrained shear strength in clays (Mersi et al. 1981). The solution of the wave propagation equation is performed in either frequency or time domain (Park and Hashash, 2004). The soil medium is divided into sub layers with absolute displacement u_j, defined at j^th sub layer interface and with shear stress τ_j, defined at the mid points of each interface. As Kramer (1996) explained, the response of soil deposit under dynamic loading is governed by the equation of motion as bellows:

(11)

The differentiation for a soil divided to N sub layers of thickness Δz and processing for the small time increment Δt is computed by using finite difference method as bellows:

(12)

(13)

: Velocity of the motion, : acceleration

The result of the combination the equations 11, 12 and 13 will get:

(14)

Moreover, it can be simplified as:

(15)

As mentioned above, for the soil surface, the shear stress is equal to zero and boundary condition for each sub layer must be satisfied. For soil, rock boundaries Joyner and Chen (1975) proposed the following equation for soil rock boundaries:

(16)

By using the equations, 15, 16 the boundary conditions are satisfied. Kramer (1996) proposed the shear for each layer as below:

(17)

As above mentioned equations show, the shear stress is computed by using current shear strain and stress-strain history . Thus, the proposed method satisfies the nonlinear and inelastic behavior of soil under earthquake excitation.

4. Data and Analysis Frame Work

The selected region is an earth dam with clay core with crest length of 445 m and width of 10 m. The dam has been built over Cheshmeh Ali River, 12 km northwest of Damghan city. About 200 billion rials have been invested in the project and 135 billion rials are needed for building irrigation and drainage network of the dam. The dam is to control floods and save surplus agricultural water to increase area of farmlands by 1,500 hectares. This dam is located in 54° 14' 31˝ East longitude and 36° 13' 50˝ North latitude. The other characteristics of this dam are given in table (1). The oldest rock outcrops in this region belongs to Barut formation and the others belong to recent sediment and deposits. The lithology of the selected site is a combination of limestone, dolomite, shale parts, dark red sandstones, volcanic portions, pink marl, white quartzite and conglomerate unit. In seismotectonic point of view, according to recent classification, the target area is near and beside the intersection of the Alborz-Kopedagh and Central Iran seismotectonic provinces.

Table 1. characteristics of Damghan earth dam Property Quantity Height 51.5 m (from the bed) and 54.5 m (from the foundation) Effective capacity 12.8 Mm3 Total capacity 21 Mm3 Lake area 1.5 Km2

The analysis program is pointed in figure 2. As indicated, at the first step the field and laboratory test were conducted and then the drilled boreholes were investigated. A total of 23 drilled boreholes with maximum depth of 125m were drilled in the selected site to investigate the geotechnical characteristics of the subsurface layers of the dam area and their results with statistical analysis are indicated in tables (2), (3), (4) and (5).

In this step, the site characterized and the condition of subsurface layers was determined. By Collected and obtained data from the previous step the idealized soil profile for the studied area was determined and the L component of Kahak earthquake with magnitude of 5.7, which happened in Iran (2007), was applied to idealized constructed soil profile of the area. To compute the surface parameters, the authors would be forced to combine several software packages with the generated GUI computer code named as “Abbas Converter V2.01” (figure3) that is written with C Sharp Computer language and it is the developed version of old ones which was introduced by Abbaszadeh Shahri et al., (2009). By the generated code, it is possible to take over the encountered problem and by this code the calculated output parameters such as motion, surface response, spectral acceleration, stress, strain and PGA profile for rigid and elastic half space bedrock were computed and compared to each other as shown in figure 4 and 5. To modify and validate the applied method in this study the damping ratio and shear stress reduction curve of the investigated area was computed and compared by several known curves. As shown in figure

6, 7, 8 and 9 the computed curves have a good agreement with others and this point can certify the validation of the study. For better understand the numerical analysis of this study was indicated in table (6).

Figure 2. Summarized analysis frame work in this study

Figure 3. Start screen of generated computer code

Table 2. drilled boreholes with their location in the investigated site borehole Coordination depth location X -25 Y -401 BH1 3047.7 2555.9 100 River bed BH1-1 3008.12 2600.66 100 River bed BH1-2 3096.67 2530.77 100 River bed BH2 3116.47 2660.90 100 Left bank BH3 3169.91 2690.93 100 Left bank BH4 3308.94 2822.88 60 Left bank BH5 3060.98 2572.14 60 River bed BH6 3022.32 2529.42 125 Right bank BH7 2632.05 2238.09 100 Right bank BH8 2944.84 2517.04 100 Right bank BH9 2874.95 2705.44 100 Left bank BH10 2963.65 2767.24 100 Left bank BH11 3255.15 2468.86 100 Right bank TUA1 2620.06 2695.00 60 Tunnel line TUA2 2951.23 2426.44 50 Tunnel line TUA3 3181.39 2282.41 60 Tunnel line TUA4 2810.07 2657.73 40 Tunnel line NS2 3385.14 3064.00 40 Near spillway NS3 3654.16 2729.94 40 Near spillway

Table 3. statistical analysis of rock quality designation of the left support of the dam borehole depth(m) RQD(Ave) Variation Standard deviation Ave min max BD5 85 38% 0% 84% 685.76 26.19 BD4 25 86% 83% 94% 28.7 5.35 BD6 25 53% 38% 63% --- --- BH2 100 74% 0% 98% 706 26.6 BH3 100 39% 0% 84% 743.8 27.3 BH4 60 66% 39% 85% 229.54 15.2 BH5 60 76% 52% 97% 173.42 13.2 BH10 100 35% 2% 83% 473.3 21.8 NS2 40 73% 22% 90% 518.6 22.8 NS3 40 72% 19% 97% 810.3 28.5 left support --- 57% 26% 88% 582.3 24.13

Table 4. statistical analysis of rock quality designation of the dam foundation borehole RQD (%) Variation Standard deviation Ave min max BH9 44 0 91 992.9763 31.5 BH(1-1) 72 43 99 322.74 18 BH1 79 50 100 272.3026 16.15 BH(1-2) 74 29 94 201.3132 14.2 BH11 58 2 95 697.87 26.42 foundation 65 25 96 497.44 22.3 Upper Miocene conglomerate 58 15 95 671.195 25.91 Elica limes 77 40 97 236.81 15.4

Table 5. statistical analysis of rock quality designation of the right support of the dam borehole RQD(%) Variation Standard deviation Ave min max BD2 51 0 95 959.4103 31 BD3 54 9 95 706.1554 26.6 BH6 66 22 100 568.59 23.85 BH8 57 10 100 926.85 30.44 BD1 52 27 94 367.32 19.2 TUA1 72 0 93 878.0227 29.63 TUA2 54 31 72 180.7 13.44 TUA3 37 3 59 427.61 20.7 TUA4 68 30 59 663.9821 26 BH7 81 24 94 305.757 17.5 BD7 90 82 99 22.1 4.7 Right support 61 18 93 552.25 23.5 Around the dam body 56 13 97 702.6 26.5

Figure 4. Comparison between elastic and rigid half space computed time history

Figure 5. Comparison between elastic and rigid half space computed response spectra

Figure 6. Comparison of obtained shear modulus curve with known ones

Figure 7. Comparison of obtained damping ratio curve with known one

Figure 8. Computed shear modulus curve for various depth of idealized soil profile

Figure 9. Computed damping ratio curve for various depth of idealized soil profile

Table 6. numerical computed results of this study for the selected site Computed on surface Input condition parameter 1.629g(0.67s) 1.586g(0.67s) Elastic response 1.634g(0.67s) 1.588g(0.67s) Rigid 0.1542%(59s) 0.000866%(59s) Elastic Strain 0.1664%(59s) 0.000875%(59s) Rigid 10.18(13.82Hz) ----- Elastic Amplification factor 15.65(13.74Hz) ----- Rigid 1.609g(0.67s) ----- Elastic Spectral acceleration 1.604g(0.66s) ----- Rigid -0.429g(59s) -0.389g(59s) Elastic motion -0.431g(59s) -0.39g(59s) Rigid 1.07(59s) 0.654(59s) Elastic stress 1.111(59s) 0.66(59s) Rigid 0.589(0.63Hz) 0.63(0.586Hz) Elastic Fourier amplitude 0.59(0.63Hz) 0.586(0.63Hz) Rigid 66.5(11.99Hz) 2.106(11.99Hz) Elastic Fourier amplitude ratio 80.1(11.99Hz) 22.48(26.16Hz) Rigid

5. Conclusions and Discussion

This study tried to follow in conducting a meaningful nonlinear site response and amplification study. After evaluating the accelerograms at the bedrock, the ground response analysis at the surface, in terms of time history and response spectra, has been obtained by nonlinear standard hyperbolic model and authors have been trying to find a practical and appropriate solution for ground response analysis under earthquake forces for the selected site. For achieve to this aim, a user graphical interface software was produced and developed by authors and on base of this generated code, a new practical reliability geotechnical based method procedure has been presented. To validate of this code and to prove the accuracy and capability of the proposed procedure, a case study on ground response analysis of a site in Damghan earth dam in Semnan province of Iran, during the Kahak earthquake (2007) is executed. This study shows that this generated computer code can be applied as a strong and reliable tool for future. The obtained results in this study shows a good agreement by known applicable procedures where as in modulus degradation can recognized between back calculation and some of previous laboratory test results for soils in the site of dam. Distribution of maximum acceleration along the depth and spectrum ratios has proved that rigid half space bedrocks compute larger peak accelerations. The PGA value at the ground surface obtained from the generated computer code can use to prepare the PGA map of selected site. They are not distributed uniformly due to variation in the soil profile at various locations. More that this PGA is comparable to obtained peak horizontal acceleration values using SPT data and the shape of variation of peak acceleration with depth are similar to the SPT data. The calculated amplification factor ranged in elastic and rigid condition can be used to prepare the amplification map of Damghan earth dam.