# Broad-band simulation of M7.2 earthquake on the North Tehran fault, considering non-linear soil effects

Broad-band simulation of M7.2 earthquake on the North Tehran fault, considering non-linear soil... Summary The North Tehran fault (NTF) is known to be one of the most drastic sources of seismic hazard on the city of Tehran. In this study, we provide broad-band (0–10 Hz) ground motions for the city as a consequence of probable M7.2 earthquake on the NTF. Low-frequency motions (0–2 Hz) are provided from spectral element dynamic simulation of 17 scenario models. High-frequency (2–10 Hz) motions are calculated with a physics-based method based on S-to-S backscattering theory. Broad-band ground motions at the bedrock level show amplifications, both at low and high frequencies, due to the existence of deep Tehran basin in the vicinity of the NTF. By employing soil profiles obtained from regional studies, effect of shallow soil layers on broad-band ground motions is investigated by both linear and non-linear analyses. While linear soil response overestimate ground motion prediction equations, non-linear response predicts plausible results within one standard deviation of empirical relationships. Average Peak Ground Accelerations (PGAs) at the northern, central and southern parts of the city are estimated about 0.93, 0.59 and 0.4 g, respectively. Increased damping caused by non-linear soil behaviour, reduces the soil linear responses considerably, in particular at frequencies above 3 Hz. Non-linear deamplification reduces linear spectral accelerations up to 63 per cent at stations above soft thick sediments. By performing more general analyses, which exclude source-to-site effects on stations, a correction function is proposed for typical site classes of Tehran. Parameters for the function which reduces linear soil response in order to take into account non-linear soil deamplification are provided for various frequencies in the range of engineering interest. In addition to fully non-linear analyses, equivalent-linear calculations were also conducted which their comparison revealed appropriateness of the method for large peaks and low frequencies, but its shortage for small to medium peaks and motions with higher than 3 Hz frequencies. Earthquake ground motions, Site effects, Wave propagation 1 INTRODUCTION The eastern segment of the North Tehran fault (NTF), at the southern boundary of the Alborz range and North of the city of Tehran, can produce earthquakes with magnitudes M7.2–7.4 (Anderson 1997; Ritz et al. 2012). The return period for such earthquakes is estimated about 1000–1150 yr and the last earthquake which is historically corresponded to the NTF activity has occurred about 900 yr ago (Berberian & Yeats 1999). Considering the old buildings at some parts of Tehran and the highly populated city, occurrence of such event can be catastrophic. Hence, attaining a quantitative estimate of the ground motions at various parts of the city due to the next M7.2 earthquake is extremely necessary and require comprehensive researches. Results of such studies can be used in risk assessment and municipal decision makings of the Tehran metropolitan. On the other hand, ground motion prediction equations suffer from increasing uncertainties at near-fault regions and their predictions are not reliable for short distances (Dalguer & Mai 2012). In a previous paper (Majidinejad et al. 2017), we presented results for 30 dynamic rupture scenarios of an M7.2 earthquake on the NTF for frequencies up to 2 Hz. The simulation also included the Tehran basin's velocity structure. Results revealed the significant effects of the surface waves on the ground motions of the city. In the following study, results are extended up to 10 Hz in order to consider the whole range of the engineering frequencies. Broad-band ground motions are obtained by combining the low frequencies (0–2 Hz) from dynamic simulations with high frequencies (2–10 Hz) from a scattering method. Site local effects, which become considerable at high frequencies, are also included for sites on a cross-basin profile. These effects are investigated by linear, equivalent-linear and fully non-linear soil models and effect of non-linearity on soil response is assessed. 2 BACKGROUND In the previous paper (Majidinejad et al. 2017), we used a 2-D model to simulate the spontaneous dynamic rupture of the NTF and propagation of the seismic waves in the heterogeneous medium of the Tehran basin. Results of that study provided deterministic aspects of the next Tehran M7.2 earthquake. However, considering huge computational cost, our limited knowledge from the earth's heterogeneous crust and also fault's complex behaviour, results of those dynamic simulations were limited to 2 Hz. All broad-band simulations are encountered with such limitations on the available computational resources and lack of knowledge about geophysical properties. The limit on the upper frequency of deterministic methods makes them unable to predict the high-frequency ground motion parameters. On the other hand, empirical-stochastic methods are designed to produce ground motion parameters similar to previously observed earthquakes (Irikura 1986). Several methods have been introduced for generation of high-frequency component of ground motions. Some of them use stochastic synthetics which convolve with appropriate spectrum of the source (Boore 1983, 2003) and some other use wave scattering theory to simulate high-frequency motions (Zeng et al. 1994). However, all these methods are unable to reproduce acceptable low-frequency waveforms. To overcome these limitations, combined broad-band methods emerged. These methods are based on combination of deterministic low-frequency with stochastic high-frequency motions (Mai & Beroza 2003). Taking into account non-linear behaviour of soil in deterministic low-frequency simulations increases the computational costs and time significantly. Moreover, previous studies have demonstrated that the low-frequency ground motions are not affected considerably by non-linear response of the shallow soil layers. In contrast, it is proven that non-linear soil effects are crucial for high-frequency contents, especially during large earthquakes (Roten et al. 2012). To consider the soil non-linearities in microzonation and geotechnical engineering studies, it is a well-known approach to employ a 1-D non-linear soil column model as a representation of shallow soil layers and calculate its response due to seismic excitation at depth (Hashash & Park 2001; Hartzell et al. 2002). Moreover, since the computations for non-linear behaviour of the soil is time consuming and requires several iterations for each time increment, an equivalent-linear analysis can be used instead. In this method (Seed & Idriss 1969), dynamic soil parameters are kept constant over whole duration of seismic excitation at each iteration of the repetitive linear soil response calculations. The maximum shear strain occurred at each point during excitation is used to determine the appropriate shear modulus and damping ratio for the next iteration (i.e. linear solution). Equivalent-linear method is primarily designed for small strain levels in which soil complex non-linearities are not prevailed. It is assumed that shallow layers of soil, which are not exactly in the vicinity of the fault, meet the criteria for using this method in calculations of the earthquake ground motions as it is used widely by geotechnical engineers (Towhata 2008). The single previous study on the broad-band simulation of the NTF’s probable earthquake was performed by Zafarani et al. (2012, 2013). They produced PGA and PGV maps of the Tehran by combining finite-fault method (Motazedian & Atkinson 2005) results for low-frequency motions and Boore's (2003) stochastic method synthetics for high-frequency component. They considered linear soil amplifications by implementing generic hard rock amplification factors (for the whole region of the study) proposed by Boore & Joyner (1997). They concluded that the highest Peak Ground Accelerations (PGAs) occur on the hangingwall side of the fault with a value approximately equal to 700 cm s−2. The recent study performed by Majidinejad et al. (2017) predicted the highest mean SA values in the frequency range of 0–2 Hz occur on the footwall side of the NTF and that there is a significant increase in the values over the deepest parts of the Tehran basin, due to the constructive interference of body shear and surface Rayleigh waves. However, restriction of these results to the upper limit of 2 Hz, prevents a comparison between them and Zafarani et al.’s (2013) broad-band PGAs. The main goal of this study is to provide broad-band synthetic ground motions for Tehran based on the low-frequency results of Majidinejad et al. (2017) and also to estimate the non-linear soil effects on ground motion intensities using soil profiles in a cross-section of the Tehran basin. 2.1 Low frequency (0–2 Hz) The low-frequency ground motions are obtained by 17 scenarios of spontaneous rupture propagation. The NTF is considered as a purely reverse fault with a 75° dip angle. Each rupture scenario is generated by a self-similar initial stress distribution on the fault plane that produces kinematic fault properties such as rupture velocity and length, slip rise time and final slip distribution corresponding to an M7.2 earthquake. Another important aspect which is considered by the low-frequency models is the deep sedimentary basin beneath Tehran. Fault-normal and vertical components of the low-frequency ground motions reveal remarkable directivity effects of the updip rupture propagation. LF ground motion components are amplified significantly due to the low-velocity layers of the basin on the footwall side, near downtown. Results also show that the formation of Rayleigh waves at the edges of the basin changes the frequency content and duration of the LF components (Majidinejad et al. 2017). 2.2 High frequency (2–10 Hz) In this study, we use the ‘Broad-Band toolbox’ code, developed by Mai et al. (2010), to generate and combine the HF components with the LF motions. The HF motions are produced based on a method proposed by Zeng et al. (1991) and Zeng (1993) for point sources and developed by Mena et al. (2010) for extended faults. The background theory of the method will be shortly described here. As the seismic waves travel from source to stations, they scatter by various distributed heterogeneities in the medium. Such phenomenon causes dispersion of the seismic wave's energy and reduces energy density of the directly arrived waves to site. However, in this case, the reduction of the energy density is different from other dissipative mechanisms (such as viscoelastic or plastic behaviour of the material) and the energy is transferred to secondary reflected waves which appear as coda in time histories. Zeng (1993) used S-to-S backscattering theory and proposed an equation for energy partitioning of a spherical wave into direct, first- and second-scattered and multiple (orders higher than 2) scattered waves. Eq. (1) represents the multiple scattered portion of the wave's energy at the time t due to an impulse excitation in the reference point of an infinite 3-D medium.   \begin{eqnarray} &&{ {E_{{\rm{multiple}}}}(\vec{r},t)}\nonumber\\ &&= \mathop \smallint \limits_{ - \infty }^{ + \infty } \frac{{{e^{i{\rm{\Omega }}}}}}{{2\pi }}d{\rm{\Omega }} \times \mathop \smallint \limits_0^\infty \frac{{{{\left( {\frac{{{\eta _s}}}{k}} \right)}^3} \times {{\left[ {{\rm{ta}}{{\rm{n}}^{ - 1}}\left( {\frac{k}{{\eta + {{{i{\rm{\Omega }}}}/{v}}}}} \right)} \right]}^4} \times {\rm{sin}}\left( {kr} \right)}}{{2{\pi ^2}vr\left[ {1 - \frac{{{\eta _s}}}{k}{\rm{ta}}{{\rm{n}}^{ - 1}}\left( {\frac{k}{{\eta + {{{i{\rm{\Omega}}}}/{v}}}}} \right)} \right]}}dk \end{eqnarray} (1) In this equation, ‘r’ is the source-to-receiver distance, ηs and η are the scattering and total shear wave attenuation coefficients and Vs is the shear wave speed. In the method used for producing the HF motions, the scattering Green's functions are generated by multiplying random wavelets (with average energy of unity) to the energy envelope calculated by eq. (1) for each station (Mai et al. 2010). Mena et al. (2010) developed the method from point source to an extended model by introducing the extended fault as a composition of small subfaults. They calculated the contribution of subfaults to the overall ground motion by taking into account the delay of slip initiation on each subfault due to its distance from hypocentre and also delay due to the subfault to station distance. The HF motion attributed to each subfault is obtained by convolution of the subfault's source time function (STF) to the scattering Green's function. The amplitude of the STF for each subfault is determined by scaling laws and with respect to the seismic moment on that subfault (Irikura 1986). By this method, synthetic high-frequency ground motions that take into account physical and random features of the fault and medium to some extent are produced. However, it should be noted that since eq. (1) is derived for an infinite medium, only body waves are measured and therefore, the effect of surface scattered waves is not reflected in the HF component. The attenuation and scattering coefficients which represent density and size of heterogeneities in the medium and control the coda envelope are set to 0.06 and 0.03 km − 1, respectively, with regard to Mena et al. (2010). We use the frequency-dependent attenuation model proposed by Farahani et al. (2012) for the Tehran basin, as used for the LF simulation in Majidinejad et al. (2017). In order to model the HF site-specific attenuation, kappa coefficient is set to κ = 0.04 s, according to Hassani et al. (2011) results for the region. Dreger et al.’s (2007) STF is used for calculation of high-frequency motions. This slip rate function is compatible with the slip functions generated by dynamic spontaneous rupture models, as it can provide the overall shape and rapid rise time of the dynamic models appropriately. A comprehensive study by Mena et al. (2010) showed that in comparison with traditional boxcar and triangular STFs and also dynamically consistent Yoffe (Tiniti et al. 2005) and Liu et al. (2006)’s STFs, Dreger's STF is the most appropriate, since it has a smooth spectrum and does not contain any spectral holes or abrupt spectral changes at low frequencies. 2.3 Broad band (0–10 Hz) The Broad-Band toolbox combines the LF and HF components of the motion with the method proposed by Mai & Beroza (2003), which provides a solution to the inconsistency of the energy spectrums of LF and HF contents near matching frequency (Thráinsson et al. 2000). The main assumption of the method is that the HF content of motion is not completely stochastic and contains some physical aspects of the rupture propagation and emitted seismic waves, such as arrival phases of P and S waves, magnitude of the motions, etc. This criterion is met by HF ground motions generated by Zeng et al. (1995) and Mena et al. (2010) method. Mai and Beroza's (2003) method searches for the matching frequency within a bandwidth such that the phase spectrum of the HF component matches optimally with the realistic LF phase spectrum. Then amplitude spectrum of the HF component is scaled by a close to unity factor to the LF amplitude at the matching frequency. Since physical characteristics of the earthquake are considered during calculation of the HF component of motions, the scale factor used for fitting the amplitude spectrums is often close to unity. By this method, the matching frequency may alter for each station and direction of motion. In this study, the bandwidth for the matching frequency search is set to 1.8–2 Hz and the maximum ground motion frequency is limited to 20 Hz. However, the resulting broad-band ground motions are low-pass filtered up to 10 Hz in order to facilitate the site effect calculations and avoid the numerical instabilities. Synthetic broad-band ground motions are provided at stations with 1 km distance strides on a north–south profile of the city. 2.4 Tehran alluviums Tehran is spread over three stratigraphic zones. The Southern Alborz Mountains foothills, north of the city, comprise the first zone and are mainly composed of shale, Eocene tuffs, andesite, basalt and pyroclastic from Palaeozoic, Mesozoic and Tertiary eras (Emami et al. 1993). The second zone, which includes Sepayeh Mountains at the east and Bibisharbano Mountains to the south of the city, is made up of Triassic period materials such as limestone, marl limestone and dolomite (Shafiee & Azadi 2007). The Tehran alluvial plain, as the third stratigraphic zone, is located above the Tehran basin and consists of Quaternary formations which are mostly accumulated by erosion and redeposition of former sediments. The alluvium spreads from foothills at the north to the low farmlands at the south of Tehran. It is composed of varying thickness layers with materials of different types and size (Rieben 1955; Jafari et al. 2002). Pedrami (1981) categorized the Tehran plain alluvium to five geological units, including units A and Bn in the north, unit Bs in the south, unit C in the north, west and centre and unit D in centre and south (Fig. 1). Jafari et al. (2001) collected data from 700 drilled boreholes over the study area to investigate the material types of the geological units. Results identified unit A, the oldest alluvial layer, mainly consists of gravel, sand and silt. They proposed unit B to be a conglomerate mixture of cobbles, gravels, sand and silt and also dolomite and limestone at some southern parts. Unit C is characterized by a mixture of pebble, sand, silt and clay and the youngest layer, unit D, is composed of soft silty-clay materials. Shafiee & Azadi (2007) used data from 188 seismic profiles, distributed throughout the city, to classify the geological units of the Tehran plain according to their shear wave velocity. They reported average Vs, 30 for the Tehran plain as in Table 1. Figure 1. View largeDownload slide Map of geological units of Tehran (Shafiee & Azadi 2007). The thick dashed black line represents top edge of the NTF’s rupture plane. The thick black line stands for the N-S profile of the study. Velocity layers of deep basin of Tehran in the N-S profile are shown in the left. Figure 1. View largeDownload slide Map of geological units of Tehran (Shafiee & Azadi 2007). The thick dashed black line represents top edge of the NTF’s rupture plane. The thick black line stands for the N-S profile of the study. Velocity layers of deep basin of Tehran in the N-S profile are shown in the left. Table 1. Characteristics of geological units of the Tehran plain. Geological unit  Constituting materials  $${\bar{V}_{s,30}}({\rm{m}}\,{{\rm{s}}^{ - 1}})$$  PI (%)  A  Conglomerate with silt–sand–gravel and silt–clay mixtures  830  0  B  Cobble, boulder, gravel and sand  770  0  C  Gravel, sand, silt and clay  400  12  D  Silt and clay  280  25  Geological unit  Constituting materials  $${\bar{V}_{s,30}}({\rm{m}}\,{{\rm{s}}^{ - 1}})$$  PI (%)  A  Conglomerate with silt–sand–gravel and silt–clay mixtures  830  0  B  Cobble, boulder, gravel and sand  770  0  C  Gravel, sand, silt and clay  400  12  D  Silt and clay  280  25  View Large 2.5 Soil properties One of the deterministic parameters in non-linear soil behaviour is the plasticity index (PI). In this study, the PI for different soil types are assigned based on the geological units of Tehran alluviums and also based on the amount of clay content in soil composition (Hartzell et al. 2002; Roten et al. 2012). The fine-grained soil in the southern parts of the Tehran consists of clay and silt. Jafari et al. (2002) performed a comprehensive investigation on alluviums of the southern Tehran by evaluating dynamic properties of the fine-grained soils using dynamic laboratory tests and field geoseismic investigations. The laboratory tests included torsional resonant column tests in a shear strain range of 10 − 6–10 − 4and stress-controlled cyclic triaxial tests in a wider range of shear strains (10 − 5–10 − 2) to assess the dynamic properties of the soils collected from sites at the southern parts of the city. They categorized the fine-grained soil in the area to three types with respect to the PI of the samples. Silt specimens were characterized by very low PI (PI < 7); silty clay samples were identified by low index (7 < PI < 15) and medium PIs (15 < PI < 30) were assigned to clay samples. They concluded that an increase of the confining pressure on the soil samples leads to more shear modulus ratios ($${G/{G_0}}$$) and lower damping ratio, however, this effect reduces as the soil plasticity increases. Finally, they have compared their results by Vucetic & Dobry (1991) models and concluded that for the medium PI clay, the results agree very well with 15 < PI < 30 in the model. Ghorbani et al. (2014) proposed the value of PI for the fine silty clay soil of Tehran as 12 per cent. Therefore, based on the main contents of the Tehran geological units, we categorized the Tehran soil to three types of gravel and sand, silty clay and clay and assumed their PIs as (Table 1). In order to assign the appropriate soil properties at each site, we used the reports of JICA (2000) and Keshavarz Bakhshayesh 2003 as well as the geological map of Shafiee & Azadi (2007) to determine the soil type at each depth due to the alluvium profiles and regional boreholes. Since the younger Alluvium of the southern Tehran are composed of remarkable clay contents, we used a PI of 25 per cent for unit D. Geological unit C has lower clay and higher silt, so it is reflected in the lower PI value of 12 per cent. Finally, the oldest A and B units are mostly composed of sand and gravel, therefore the negligible clay content is represented by PI = 0 per cent. Tehran alluviums have experienced many loading and unloading cycles such as fluctuations of the water table and desiccation during deposition, deposition of upper layers, aging effects and also cementation and upward motion of lower deposits due to Quaternary faults movements. So, due to the geological history of Tehran sedimentary basin, it is plausible to consider overconsolidated soil in the region. Several studies are established to estimate the overconsolidation ratio (OCR) of the Tehran alluvium (Jafari et al. 2000). These studies suggest that south of Tehran soils are generally characterized by a low OVR, approximately about 1–1.5. On the other hand, after numerous laboratory investigations, Darendeli (2001) stated that the OCR in comparison with other factors has minor effect on the dynamic properties of soil samples. Therefore, we used the OCR = 1 for all types of the soils in Tehran, both in calculation of at-rest lateral earth pressure coefficient (K0) and Darendeli's equations. Several models have been proposed by researchers to predict the soil non-linear behaviour (Vucetic & Dobry 1991; Darendeli 2001; Amir-Faryar & Aggour 2016). However, regional studies on the non-linear behaviour of various soil types, around the world, show that some modifications to the parameters used in the proposed models are necessary in order to provide better fit between predictions and the local experimental results (e.g. Bay & Sasanakul 2005). We used Darendeli's (2001) model which takes into account the consolidation state, effective stress and PI in calculation of soil dynamic behaviour. The model proposed by Darendeli (2001) also requires the number of cycles (N) and the frequency of excitation ( f). Jafari et al. (2002) set these parameters to N = 10 and f = 1 Hz with respect to the Standard (2003) ASTMD 3999 recommendation and in this study we use similar values for the parameters of the Darendeli's model (following Roten et al. 2012). Selected parameters are verified with the experimental data obtained by Jafari et al. (2002) for the Tehran region (Fig. 2). Figure 2. View largeDownload slide Comparison of Darendeli's (2001) model and experimental results for the Tehran soil (Jafari et al. 2002) for two typical soil samples with medium and low plasticity indices. Figure 2. View largeDownload slide Comparison of Darendeli's (2001) model and experimental results for the Tehran soil (Jafari et al. 2002) for two typical soil samples with medium and low plasticity indices. Due to the computational limits, the minimum shear wave velocity in the LF spectral element model was constrained to 600 m s−1 on soft quaternary and 800 m s−1 on Tertiary deposits (Majidinejad et al. 2017). The velocity model used for shallow soil profiles is proposed based on reflection and borehole tests (JICA 2000) and is consistent with the velocity model used for the LF 2-D simulations. JICA (2000) proposed 40 shallow profile models over Tehran, categorized according to shear wave velocity, material type, thickness of the layers and the depth to the geotechnical bedrock. The maximum velocity of the JICA profiles at deepest layers are 550 and 800 m s−1 for quaternary and Tertiary deposits, respectively, which coincide with the minimum velocities in 2-D spectral element model uppermost layers. 2.6 Analysis model In order to perform the site response analysis with linear, equivalent-linear and fully non-linear models, we used the code DEEPSOIL, developed by Hashash et al. (2012), which utilizes multiple lumped mass methods for simulation of wave propagation through the soil layers. Since the soil dynamic parameters such as stiffness ($${G/{G_0}}$$) and damping (ξ) are obtained under harmonic loading tests, hence, in the case of equivalent-linear analyses, direct implementation of the maximum strain during seismic excitation overestimates the stiffness reduction and hysteresis attenuation level. Therefore, there is a common correction method which applies a fraction of maximum strain (effective strain) to update $${G/{G_0}}$$ and ξ values from experimental curves. This ratio should reflect the non-linearity level of the model during an earthquake with magnitude M, thus an effective shear strain ratio (SSR) defined as $${\rm{SSR}} = {{( {M - 1} )} /{10}}$$ would be proper (Idriss & Sun 1992). The fully non-linear analyses employ the hyperbolic modified (Kondner & Zelasko 1964) model for the backbone curve formulation (Matasović & Vucetic 1995). The model is extended to pressure-dependent in DEEPSOIL by making the shear strength as a function of the confining pressure (Hashash & Park 2001). Masing's (1926) model which is widely used for constructing cyclic loops is applied in DEEPSOIL for representing the cyclic loading/unloading behaviour of the soil (Hashash et al. 2012). Implementing Masing rules force the cyclic loops to follow the backbone curve but overestimate the damping ratio at large strains (Stewart et al. 2014). Phillips & Hashash (2009) proposed a modification to apply a pinching manner to the cyclic loops of the soil at large strains. This method is implemented in DEEPSOIL (Stewart et al. 2014). Roten et al. (2012) and Hartzell et al. (2002) showed that even under large earthquake excitations, the layers below ∼90 m depth remain elastic. On the other hand, due to the limited depth extent of some profiles of the JICA study, we extend our profiles with a velocity of 600 m s−1 (or 800 m s−1 correspondingly) to a depth of 100 m for these profiles, following Roten et al. (2012) and Hartzell et al. (2002). Also an elastic half-space is defined as the boundary condition beneath the soil column to prevent reflection of waves from bottom and specifically formation of unrealistic resonances at the natural frequencies of soil columns (Roten et al. 2012). We introduced a shear wave velocity of 800–1200 m s−1 for the bedrock to represent the average effect of shallower layers of the deep basement. Ground motions at the bedrock level are obtained from combination of the LF motions at the free surface of the spectral element simulation with the HF motions calculated for the same fault and velocity structure model. In fact, due to its negligible effect on ground motions in comparison with crustal thickness, we have ignored the thickness of the shallow soils (with maximum depth of 50–100 m) in the first step of ground motion calculations (i.e. generation of LF and HF motions) as done by Hartzell et al. (2002). The input ground motions in DEEPSOIL are defined as ‘outcrop’ motions since they are obtained from the free surface of the spectral element model. The ‘outcrop’ motion is the motion that occurs at a free surface and is simply twice the upward propagating wave motion. Using ‘outcrop’ input ground motions is consistent with the implementation of elastic half-space as the boundary condition of soil column. The water table depth in Tehran decreases from north to the south, consistent with the topography and hydraulic gradient of the underground water. At the elevated parts of the city, which are characterized by hard and dense soils and rocks, shallow ground water is not available. However, at the central and southern parts of the city ground water depth changes from 50 to 15 m as directing southwards. Oscillations in pore-water pressure change dynamic soil properties at each time step even at small strains, for example maximum shear modulus and maximum shear stress. Consequently, implementation of the effective stress in calculations is too time consuming against total stress calculations, specifically for studies with large amounts of models. Moreover, seepage analysis parameters such as permeability coefficients are not well constrained for Tehran soils, therefore, following Roten et al. (2012) and Hartzell et al. (2002), we use the total stress approach in calculations. Annual measurements of the ground water level, attained by IWRMC (Iran Water Resource Management Company), are used for total stress calculations in this study (Keshavarz Bakhshayesh 2003). 3 RESULTS AND DISCUSSION We applied 50 s low-frequency time histories to generate broad-band synthetics at the bedrock level and then studied effects of local soil response with linear, equivalent-linear and non-linear soil column models. Linear soil response is compared with non-linear response for each station to evaluate the effect of soil non-linearities on the ground motions. 3.1 Broad-band synthetics at bedrock Low-frequency synthetic ground motions are provided for all of the 42 stations (with 1 km intervals) along the N-S profile obtained from the set of 17 low-frequency simulations in the previous study (Majidinejad et al. 2017). Each of these LF synthetics are combined with seven high-frequency scatterograms to achieve 119 broad-band 0–10 Hz ground motions. Fig. 3 depicts the LF and broad-band velocity and acceleration time histories and amplitude spectra of a typical scenario model for two stations of the footwall side along the N-S profile. One of the stations is located near the fault trace and out of the basin (Rx = + 2 km) and the other is located further above the deep parts of basin (Rx = + 14 km). As it is observable, the amplitudes of the LF model fall at frequencies above 2 Hz, while the broad-band synthetic amplitudes are considerable up to 10 Hz. As it can be seen, the general form of the broad-band time histories are consistent with real observations in the sense that ground motions first begin with HF motions and then continue with arrival of the main LF shear waves (Hartzell et al. 2004). Late surface Rayleigh waves with large amplitudes are visible at the station above the deep Tehran basin. By achieving broad-band synthetics, we are now able to calculate the HF spectral amplitudes and PGA of the expected ground motions in Tehran for an M7.2 earthquake induced by the NTF. As comprehensively discussed in Majidinejad et al. (2017), we assume that the fault-normal component is the predominant component of the ground motions due to the NTF earthquake and contains major portion of the seismic wave's energy taking into account the fault mechanism and 2-D geometry of the fault and basin. Figure 3. View largeDownload slide LF and broad-band ground motion time histories and amplitude spectra at two typical stations. Figure 3. View largeDownload slide LF and broad-band ground motion time histories and amplitude spectra at two typical stations. Fig. 4 shows the comparison of PGA and response spectra values from broad-band simulations at the Bedrock level with the Boore & Atkinson (2008) GMPEs. Results demonstrate very good agreement with the empirical relationships on the average values. Almost all the data points, except the ones corresponding to the deepest parts of the basin, fall within one standard deviation of the empirical relationships. Differences between average values of the synthetic ground motions parameters and empirical values are less than 50 per cent of the empirical standard deviations at all bins. It is noteworthy to point out that the inclusion of the high-frequency contents in broad-band synthetics improves the match between synthetics and empirical relationships, even at low-frequency SAs. Amplification of the SA values due to the deep sedimentary basin of Tehran is also present at high frequencies and HF SAs follow a similar trend to the trend of 1.5–2 Hz SAs in the LF models. This can be attributed to the scaling of the Fourier amplitudes of the HF motions to the LF amplitudes during generation of BB ground motions. However, there are more spatial variations in the HF SA values with respect to the LF SA values due to the stochastic manner of HF contents (Majidinejad et al. 2017). Figure 4. View largeDownload slide Ground motion parameters obtained from broad-band simulation at the geotechnical bedrock level compared with the empirical equation of Boore & Atkinson (2008). Black bars represent for simulation results, including the standard deviation (error bars), binned into logarithmically spaced distance ranges. The thick black line corresponds to the mean values of empirical GMPE; the dashed line show mean ± σ of the empirical GMPE correspondingly. Figure 4. View largeDownload slide Ground motion parameters obtained from broad-band simulation at the geotechnical bedrock level compared with the empirical equation of Boore & Atkinson (2008). Black bars represent for simulation results, including the standard deviation (error bars), binned into logarithmically spaced distance ranges. The thick black line corresponds to the mean values of empirical GMPE; the dashed line show mean ± σ of the empirical GMPE correspondingly. In order to take into account the soil column behaviour, for each LF model, we selected the broad-band time history that produces the nearest responses to the mean values obtained from all correspondent LF models. Therefore, we attained 17 final broad-band models—with dissimilar LF and HF contents—which are then used for investigation of the local site effects at 30 stations situated in the municipal regions of Tehran on the N-S profile (i.e. a total number of 510 simulations for each of the linear and non-linear soil columns analyses). 3.2 Linear results Soil column models are set such that minimum resolvable frequency of 30 Hz is available for all models, with layers defined at 2 m interval depths. Fig. 5 compares the average 3 and 5 Hz SAs, obtained by linear amplification of the broad-band synthetic ground motions, with the empirical GMPE of Boore & Atkinson (2008) in the N-S profile of city at the free-surface level. The required parameters for the empirical relationships, such as Vs, 30, Rx, etc., are evaluated at each site with regard to the values in the models. Although the average values of the broad-band synthetics parameters at the bedrock level appeared to be within one-half standard deviation of the empirical relationships, linearly amplified values of these parameters at the surface do not coincide with the empirical relationships. The mean values differ more than one standard deviation from the empirical observations at many parts of the N-S profile, in particular over the deep parts of the basin. This inconsistency should be attributed to the insufficient accuracy of the linear modeling in representing soil behaviour during severe ground motions. While central and southern parts of the city experience considerable amplifications because of shallow soil layers, the northern parts are not affected significantly by local site effects, due to the shallow bedrock depth on the hangingwall side (Fig. 5). Figure 5. View largeDownload slide Average spectral accelerations at the ground surface level of the stations along the N-S profile obtained by linear, equivalent-linear and fully non-linear analyses. Empirical values for each station on the profile is also calculated based on site specific values and presented for comparison. Figure 5. View largeDownload slide Average spectral accelerations at the ground surface level of the stations along the N-S profile obtained by linear, equivalent-linear and fully non-linear analyses. Empirical values for each station on the profile is also calculated based on site specific values and presented for comparison. Maximum linear SA values occur on the footwall side of the fault and reach up to ∼2.4 and ∼3.7 g on average at 3 and 5 Hz, respectively. On the hangingwall side, the maximum linear SAs at these frequencies are about ∼1.1 and ∼1.6 g which occur in the vicinity of the fault. A notable remark is that, in agreement with Dalguer & Mai (2012) and in contrast with empirical relationships, the maximum values of responses do not occur just in the vicinity of the fault, but in a distance 2–3 km from its trace. Also, it can be seen that at stations very close to the fault, results are lower than empirical predictions. 3.3 Non-linear results Many energy dissipation mechanisms occur at the shallow layers of soil during large earthquake ground motions, therefore, a non-linear model of the soil which takes into account soil's stiffness degradation and hysteresis damping may be capable of reducing the unrealistic high values of the linear amplification. In order to evaluate the soil non-linear effects on the HF ground motions, we performed both full non-linear and equivalent-linear analyses. The reference strain and damping ratios are defined at 2 m intervals of the vertical profiles. We also assumed that the soil layers with shear wave velocities above 750 m s−1 do not undergo non-linear strains, due to their stiff and dense soil compositions. Linear, non-linear and equivalent-linear results of a scenario earthquake are compared in Fig. 6 for a station at distance of 13 km from the fault trace, on the footwall side and above the basin. As it is clear in Fig. 6, linear soil response amplifies the ground motions significantly which leads to large apexes in time history, while, these peaks are lowered and also smoothed by the non-linear behaviour of the soil column. It is noteworthy that the reductions due to the soil non-linearities are more pronounced at the high amplitude peaks, which are followed by large inelastic strains. Comparison of the non-linear and equivalent-linear results shows that although equivalent-linear method predicts the deamplifications at large peaks properly, it also reduces the modest and small apexes of the time history which results in underestimation of medium peaks and a somewhat smoothed time history. This observation clearly associates with the assumption of the equivalent-linear method, which utilizes the stiffness and damping ratio calculated for the maximum peak at all time steps. In this case, the equivalent-linear method has underpredicted the spectral accelerations with respect to the fully non-linear analysis at frequencies above 3 Hz (Fig. 6, right). Figure 6. View largeDownload slide Effect of shallow soil layers on broad-band ground motions of a station located 13 km from fault's trace, on the footwall side. Left: time histories and right: response spectral accelerations. Figure 6. View largeDownload slide Effect of shallow soil layers on broad-band ground motions of a station located 13 km from fault's trace, on the footwall side. Left: time histories and right: response spectral accelerations. Comparison between linear and non-linear soil responses reveals considerable decrease in spectral acceleration values of the HF content and also increase in soil column resonant periods due to the reduced stiffness and elevated damping of the non-linear soil layers (Fig. 6, right). As stated by Hartzell et al. (2004) and also reported by Roten et al. (2012), the opposing effects of stiffness reduction and increase of damping in non-linear soil response, determine further amplification or deamplification of the responses with respect to the linear behaviour. Since damping effect is predominant at high frequencies (>2 Hz), the non-linear responses at these frequencies are usually below the linear values. On the other hand, at low frequencies (<1 Hz) damping effects are minor, therefore, non-linearities lead to higher amplifications at this range (Fig. 6, right). Fig. 5 also displays the average spectral accelerations in the case of non-linear soil response for comparison with ground motion prediction equations. Intense SA values above the basin are diminished in this case and average values of numerical results are now limited to the mean ± one standard deviation of the GMPEs at many stations. The maximum non-linear deamplification in response takes place at a station 18.5 km from the fault, on the footwall side. At this site, the non-linear deamplification ratio (i.e. $$1 - {{{\rm{S}}{{\rm{A}}_{{\rm{NL}}}}} / {{\rm{S}}{{\rm{A}}_{{\rm{Lin}}}}}}$$) reaches to 33, 60 and 63 per cent for frequencies of 2, 3 and 5 Hz, respectively. On the footwall side of the NTF, in a distance range of 0–8 km from the fault trace where the soil profiles depth to the geotechnical bedrock is not very high, due to the hysteresis behaviour in few thin soil layers, dissipation of energy is negligible. Hence, these stations experience low deamplifications or even higher amplifications as a result of stiffness reductions in the non-linear response (particularly at low frequencies). At distances within 8–22 km from the fault trace, depth of low-velocity layers to the geotechnical bedrock increases and consequently, at the same time with the reduction of dynamic stiffness during non-linear soil behaviour, the amount of hysteresis damping increases significantly. This highly elevated dissipation in non-linear response leads to considerable reduction of responses with respect to the linear behaviour. Although stations further than 22 km are located on weak soil layers, low bedrock motions at these stations are unable to develop significant plasticity in soil columns, therefore, these sites exhibit lower non-linear deamplifications. Figs 7(a) and (b) show the vertical profiles of maximum accelerations, averaged over scenario earthquakes, at stations located 9.5 and 18.5 km from the fault trace and for the case of linear and non-linear analyses. Both station are located on the footwall side and above the deep sedimentary basin (Fig. 1) but one has medium depth, low-velocity layers with average shear wave velocity of Vs, 30 = 400 m s − 1  at the top of the geotechnical bedrock (Rx = 9.5 km) and the other is located above thick very low-velocity layers, Vs, 30 = 220 m s − 1, down to the geotechnical bedrock (Rx = 18.5 km). Figs 7(c) and (d) depict the hysteresis loops of soil at depth of 30 m for these stations, which clearly show the extent of plasticity in soil columns. Limited plasticity at the station out of the basin explains similar linear and non-linear responses. Figs 7(e) and (f) are hysteresis loops at depth of 90 m and confirm the presumption of very limited plasticity at depths below 100 m. Figure 7. View largeDownload slide Effect of soil's plasticity development on deamplification caused by non-linear soil response at two representative stations. (a) and (b) vertical profiles of average PGA in linear and non-linear response. (c) and (d) hysteresis loops of the soil layer at depth of 30 m. (c) and (d) hysteresis loops at depth of 90 m. Figure 7. View largeDownload slide Effect of soil's plasticity development on deamplification caused by non-linear soil response at two representative stations. (a) and (b) vertical profiles of average PGA in linear and non-linear response. (c) and (d) hysteresis loops of the soil layer at depth of 30 m. (c) and (d) hysteresis loops at depth of 90 m. 3.4 Non-linear deamplification models Although we investigated just one N-S profile along Tehran, by considering JICA (2000) soil profiles and taking into account depth of the geotechnical bedrock, geological units and average shear wave velocity, it appears plausible to suggest three major site classes I, II and III for the Tehran plain as in Fig. 8. Sites located in each site class are roughly assumed to represent similar responses to the same bedrock excitations. Class I includes shallow bedrock profiles (depth <30 m) and coarse-grained materials with approximately high shear wave velocities (Vs, 30 above 450 m s−1). Class II consists of mixture of sand and clay and the bedrock depth at this class changes from 35 to 60 m. Average shear wave velocity (Vs, 30) at this class ranges between 300 and 450 m s−1. Finally, class III is characterized by deep bedrock (75m< depth) and fine clayey materials with low shear wave velocities (Vs, 30 < 300 m s−1). Figure 8. View largeDownload slide Map of the site classes proposed based on behaviour of soil profiles in JICA (2000). Thick black lines represent the municipal district's boundaries of the Tehran. Figure 8. View largeDownload slide Map of the site classes proposed based on behaviour of soil profiles in JICA (2000). Thick black lines represent the municipal district's boundaries of the Tehran. At this stage, we propose deamplification factors for the Tehran site classes, as done by Roten et al. (2012) for the Salt Lake City. Such results can be used for stations all over the city in order to calculate the non-linear soil responses with respect to the corresponding linear spectral accelerations. Since we aim to introduce models for any possible earthquakes in the city and regardless of the causative fault—in order to remove the source to site distance effects on input motions—we apply all the synthetic ground motions (510 time histories) to each of archetype soil profile models (10 models). Therefore, total number of 10 200 simulations are performed and results are categorized into three major site classes. Linear and non-linear simulation results are analysed to obtain statistical relations between two (linear and non-linear) sets at each class. A possible caveat here is the duration of ground motions for earthquakes with magnitudes smaller (or larger) than M7.2. Although the deamplification functions correct for distance, they do not account for shorter (or longer) duration that would be characteristic of smaller (larger) earthquakes, but which would result in different non-linear site response. Fig. 9 shows the non-linear versus linear spectral acceleration scatters for three site classes, considering all frequencies at the range of interest (top) and frequency of 5 Hz (bottom). As it can be obtained from Fig. 9, the maximum linear 5 Hz spectral accelerations reach up to 75.61, 67.56 and 50.81 m s−2 for site classes I, II and III, while the non-linear responses are limited to 31.38, 20.86 and 8.37 m s−2 correspondingly. It is notable that non-linear soil effects become evident for linear spectral accelerations larger than ∼2 m s−2. Also the dispersion of data points increases at larger SA values due to the more plasticity development and non-linearities. As anticipated site, classes I and III show the lowest and highest deamplifications to the normalized source-to-site distance ground motions, respectively, however, this conclusion may differ slightly for the earthquake taking place on the NTF as the stiffest site class I is the nearest and the loosest site class III is the furthest site class to the NTF. Figure 9. View largeDownload slide Scatterogram of non-linear SAs versus corresponding linear SAs for three site class of Tehran plain. Top: spectral accelerations at all frequencies. Bottom: spectral accelerations at frequency of 5 Hz. Black lines correspond to the correction functions calculated in the form of eq. (1) and dashed lines correspond to their standard deviations, respectively. Figure 9. View largeDownload slide Scatterogram of non-linear SAs versus corresponding linear SAs for three site class of Tehran plain. Top: spectral accelerations at all frequencies. Bottom: spectral accelerations at frequency of 5 Hz. Black lines correspond to the correction functions calculated in the form of eq. (1) and dashed lines correspond to their standard deviations, respectively. The overall trends of scattered data in Fig. 9 suggest a second-order relationship between the logarithm values of the linear and non-linear spectral accelerations, therefore, we use a relationship in the form of eq. (2) to represent the non-linear response as a function of linear spectral acceleration (Roten et al. 2012):   $$\log ({\rm{S}}{{\rm{A}}_{{\rm{NL}}}}) = b * \log ({\rm{S}}{{\rm{A}}_{{\rm{Lin}}}}) + a * {(\log ({\rm{S}}{{\rm{A}}_{{\rm{Lin}}}}))^2}.$$ (2) Table 2 lists the values and standard deviations calculated for regression parameters of eq. (2) for spectral accelerations at various frequencies in the range of interest and also PGA. Values in Table 2 can be interpreted by characteristics of the site classes. Site class I, which shows negligible deamplifications, is approximated by close to unity ‘b’ values and small values of ‘a’ parameter. Considerable deamplifications in site classes II and III, which are composed of soft materials, are reflected in high absolute values of ‘a’ parameter. Although at some frequencies, ‘a’ values of site class II are higher than corresponding values in site class III, parameter ‘b’ at these frequencies is meaningfully smaller for site class III. Small values of ‘b’ parameter imply premature initiation of non-linearities in site class III. Table 2. Parameters of regression equation (2) and their standard deviation.   Site class I  Site class II  Site class III  Frequency (Hz)  b  a  std-b  std-a  b  a  std-b  std-a  b  a  std-b  std-a  10  1.19  −0.16  0.248  0.135  1.07  −0.19  0.145  0.068  1.29  −0.52  0.313  0.216  9.5  1.12  −0.08  0.146  0.064  1.05  −0.17  0.170  0.087  0.87  −0.15  0.204  0.123  9  1.29  −0.21  0.173  0.091  1.05  −0.20  0.174  0.084  0.73  −0.11  0.126  0.057  8.5  1.18  −0.16  0.121  0.046  1.07  −0.20  0.192  0.093  0.73  −0.12  0.153  0.085  8  1.07  −0.14  0.146  0.052  1.14  −0.24  0.164  0.080  0.77  −0.13  0.182  0.107  7.5  0.96  −0.08  0.143  0.063  1.16  −0.27  0.163  0.080  0.79  −0.16  0.170  0.104  7  0.98  −0.10  0.122  0.057  1.07  −0.24  0.182  0.093  0.66  −0.09  0.184  0.106  6.5  1.01  −0.11  0.123  0.060  1.05  −0.23  0.200  0.105  0.66  −0.07  0.204  0.124  6  1.05  −0.15  0.159  0.072  1.02  −0.21  0.164  0.084  0.69  −0.10  0.130  0.065  5.5  1.02  −0.13  0.221  0.110  1.03  −0.22  0.158  0.078  0.61  −0.03  0.140  0.078  5  0.98  −0.11  0.210  0.101  1.00  −0.21  0.159  0.080  0.62  −0.06  0.138  0.077  4.5  1.04  −0.16  0.184  0.091  0.97  −0.19  0.151  0.076  0.63  −0.06  0.137  0.071  4  1.05  −0.15  0.164  0.087  0.98  −0.20  0.150  0.077  0.67  −0.08  0.149  0.085  3.5  1.02  −0.14  0.145  0.077  0.96  −0.17  0.160  0.088  0.69  −0.09  0.162  0.104  3  1.05  −0.15  0.149  0.073  0.94  −0.16  0.142  0.074  0.72  −0.12  0.157  0.097  2.5  1.00  −0.13  0.132  0.056  0.97  −0.17  0.158  0.083  0.80  −0.19  0.171  0.107  2  1.02  −0.14  0.160  0.076  0.98  −0.17  0.172  0.094  0.85  −0.22  0.155  0.091  All  0.95  −0.06  0.025  0.030  0.90  −0.10  0.024  0.035  0.80  −0.15  0.023  0.047  PGA  1.02  −0.12  0.178  0.082  0.91  −0.15  0.184  0.095  0.56  −0.03  0.170  0.119    Site class I  Site class II  Site class III  Frequency (Hz)  b  a  std-b  std-a  b  a  std-b  std-a  b  a  std-b  std-a  10  1.19  −0.16  0.248  0.135  1.07  −0.19  0.145  0.068  1.29  −0.52  0.313  0.216  9.5  1.12  −0.08  0.146  0.064  1.05  −0.17  0.170  0.087  0.87  −0.15  0.204  0.123  9  1.29  −0.21  0.173  0.091  1.05  −0.20  0.174  0.084  0.73  −0.11  0.126  0.057  8.5  1.18  −0.16  0.121  0.046  1.07  −0.20  0.192  0.093  0.73  −0.12  0.153  0.085  8  1.07  −0.14  0.146  0.052  1.14  −0.24  0.164  0.080  0.77  −0.13  0.182  0.107  7.5  0.96  −0.08  0.143  0.063  1.16  −0.27  0.163  0.080  0.79  −0.16  0.170  0.104  7  0.98  −0.10  0.122  0.057  1.07  −0.24  0.182  0.093  0.66  −0.09  0.184  0.106  6.5  1.01  −0.11  0.123  0.060  1.05  −0.23  0.200  0.105  0.66  −0.07  0.204  0.124  6  1.05  −0.15  0.159  0.072  1.02  −0.21  0.164  0.084  0.69  −0.10  0.130  0.065  5.5  1.02  −0.13  0.221  0.110  1.03  −0.22  0.158  0.078  0.61  −0.03  0.140  0.078  5  0.98  −0.11  0.210  0.101  1.00  −0.21  0.159  0.080  0.62  −0.06  0.138  0.077  4.5  1.04  −0.16  0.184  0.091  0.97  −0.19  0.151  0.076  0.63  −0.06  0.137  0.071  4  1.05  −0.15  0.164  0.087  0.98  −0.20  0.150  0.077  0.67  −0.08  0.149  0.085  3.5  1.02  −0.14  0.145  0.077  0.96  −0.17  0.160  0.088  0.69  −0.09  0.162  0.104  3  1.05  −0.15  0.149  0.073  0.94  −0.16  0.142  0.074  0.72  −0.12  0.157  0.097  2.5  1.00  −0.13  0.132  0.056  0.97  −0.17  0.158  0.083  0.80  −0.19  0.171  0.107  2  1.02  −0.14  0.160  0.076  0.98  −0.17  0.172  0.094  0.85  −0.22  0.155  0.091  All  0.95  −0.06  0.025  0.030  0.90  −0.10  0.024  0.035  0.80  −0.15  0.023  0.047  PGA  1.02  −0.12  0.178  0.082  0.91  −0.15  0.184  0.095  0.56  −0.03  0.170  0.119  View Large We also compared the amplification ratios from fully non-linear analyses with most recent empirical and semi-empirical relationships proposed by Seyhan & Stewart (2014), Sandikkaya et al. (2013) and Kamai et al. (2014), which are implemented in updated NGA2-west GMPEs. Kamai et al. (2014) and Seyhan & Stewart (2014) equations are obtained by using a data set composed of realistic data and simulations results. In order to estimate the non-linear soil response, they have performed equivalent-linear 1-D simulations. Sandikkaya et al.’s (2013) model is completely based on empirical data. The method used for normalizing the Kamai et al.’s (2014) amplification ratios is the same as the one described by Sandikkaya et al. (2013), such that the values obtained for different reference shear wave velocities become comparable. Simulations show more non-linearity in comparison with relationships, therefore, the amplification ratios for higher base ground motions drop faster than these predictions. Dependency of the non-linear soil effects on the bedrock motions is more notable at the 5 Hz frequency. As it is observable in Fig. 10, the level of amplification ratios for the Vs, 30 = 800 m s − 1 is almost constant with a very slight positive slope at large ground motions which can be due to the reduction of soil stiffness at small plastic strains, while non-linear damping is still negligible. Profiles with Vs, 30 = 480 m s − 1 show a somewhat constant amplification ratio at the frequency of 3 Hz and at small to moderate excitations of the bedrock at 5 Hz. Amplification ratios of profiles with low average velocity Vs, 30 < 300 m s − 1 exhibit a descending manner at both frequencies and their independency from the bedrock motion occurs at levels below ∼0.2 g. Results show that for large ground motions, the highest 3 Hz amplifications occur at profiles with medium shear wave velocity (i.e. Vs, 30 ∼ 500 m s − 1), while this happens at low shear wave velocities (i.e. Vs, 30 ∼ 300 m s − 1) for lower motions. The amplification ratio at 5 Hz is approximately the same in medium and loose profiles. The amplification ratios are generally higher than their corresponding values at 5 Hz which agree with the results obtained by Kamai et al. (2014). Our results are more consistent with the Sandikkaya et al.’s (2013) model, however, the general form of all studies match with each other. The agreement between the results and Kamai et al.’s (2014) model is also acceptable for almost all cases, except for the medium shear wave profiles at 5 Hz. The difference between the analyses results and relationships may originate from various reasons such as different non-linear curves implemented for prediction of soil non-linear behaviour, the difference between non-linear and equivalent-linear results and also assumptions made for reference shear wave velocities. Figure 10. View largeDownload slide Comparison of the non-linear amplification ratios with (semi-) empirical results of Kamai et al. (2014), Seyhan & Stewart (2014) and Sandikkaya et al. (2013). Figure 10. View largeDownload slide Comparison of the non-linear amplification ratios with (semi-) empirical results of Kamai et al. (2014), Seyhan & Stewart (2014) and Sandikkaya et al. (2013). 4 CONCLUSIONS In this study, we generated broad-band time histories by combining low-frequency synthetic motions, obtained by dynamic simulations, with stochastic physics based high-frequency motions calculated from wave scattering theory. Since the central and southern parts of the city of Tehran are made up of many old low-rise masonry buildings with high natural frequencies, providing such broad-band ground motions is an essential step to make deterministic low-frequency ground motions applicable to practical engineering tasks for the city. Results obtained from simulations at the bedrock level show good consistency with the empirical values in the sense that average values of ground motion parameters fall within a half standard deviation from mean values of the empirical relationships. As discussed by Hartzell et al. (2002) and Frankel & Stephenson (2000) existence of a deep basin in the vicinity of a reverse fault causes large fault-normal and vertical motions. Large amplitudes are associated with upward shear waves in the basin and formation of surface waves at the edges. In the next step, we investigated local site effects on ground motions by taking into account the soil profiles proposed for the Tehran plain by previous studies. Both linear and non-linear soil responses were calculated for stations along the N-S profile of the city. Hartzell et al. (2002) compared the effect of variations in source models and randomized soil profile parameters by considering 18 input synthetic ground motions (i.e. at the bedrock level) and 20 random soil profiles. They concluded that variations in source models are predominant on variability of the ground motion parameters. Therefore, considering 17 broad-band models for each station in our study may be appropriate for considering slight to moderate changes in soil profile characteristics (Toro 1997). Linear soil response models predicted large unrealistic amplification factors, ranging from 1.5 to more than 3, which do not seem realistic with respect to the empirical relations and observations (Boore & Joyner 1997). This large amplifications resulted in overestimation of the ground motion parameters with respect to empirical relationships at the surface. In contrast, taking into account the non-linear soil behaviour led to plausible results which follow the empirical relationships with one standard deviation. Average values of the horizontal PGA at the surface reach up to 1.09 g on the hangingwall side. On the footwall side, which contains most parts of the city, mean PGA value is about 0.93 g at the northern, 0.59 g at the central and 0.40 g at the southern parts. According to the results, the extent of plasticity development in soil column plays a major role in deamplification due to non-linear soil response. The development of plasticity is not only controlled by the soil properties but also by the ground motion intensity at the base of the column, therefore, soft soil columns with large linear ground motion response experience the highest deamplifications due to the non-linearities. This phenomenon explains why stations very close to the fault show lower amplification factors. It also explains why in some previous earthquakes such as Northridge 1994, main shocks were characterized by lower amplifications with respect to their aftershocks (Field et al. 1997). As a conclusion, since during moderate and small earthquakes soil columns remain elastic, high linear amplifications are not improbable and therefore in the purpose of performing microzonation studies, both the linear and non-linear behaviour of soil profiles should be considered in estimation of the local site effects. The spatial distribution of amplification factors in linear and non-linear cases are quite different, specifically at higher frequencies. Stations with the highest amplifications in linear case, experience the highest deamplifications due to non-linear soil response. Conclusions obtained here place a question on the values reported by Haghshenas (2005) for amplifications of the Tehran alluvium. They reported very large, up to 7, amplification ratios for stations above the Tehran basin. Unrealistically large values in that study are calculated based on microearthquake investigations, hence, do not take into account the soil non-linearities at all. However, it is interesting that our linear simulations also predict overall amplification factors (amplifications caused by deep sedimentary basin up to the bedrock multiplied by near-surface amplifications from bedrock to surface) up to 7.5–8.5 at deepest parts of the basin. Finally, it is notable that although equivalent-linear results underestimate the modest peaks of the motions, they can properly estimate the maximum and large peaks of the motion with very effective reduction in computational time cost. Results demonstrate the need for more accurate studies on properties of the soil profiles of the city. We assigned soil non-linear parameters with respect to the alluvium's structure, composition and loading history. In this regard, more boreholes with larger depths should be investigated over the city to provide a complete database for microzonation studies. Until that time, we proposed a model for three site classes in Tehran to estimate the non-linear soil column response from linear values. In order to propose a more applicable and comprehensive model, we eliminated source-to-site distance effect by considering all time histories for all soil profiles. The correction equation and corresponding parameters for frequencies in the range of 2–10 are proposed in this paper. Broad-band motions calculated here, can be used for estimation of damages and required seismic rehabilitation plans in the city due to the probable M7.2 earthquake of the NTF. The final goal of this research is to propose SA maps that include soil effects for the city of Tehran and also providing a specific synthetic-based GMPE for the Tehran region by using simulation results which will be presented in a forthcoming paper. 4.1 Data and resources The code ‘Broad-Band toolbox’ is developed by Martin Mai and Walter Imperatori at ETH Zurich and can be accessed from (https://ces.kaust.edu.sa/Pages/Software.aspx, last accessed 2017 January). Attenuation curves were plotted by OpenSHA attenuation relationship plotter (http://www.opensha.org, last accessed 2017 January). Non-linear soil analyses are performed by the academic software DEEPSOIL, developed by Youssef Hashash at University of Illinois at Urbana-Champaign and is accessible from http://deepsoil.cee.illinois.edu/. Acknowledgements The authors are grateful to Martin Mai for his kindly providing ‘Broad-Band toolbox’ software. We are also very grateful to Dr Eiichi Fukuyama and two anonymous reviewers for their insightful and constructive comments, which significantly improved the manuscript. H. Zafarani thanks the continuing support of the International Institute of Earthquake Engineering and Seismology during this research, in the framework of the Probabilities of Earthquake Ruptures in Iran (PERSIA) project. REFERENCES Amir-Faryar B., Aggour M.S., 2016. Effect of fibre inclusion on dynamic properties of clay, Geomech. Geoeng. , 11( 2), 104– 113. https://doi.org/10.1080/17486025.2015.1029013 Google Scholar CrossRef Search ADS   Anderson J.G., 1997. Benefits of scenario ground motion maps, Eng. Geol. , 48( 1–2), 43– 57. https://doi.org/10.1016/S0013-7952(97)81913-8 Google Scholar CrossRef Search ADS   Bay J.A., Sasanakul I., 2005. Measurement of nonlinear dynamic properties of bonneville silty-clays: NIW program, Department of Civil and Environmental Engineering, Utah State University, Logan, Utah, Award 04HQGR0055. Berberian M., Yeats R.S., 1999. Patterns of historical earthquake rupture in the Iranian Plateau, Bull. seism. Soc. Am. , 89, 120– 139. Boore D.M., 1983. Stochastic simulation of high-frequency ground motions based on seismological models of the radiated spectra, Bull. seism. Soc. Am. , 73, 1865– 1894. Boore D.M., 2003. Simulation of ground motion using the stochastic method, Pure appl. geophys. , 160( 3), 635– 676. https://doi.org/10.1007/PL00012553 Google Scholar CrossRef Search ADS   Boore D.M., Atkinson G.M., 2008. Ground-motion prediction equations for the average horizontal component of PGA, PGV, and 5%-samped PSA at spectral periods between 0.01 s and 10.0 s, Earthq. Spectra , 24( 1), 99– 138. https://doi.org/10.1193/1.2830434 Google Scholar CrossRef Search ADS   Boore D.M., Joyner W.B., 1997. Site amplifications for generic rock sites, Bull. seism. Soc. Am. , 87, 327– 341. Dalguer L., Mai P., 2012. Prediction of near-source ground motion exceeding 1 g at low frequencies (<2 Hz) from Mw ∼ 6.5 deterministic physics-based dynamic rupture simulations, in Proceedings of the 15th World Conference on Earthquake Engineering (15WCEE), Lisbon, Portugal. Darendeli M.B., 2001. Development of a New Family of Normalized Modulus Reduction and Material Damping Curves, PhD thesis , University of Texas. Dreger D., Tinti E., Cirella A., 2007. In Seismological Society of America 2007 Annual Meeting , 2007 April 11–13, Waikoloa, Hawaii. Emami M.H., Amini B., Jamshidi K.H., Afsharyanzadeh A.M., 1993. Geology Map of Tehran (1:100,000 scale) , Geological Survey of Iran. Farahani J.V., Zaré M., Cichowicz A., 2012. Attenuation of high-frequency P and S waves in south and southeast Tehran using blast data, Soil Dyn. Earthq. Eng. , 40, 99– 108. https://doi.org/10.1016/j.soildyn.2012.03.005 Google Scholar CrossRef Search ADS   Field E.H., Johnson P.A., Beresnev I.A., Zeng Y., 1997. Nonlinear ground-motion amplification by sediments during the 1994 Northridge earthquake, Nature , 390( 6660), 599– 602. https://doi.org/10.1038/37586 Google Scholar CrossRef Search ADS   Frankel A., Stephenson W., 2000. Three-dimensional simulations of ground motions in the Seattle region for earthquakes in the Seattle fault zone, Bull. seism. Soc. Am. , 90( 5), 1251– 1267. https://doi.org/10.1785/0119990159 Google Scholar CrossRef Search ADS   Ghorbani A., Kalantari F., Zohori M., 2014. Estimation of the drained shear strength parameters of fine grained soil of southern Tehran using nonlinear regression model, J. Eng. Geol. , 7, 1869– 1900. Haghshenas E., 2005. Conditions Géotechniques et aléa Sismique local à Téhéran 2005, PhD thesis, Joseph Fourier University (Grenoble I), Grenoble France. Hartzell S., Bonilla L., Williams R.A., 2004. Prediction of nonlinear soil effects, Bull. seism. Soc. Am. , 94( 5), 1609– 1629. https://doi.org/10.1785/012003256 Google Scholar CrossRef Search ADS   Hartzell S., Leeds A., Frankel A., Williams R.A., Odum J., Stephenson W., Silva W., 2002. Simulation of broadband ground motion including nonlinear soil effects for a magnitude 6.5 earthquake on the Seattle fault, Seattle, Washington, Bull. seism. Soc. Am. , 92( 2), 831– 853. https://doi.org/10.1785/0120010114 Google Scholar CrossRef Search ADS   Hashash Y., Groholski D., Phillips C., Park D., Musgrove M., 2012. DEEPSOIL 5.1., User Manual and Tutorial , p. 107, University of Illinois, Urbana, IL, USA. Hashash Y.M., Park D., 2001. Non-linear one-dimensional seismic ground motion propagation in the Mississippi embayment, Eng. Geol. , 62( 1–3), 185– 206. https://doi.org/10.1016/S0013-7952(01)00061-8 Google Scholar CrossRef Search ADS   Hassani B., Zafarani H., Farjoodi J., Ansari A., 2011. Estimation of site amplification, attenuation and source spectra of S-waves in the East-Central Iran, Soil Dyn. Earthq. Eng. , 31( 10), 1397– 1413. https://doi.org/10.1016/j.soildyn.2011.05.017 Google Scholar CrossRef Search ADS   Idriss I., Sun J.I., 1992. SHAKE91: A Computer Program for Conducting Equivalent Linear Seismic Response Analyses of Horizontally Layered Soil Deposits , Center for Geotechnical Modeling, Department of Civil and Environmental Engineering, University of California, Davis, CA. Irikura K., 1986. Prediction of strong acceleration motion using empirical Green's function, in Proc. 7th Japan Earthq. Eng. Symp , pp. 151– 156, IEEE. Jafari M., Razmkhah A., Sohrabi A., Keshavarz M., Pourazin K., 2001. Complementary seismic microzonation studies for South of Tehran, International Institute of Earthquake Engineering and Seismology, Tehran, Iran , p. 126. Jafari M., Saniei H., Haidari M., 2000. Structural classification of clay soils and its application in classifying Tehran City clays, Geosciences , 9, 2– 17. Jafari M.K., Shafiee A., Razmkhah A., 2002. Dynamic properties of fine grained soils in south of Tehran, J. Seismol. Earthq. Eng. , 4, 25– 35. Jica C., 2000. The study on seismic microzoning of the Greater Tehran Area in the Islamic Republic of Iran, Pacific Consultants Int. Rep., OYO Cooperation, Japan . Kamai R., Abrahamson N., Silva W.J., 2014. Nonlinear horizontal site amplification for constraining the NGA-West2 GMPEs, Earthq. Spectra , 30( 3), 1223– 1240. https://doi.org/10.1193/070113EQS187M Google Scholar CrossRef Search ADS   Keshavarz Bakhshayesh M., 2003. The effect of southern Tehran faults on alluvial deposits, in Proceeding 4th International Conference on Sismology and Earthquake Engineering (SEE 4), Tehran 12–14. Kondner R.L., Zelasko J.S., 1964. Void ratio effects on the hyperbolic stress-strain response of a sand, in Laboratory Shear Testing of Soils , ASTM International. Google Scholar CrossRef Search ADS   Liu P., Archuleta R.A., Hartzell S.H., 2006. Prediction of broadband ground-motion time histories: hybrid low/high- frequency method with correlated random source parameters, Bull. seism. Soc. Am. , 96( 6), 2118– 2130. https://doi.org/10.1785/0120060036 Google Scholar CrossRef Search ADS   Mai P.M., Beroza G., 2003. A hybrid method for calculating near-source, broadband seismograms: application to strong motion prediction, Phys. Earth planet. Inter. , 137( 1–4), 183– 199. https://doi.org/10.1016/S0031-9201(03)00014-1 Mai P.M., Imperatori W., Olsen K.B., 2010. Hybrid broadband ground-motion simulations: combining long-period deterministic synthetics with high-frequency multiple S-to-S backscattering, Bull. seism. Soc. Am. , 100( 5A), 2124– 2142. https://doi.org/10.1785/0120080194 Google Scholar CrossRef Search ADS   Majidinejad A., Zafarani H., Vahdani S., 2017. Dynamic simulation of ground motions from scenario earthquakes on the North Tehran Fault, Geophys. J. Int. , 209( 1), 434– 452. https://doi.org/10.1093/gji/ggx017 Google Scholar CrossRef Search ADS   Masing G., 1926. Eigenspannungen und verfestigung beim messing, in Proceedings of the 2nd International Congress of Applied Mechanics , pp. 332– 335, Zürich. Matasović N., Vucetic M., 1995. Generalized cyclic-degradation-pore-pressure generation model for clays, J. Geotech. Eng. , 121( 1), 33– 42. https://doi.org/10.1061/(ASCE)0733-9410(1995)121:1(33) Google Scholar CrossRef Search ADS   Mena B., Mai P.M., Olsen K.B., Purvance M.D., Brune J.N., 2010. Hybrid broadband ground-motion simulation using scattering green's functions: application to large-magnitude events, Bull. seism. Soc. Am. , 100( 5A), 2143– 2162. https://doi.org/10.1785/0120080318 Google Scholar CrossRef Search ADS   Motazedian D., Atkinson G.M. 2005. Stochastic finite-fault modeling based on a dynamic corner frequency, Bull. seism. Soc. Am. , 95( 3), 995– 1010. https://doi.org/10.1785/0120030207 Google Scholar CrossRef Search ADS   Pedrami M., 1981. Pasadenian orogeny and geology of last 700,000 years of Iran, Int. Rep., Geological Survey of Iran (in Persian) , 273. Phillips C., Hashash Y.M., 2009. Damping formulation for nonlinear 1D site response analyses, Soil Dyn. Earthq. Eng. , 29( 7), 1143– 1158. https://doi.org/10.1016/j.soildyn.2009.01.004 Google Scholar CrossRef Search ADS   Rieben H., 1955. The geology of the Teheran plain, Am. J. Sci. , 253( 11), 617– 639. https://doi.org/10.2475/ajs.253.11.617 Google Scholar CrossRef Search ADS   Ritz J.F. et al.  , 2012. Paleoearthquakes of the past 30,000 years along the North Tehran Fault (Iran), J. geophys. Res. , 117, doi:10.1029/2012JB009147. Roten D., Olsen K., Pechmann J., 2012. 3D simulations of M 7 Earthquakes on the wasatch fault, utah, Part II: broadband (0–10 Hz) ground motions and nonlinear soil behavior, Bull. seism. Soc. Am. , 102( 5), 2008– 2030. https://doi.org/10.1785/0120110286 Google Scholar CrossRef Search ADS   Sandikkaya M.A., Akkar S., Bard P., 2013. A nonlinear site-amplification model for the next pan-european ground-motion prediction equations, Bull. seism. Soc. Am. , 103( 1), 19– 32. https://doi.org/10.1785/0120120008 Google Scholar CrossRef Search ADS   Seed H.B., Idriss I.M., 1969. Influence of soil conditions on ground motions during earthquakes, J. Soil Mech. Found. Div. , 95, 99– 138. Seyhan E., Stewart J.P., 2014. Semi-empirical nonlinear site amplification from NGA-West2 data and simulations, Earthq. Spectra , 30( 3), 1241– 1256. https://doi.org/10.1193/063013EQS181M Google Scholar CrossRef Search ADS   Shafiee A., Azadi A., 2007. Shear-wave velocity characteristics of geological units throughout Tehran City, Iran, J. Asian Earth Sci. , 29( 1), 105– 115. https://doi.org/10.1016/j.jseaes.2006.02.005 Google Scholar CrossRef Search ADS   Standard A., 2003. D3999-91. ‘‘Standard test methods for the determination of the modulus and damping properties of soils using the cyclic triaxial apparatus’’, in Annual Book of ASTM Standard . Stewart J.P., Afshari K., Hashash Y.M., 2014. Guidelines for performing hazard-consistent one-dimensional ground response analysis for ground motion prediction, in PEER Rep. , p. 16. Thráinsson H., Kiremidjian A.S., Winterstein S.R., 2000. Modeling of Earthquake Ground Motion in the Frequency Domain , John A. Blume Earthquake Engineering Center. Tiniti E., Fukuyama E., Piatanesi A., Cocco M., 2005. A kinematic source-time function compatible with earthquake dynamics, Bull. seism. Soc. Am. , 95( 4), 1211– 1223. https://doi.org/10.1785/0120040177 Google Scholar CrossRef Search ADS   Toro G., 1997. Probabilistic Models of Shear-Wave Velocity Profiles at the Savannah River Site , Aiken, South Carolina, Risk Engineering, Inc. Tech. Rep., Boulder, CO, p. 72. Towhata I., 2008. Geotechnical earthquake engineering, in Springer Series in Geomechanics and Geoengineering , Springer, Berlin. doi:10.1007/978-3-540-35783-4. Google Scholar CrossRef Search ADS   Vucetic M., Dobry R., 1991. Effect of soil plasticity on cyclic response, J. Geotech. Eng. , 117, 89– 107. https://doi.org/10.1061/(ASCE)0733-9410(1991)117:1(89) Google Scholar CrossRef Search ADS   Zafarani H., Vahidifard H., Ansari A., 2013. Prediction of broadband ground-motion time histories: the case of Tehran, Iran, Earthq. Spectra , 29, 633– 660. https://doi.org/10.1193/1.4000150 Google Scholar CrossRef Search ADS   Zafarani H., Vahidifard H., Ansari A., 2012. Sensitivity of ground-motion scenarios to earthquake source parameters in the Tehran metropolitan area, Iran, Soil Dyn. Earthq. Eng. , 43, 342– 354 https://doi.org/10.1016/j.soildyn.2012.07.007 Google Scholar CrossRef Search ADS   Zeng Y., 1993. Theory of scattered P-and S-wave energy in a random isotropic scattering medium, Bull. seism. Soc. Am. , 83, 1264– 1276. Zeng Y., Anderson J.G., Su F., 1995. Subevent rake and random scattering effects in realistic strong ground motion simulation, Geophys. Res. Lett. , 22, 17– 20. https://doi.org/10.1029/94GL02798 Google Scholar CrossRef Search ADS   Zeng Y., Anderson J.G., Yu G., 1994. A composite source model for computing realistic synthetic strong ground motions, Geophys. Res. Lett. , 21, 725– 728. https://doi.org/10.1029/94GL00367 Google Scholar CrossRef Search ADS   Zeng Y., Su F., Aki K., 1991. Scattering wave energy propagation in a random isotropic scattering medium: 1. Theory, J. geophys. Res. , 96, 607– 619. https://doi.org/10.1029/90JB02012 Google Scholar CrossRef Search ADS   © The Author(s) 2018. Published by Oxford University Press on behalf of The Royal Astronomical Society. http://www.deepdyve.com/assets/images/DeepDyve-Logo-lg.png Geophysical Journal International Oxford University Press

# Broad-band simulation of M7.2 earthquake on the North Tehran fault, considering non-linear soil effects

, Volume 213 (2) – May 1, 2018
15 pages

Publisher
The Royal Astronomical Society
ISSN
0956-540X
eISSN
1365-246X
D.O.I.
10.1093/gji/ggy038
Publisher site
See Article on Publisher Site

### Abstract

Summary The North Tehran fault (NTF) is known to be one of the most drastic sources of seismic hazard on the city of Tehran. In this study, we provide broad-band (0–10 Hz) ground motions for the city as a consequence of probable M7.2 earthquake on the NTF. Low-frequency motions (0–2 Hz) are provided from spectral element dynamic simulation of 17 scenario models. High-frequency (2–10 Hz) motions are calculated with a physics-based method based on S-to-S backscattering theory. Broad-band ground motions at the bedrock level show amplifications, both at low and high frequencies, due to the existence of deep Tehran basin in the vicinity of the NTF. By employing soil profiles obtained from regional studies, effect of shallow soil layers on broad-band ground motions is investigated by both linear and non-linear analyses. While linear soil response overestimate ground motion prediction equations, non-linear response predicts plausible results within one standard deviation of empirical relationships. Average Peak Ground Accelerations (PGAs) at the northern, central and southern parts of the city are estimated about 0.93, 0.59 and 0.4 g, respectively. Increased damping caused by non-linear soil behaviour, reduces the soil linear responses considerably, in particular at frequencies above 3 Hz. Non-linear deamplification reduces linear spectral accelerations up to 63 per cent at stations above soft thick sediments. By performing more general analyses, which exclude source-to-site effects on stations, a correction function is proposed for typical site classes of Tehran. Parameters for the function which reduces linear soil response in order to take into account non-linear soil deamplification are provided for various frequencies in the range of engineering interest. In addition to fully non-linear analyses, equivalent-linear calculations were also conducted which their comparison revealed appropriateness of the method for large peaks and low frequencies, but its shortage for small to medium peaks and motions with higher than 3 Hz frequencies. Earthquake ground motions, Site effects, Wave propagation 1 INTRODUCTION The eastern segment of the North Tehran fault (NTF), at the southern boundary of the Alborz range and North of the city of Tehran, can produce earthquakes with magnitudes M7.2–7.4 (Anderson 1997; Ritz et al. 2012). The return period for such earthquakes is estimated about 1000–1150 yr and the last earthquake which is historically corresponded to the NTF activity has occurred about 900 yr ago (Berberian & Yeats 1999). Considering the old buildings at some parts of Tehran and the highly populated city, occurrence of such event can be catastrophic. Hence, attaining a quantitative estimate of the ground motions at various parts of the city due to the next M7.2 earthquake is extremely necessary and require comprehensive researches. Results of such studies can be used in risk assessment and municipal decision makings of the Tehran metropolitan. On the other hand, ground motion prediction equations suffer from increasing uncertainties at near-fault regions and their predictions are not reliable for short distances (Dalguer & Mai 2012). In a previous paper (Majidinejad et al. 2017), we presented results for 30 dynamic rupture scenarios of an M7.2 earthquake on the NTF for frequencies up to 2 Hz. The simulation also included the Tehran basin's velocity structure. Results revealed the significant effects of the surface waves on the ground motions of the city. In the following study, results are extended up to 10 Hz in order to consider the whole range of the engineering frequencies. Broad-band ground motions are obtained by combining the low frequencies (0–2 Hz) from dynamic simulations with high frequencies (2–10 Hz) from a scattering method. Site local effects, which become considerable at high frequencies, are also included for sites on a cross-basin profile. These effects are investigated by linear, equivalent-linear and fully non-linear soil models and effect of non-linearity on soil response is assessed. 2 BACKGROUND In the previous paper (Majidinejad et al. 2017), we used a 2-D model to simulate the spontaneous dynamic rupture of the NTF and propagation of the seismic waves in the heterogeneous medium of the Tehran basin. Results of that study provided deterministic aspects of the next Tehran M7.2 earthquake. However, considering huge computational cost, our limited knowledge from the earth's heterogeneous crust and also fault's complex behaviour, results of those dynamic simulations were limited to 2 Hz. All broad-band simulations are encountered with such limitations on the available computational resources and lack of knowledge about geophysical properties. The limit on the upper frequency of deterministic methods makes them unable to predict the high-frequency ground motion parameters. On the other hand, empirical-stochastic methods are designed to produce ground motion parameters similar to previously observed earthquakes (Irikura 1986). Several methods have been introduced for generation of high-frequency component of ground motions. Some of them use stochastic synthetics which convolve with appropriate spectrum of the source (Boore 1983, 2003) and some other use wave scattering theory to simulate high-frequency motions (Zeng et al. 1994). However, all these methods are unable to reproduce acceptable low-frequency waveforms. To overcome these limitations, combined broad-band methods emerged. These methods are based on combination of deterministic low-frequency with stochastic high-frequency motions (Mai & Beroza 2003). Taking into account non-linear behaviour of soil in deterministic low-frequency simulations increases the computational costs and time significantly. Moreover, previous studies have demonstrated that the low-frequency ground motions are not affected considerably by non-linear response of the shallow soil layers. In contrast, it is proven that non-linear soil effects are crucial for high-frequency contents, especially during large earthquakes (Roten et al. 2012). To consider the soil non-linearities in microzonation and geotechnical engineering studies, it is a well-known approach to employ a 1-D non-linear soil column model as a representation of shallow soil layers and calculate its response due to seismic excitation at depth (Hashash & Park 2001; Hartzell et al. 2002). Moreover, since the computations for non-linear behaviour of the soil is time consuming and requires several iterations for each time increment, an equivalent-linear analysis can be used instead. In this method (Seed & Idriss 1969), dynamic soil parameters are kept constant over whole duration of seismic excitation at each iteration of the repetitive linear soil response calculations. The maximum shear strain occurred at each point during excitation is used to determine the appropriate shear modulus and damping ratio for the next iteration (i.e. linear solution). Equivalent-linear method is primarily designed for small strain levels in which soil complex non-linearities are not prevailed. It is assumed that shallow layers of soil, which are not exactly in the vicinity of the fault, meet the criteria for using this method in calculations of the earthquake ground motions as it is used widely by geotechnical engineers (Towhata 2008). The single previous study on the broad-band simulation of the NTF’s probable earthquake was performed by Zafarani et al. (2012, 2013). They produced PGA and PGV maps of the Tehran by combining finite-fault method (Motazedian & Atkinson 2005) results for low-frequency motions and Boore's (2003) stochastic method synthetics for high-frequency component. They considered linear soil amplifications by implementing generic hard rock amplification factors (for the whole region of the study) proposed by Boore & Joyner (1997). They concluded that the highest Peak Ground Accelerations (PGAs) occur on the hangingwall side of the fault with a value approximately equal to 700 cm s−2. The recent study performed by Majidinejad et al. (2017) predicted the highest mean SA values in the frequency range of 0–2 Hz occur on the footwall side of the NTF and that there is a significant increase in the values over the deepest parts of the Tehran basin, due to the constructive interference of body shear and surface Rayleigh waves. However, restriction of these results to the upper limit of 2 Hz, prevents a comparison between them and Zafarani et al.’s (2013) broad-band PGAs. The main goal of this study is to provide broad-band synthetic ground motions for Tehran based on the low-frequency results of Majidinejad et al. (2017) and also to estimate the non-linear soil effects on ground motion intensities using soil profiles in a cross-section of the Tehran basin. 2.1 Low frequency (0–2 Hz) The low-frequency ground motions are obtained by 17 scenarios of spontaneous rupture propagation. The NTF is considered as a purely reverse fault with a 75° dip angle. Each rupture scenario is generated by a self-similar initial stress distribution on the fault plane that produces kinematic fault properties such as rupture velocity and length, slip rise time and final slip distribution corresponding to an M7.2 earthquake. Another important aspect which is considered by the low-frequency models is the deep sedimentary basin beneath Tehran. Fault-normal and vertical components of the low-frequency ground motions reveal remarkable directivity effects of the updip rupture propagation. LF ground motion components are amplified significantly due to the low-velocity layers of the basin on the footwall side, near downtown. Results also show that the formation of Rayleigh waves at the edges of the basin changes the frequency content and duration of the LF components (Majidinejad et al. 2017). 2.2 High frequency (2–10 Hz) In this study, we use the ‘Broad-Band toolbox’ code, developed by Mai et al. (2010), to generate and combine the HF components with the LF motions. The HF motions are produced based on a method proposed by Zeng et al. (1991) and Zeng (1993) for point sources and developed by Mena et al. (2010) for extended faults. The background theory of the method will be shortly described here. As the seismic waves travel from source to stations, they scatter by various distributed heterogeneities in the medium. Such phenomenon causes dispersion of the seismic wave's energy and reduces energy density of the directly arrived waves to site. However, in this case, the reduction of the energy density is different from other dissipative mechanisms (such as viscoelastic or plastic behaviour of the material) and the energy is transferred to secondary reflected waves which appear as coda in time histories. Zeng (1993) used S-to-S backscattering theory and proposed an equation for energy partitioning of a spherical wave into direct, first- and second-scattered and multiple (orders higher than 2) scattered waves. Eq. (1) represents the multiple scattered portion of the wave's energy at the time t due to an impulse excitation in the reference point of an infinite 3-D medium.   \begin{eqnarray} &&{ {E_{{\rm{multiple}}}}(\vec{r},t)}\nonumber\\ &&= \mathop \smallint \limits_{ - \infty }^{ + \infty } \frac{{{e^{i{\rm{\Omega }}}}}}{{2\pi }}d{\rm{\Omega }} \times \mathop \smallint \limits_0^\infty \frac{{{{\left( {\frac{{{\eta _s}}}{k}} \right)}^3} \times {{\left[ {{\rm{ta}}{{\rm{n}}^{ - 1}}\left( {\frac{k}{{\eta + {{{i{\rm{\Omega }}}}/{v}}}}} \right)} \right]}^4} \times {\rm{sin}}\left( {kr} \right)}}{{2{\pi ^2}vr\left[ {1 - \frac{{{\eta _s}}}{k}{\rm{ta}}{{\rm{n}}^{ - 1}}\left( {\frac{k}{{\eta + {{{i{\rm{\Omega}}}}/{v}}}}} \right)} \right]}}dk \end{eqnarray} (1) In this equation, ‘r’ is the source-to-receiver distance, ηs and η are the scattering and total shear wave attenuation coefficients and Vs is the shear wave speed. In the method used for producing the HF motions, the scattering Green's functions are generated by multiplying random wavelets (with average energy of unity) to the energy envelope calculated by eq. (1) for each station (Mai et al. 2010). Mena et al. (2010) developed the method from point source to an extended model by introducing the extended fault as a composition of small subfaults. They calculated the contribution of subfaults to the overall ground motion by taking into account the delay of slip initiation on each subfault due to its distance from hypocentre and also delay due to the subfault to station distance. The HF motion attributed to each subfault is obtained by convolution of the subfault's source time function (STF) to the scattering Green's function. The amplitude of the STF for each subfault is determined by scaling laws and with respect to the seismic moment on that subfault (Irikura 1986). By this method, synthetic high-frequency ground motions that take into account physical and random features of the fault and medium to some extent are produced. However, it should be noted that since eq. (1) is derived for an infinite medium, only body waves are measured and therefore, the effect of surface scattered waves is not reflected in the HF component. The attenuation and scattering coefficients which represent density and size of heterogeneities in the medium and control the coda envelope are set to 0.06 and 0.03 km − 1, respectively, with regard to Mena et al. (2010). We use the frequency-dependent attenuation model proposed by Farahani et al. (2012) for the Tehran basin, as used for the LF simulation in Majidinejad et al. (2017). In order to model the HF site-specific attenuation, kappa coefficient is set to κ = 0.04 s, according to Hassani et al. (2011) results for the region. Dreger et al.’s (2007) STF is used for calculation of high-frequency motions. This slip rate function is compatible with the slip functions generated by dynamic spontaneous rupture models, as it can provide the overall shape and rapid rise time of the dynamic models appropriately. A comprehensive study by Mena et al. (2010) showed that in comparison with traditional boxcar and triangular STFs and also dynamically consistent Yoffe (Tiniti et al. 2005) and Liu et al. (2006)’s STFs, Dreger's STF is the most appropriate, since it has a smooth spectrum and does not contain any spectral holes or abrupt spectral changes at low frequencies. 2.3 Broad band (0–10 Hz) The Broad-Band toolbox combines the LF and HF components of the motion with the method proposed by Mai & Beroza (2003), which provides a solution to the inconsistency of the energy spectrums of LF and HF contents near matching frequency (Thráinsson et al. 2000). The main assumption of the method is that the HF content of motion is not completely stochastic and contains some physical aspects of the rupture propagation and emitted seismic waves, such as arrival phases of P and S waves, magnitude of the motions, etc. This criterion is met by HF ground motions generated by Zeng et al. (1995) and Mena et al. (2010) method. Mai and Beroza's (2003) method searches for the matching frequency within a bandwidth such that the phase spectrum of the HF component matches optimally with the realistic LF phase spectrum. Then amplitude spectrum of the HF component is scaled by a close to unity factor to the LF amplitude at the matching frequency. Since physical characteristics of the earthquake are considered during calculation of the HF component of motions, the scale factor used for fitting the amplitude spectrums is often close to unity. By this method, the matching frequency may alter for each station and direction of motion. In this study, the bandwidth for the matching frequency search is set to 1.8–2 Hz and the maximum ground motion frequency is limited to 20 Hz. However, the resulting broad-band ground motions are low-pass filtered up to 10 Hz in order to facilitate the site effect calculations and avoid the numerical instabilities. Synthetic broad-band ground motions are provided at stations with 1 km distance strides on a north–south profile of the city. 2.4 Tehran alluviums Tehran is spread over three stratigraphic zones. The Southern Alborz Mountains foothills, north of the city, comprise the first zone and are mainly composed of shale, Eocene tuffs, andesite, basalt and pyroclastic from Palaeozoic, Mesozoic and Tertiary eras (Emami et al. 1993). The second zone, which includes Sepayeh Mountains at the east and Bibisharbano Mountains to the south of the city, is made up of Triassic period materials such as limestone, marl limestone and dolomite (Shafiee & Azadi 2007). The Tehran alluvial plain, as the third stratigraphic zone, is located above the Tehran basin and consists of Quaternary formations which are mostly accumulated by erosion and redeposition of former sediments. The alluvium spreads from foothills at the north to the low farmlands at the south of Tehran. It is composed of varying thickness layers with materials of different types and size (Rieben 1955; Jafari et al. 2002). Pedrami (1981) categorized the Tehran plain alluvium to five geological units, including units A and Bn in the north, unit Bs in the south, unit C in the north, west and centre and unit D in centre and south (Fig. 1). Jafari et al. (2001) collected data from 700 drilled boreholes over the study area to investigate the material types of the geological units. Results identified unit A, the oldest alluvial layer, mainly consists of gravel, sand and silt. They proposed unit B to be a conglomerate mixture of cobbles, gravels, sand and silt and also dolomite and limestone at some southern parts. Unit C is characterized by a mixture of pebble, sand, silt and clay and the youngest layer, unit D, is composed of soft silty-clay materials. Shafiee & Azadi (2007) used data from 188 seismic profiles, distributed throughout the city, to classify the geological units of the Tehran plain according to their shear wave velocity. They reported average Vs, 30 for the Tehran plain as in Table 1. Figure 1. View largeDownload slide Map of geological units of Tehran (Shafiee & Azadi 2007). The thick dashed black line represents top edge of the NTF’s rupture plane. The thick black line stands for the N-S profile of the study. Velocity layers of deep basin of Tehran in the N-S profile are shown in the left. Figure 1. View largeDownload slide Map of geological units of Tehran (Shafiee & Azadi 2007). The thick dashed black line represents top edge of the NTF’s rupture plane. The thick black line stands for the N-S profile of the study. Velocity layers of deep basin of Tehran in the N-S profile are shown in the left. Table 1. Characteristics of geological units of the Tehran plain. Geological unit  Constituting materials  $${\bar{V}_{s,30}}({\rm{m}}\,{{\rm{s}}^{ - 1}})$$  PI (%)  A  Conglomerate with silt–sand–gravel and silt–clay mixtures  830  0  B  Cobble, boulder, gravel and sand  770  0  C  Gravel, sand, silt and clay  400  12  D  Silt and clay  280  25  Geological unit  Constituting materials  $${\bar{V}_{s,30}}({\rm{m}}\,{{\rm{s}}^{ - 1}})$$  PI (%)  A  Conglomerate with silt–sand–gravel and silt–clay mixtures  830  0  B  Cobble, boulder, gravel and sand  770  0  C  Gravel, sand, silt and clay  400  12  D  Silt and clay  280  25  View Large 2.5 Soil properties One of the deterministic parameters in non-linear soil behaviour is the plasticity index (PI). In this study, the PI for different soil types are assigned based on the geological units of Tehran alluviums and also based on the amount of clay content in soil composition (Hartzell et al. 2002; Roten et al. 2012). The fine-grained soil in the southern parts of the Tehran consists of clay and silt. Jafari et al. (2002) performed a comprehensive investigation on alluviums of the southern Tehran by evaluating dynamic properties of the fine-grained soils using dynamic laboratory tests and field geoseismic investigations. The laboratory tests included torsional resonant column tests in a shear strain range of 10 − 6–10 − 4and stress-controlled cyclic triaxial tests in a wider range of shear strains (10 − 5–10 − 2) to assess the dynamic properties of the soils collected from sites at the southern parts of the city. They categorized the fine-grained soil in the area to three types with respect to the PI of the samples. Silt specimens were characterized by very low PI (PI < 7); silty clay samples were identified by low index (7 < PI < 15) and medium PIs (15 < PI < 30) were assigned to clay samples. They concluded that an increase of the confining pressure on the soil samples leads to more shear modulus ratios ($${G/{G_0}}$$) and lower damping ratio, however, this effect reduces as the soil plasticity increases. Finally, they have compared their results by Vucetic & Dobry (1991) models and concluded that for the medium PI clay, the results agree very well with 15 < PI < 30 in the model. Ghorbani et al. (2014) proposed the value of PI for the fine silty clay soil of Tehran as 12 per cent. Therefore, based on the main contents of the Tehran geological units, we categorized the Tehran soil to three types of gravel and sand, silty clay and clay and assumed their PIs as (Table 1). In order to assign the appropriate soil properties at each site, we used the reports of JICA (2000) and Keshavarz Bakhshayesh 2003 as well as the geological map of Shafiee & Azadi (2007) to determine the soil type at each depth due to the alluvium profiles and regional boreholes. Since the younger Alluvium of the southern Tehran are composed of remarkable clay contents, we used a PI of 25 per cent for unit D. Geological unit C has lower clay and higher silt, so it is reflected in the lower PI value of 12 per cent. Finally, the oldest A and B units are mostly composed of sand and gravel, therefore the negligible clay content is represented by PI = 0 per cent. Tehran alluviums have experienced many loading and unloading cycles such as fluctuations of the water table and desiccation during deposition, deposition of upper layers, aging effects and also cementation and upward motion of lower deposits due to Quaternary faults movements. So, due to the geological history of Tehran sedimentary basin, it is plausible to consider overconsolidated soil in the region. Several studies are established to estimate the overconsolidation ratio (OCR) of the Tehran alluvium (Jafari et al. 2000). These studies suggest that south of Tehran soils are generally characterized by a low OVR, approximately about 1–1.5. On the other hand, after numerous laboratory investigations, Darendeli (2001) stated that the OCR in comparison with other factors has minor effect on the dynamic properties of soil samples. Therefore, we used the OCR = 1 for all types of the soils in Tehran, both in calculation of at-rest lateral earth pressure coefficient (K0) and Darendeli's equations. Several models have been proposed by researchers to predict the soil non-linear behaviour (Vucetic & Dobry 1991; Darendeli 2001; Amir-Faryar & Aggour 2016). However, regional studies on the non-linear behaviour of various soil types, around the world, show that some modifications to the parameters used in the proposed models are necessary in order to provide better fit between predictions and the local experimental results (e.g. Bay & Sasanakul 2005). We used Darendeli's (2001) model which takes into account the consolidation state, effective stress and PI in calculation of soil dynamic behaviour. The model proposed by Darendeli (2001) also requires the number of cycles (N) and the frequency of excitation ( f). Jafari et al. (2002) set these parameters to N = 10 and f = 1 Hz with respect to the Standard (2003) ASTMD 3999 recommendation and in this study we use similar values for the parameters of the Darendeli's model (following Roten et al. 2012). Selected parameters are verified with the experimental data obtained by Jafari et al. (2002) for the Tehran region (Fig. 2). Figure 2. View largeDownload slide Comparison of Darendeli's (2001) model and experimental results for the Tehran soil (Jafari et al. 2002) for two typical soil samples with medium and low plasticity indices. Figure 2. View largeDownload slide Comparison of Darendeli's (2001) model and experimental results for the Tehran soil (Jafari et al. 2002) for two typical soil samples with medium and low plasticity indices. Due to the computational limits, the minimum shear wave velocity in the LF spectral element model was constrained to 600 m s−1 on soft quaternary and 800 m s−1 on Tertiary deposits (Majidinejad et al. 2017). The velocity model used for shallow soil profiles is proposed based on reflection and borehole tests (JICA 2000) and is consistent with the velocity model used for the LF 2-D simulations. JICA (2000) proposed 40 shallow profile models over Tehran, categorized according to shear wave velocity, material type, thickness of the layers and the depth to the geotechnical bedrock. The maximum velocity of the JICA profiles at deepest layers are 550 and 800 m s−1 for quaternary and Tertiary deposits, respectively, which coincide with the minimum velocities in 2-D spectral element model uppermost layers. 2.6 Analysis model In order to perform the site response analysis with linear, equivalent-linear and fully non-linear models, we used the code DEEPSOIL, developed by Hashash et al. (2012), which utilizes multiple lumped mass methods for simulation of wave propagation through the soil layers. Since the soil dynamic parameters such as stiffness ($${G/{G_0}}$$) and damping (ξ) are obtained under harmonic loading tests, hence, in the case of equivalent-linear analyses, direct implementation of the maximum strain during seismic excitation overestimates the stiffness reduction and hysteresis attenuation level. Therefore, there is a common correction method which applies a fraction of maximum strain (effective strain) to update $${G/{G_0}}$$ and ξ values from experimental curves. This ratio should reflect the non-linearity level of the model during an earthquake with magnitude M, thus an effective shear strain ratio (SSR) defined as $${\rm{SSR}} = {{( {M - 1} )} /{10}}$$ would be proper (Idriss & Sun 1992). The fully non-linear analyses employ the hyperbolic modified (Kondner & Zelasko 1964) model for the backbone curve formulation (Matasović & Vucetic 1995). The model is extended to pressure-dependent in DEEPSOIL by making the shear strength as a function of the confining pressure (Hashash & Park 2001). Masing's (1926) model which is widely used for constructing cyclic loops is applied in DEEPSOIL for representing the cyclic loading/unloading behaviour of the soil (Hashash et al. 2012). Implementing Masing rules force the cyclic loops to follow the backbone curve but overestimate the damping ratio at large strains (Stewart et al. 2014). Phillips & Hashash (2009) proposed a modification to apply a pinching manner to the cyclic loops of the soil at large strains. This method is implemented in DEEPSOIL (Stewart et al. 2014). Roten et al. (2012) and Hartzell et al. (2002) showed that even under large earthquake excitations, the layers below ∼90 m depth remain elastic. On the other hand, due to the limited depth extent of some profiles of the JICA study, we extend our profiles with a velocity of 600 m s−1 (or 800 m s−1 correspondingly) to a depth of 100 m for these profiles, following Roten et al. (2012) and Hartzell et al. (2002). Also an elastic half-space is defined as the boundary condition beneath the soil column to prevent reflection of waves from bottom and specifically formation of unrealistic resonances at the natural frequencies of soil columns (Roten et al. 2012). We introduced a shear wave velocity of 800–1200 m s−1 for the bedrock to represent the average effect of shallower layers of the deep basement. Ground motions at the bedrock level are obtained from combination of the LF motions at the free surface of the spectral element simulation with the HF motions calculated for the same fault and velocity structure model. In fact, due to its negligible effect on ground motions in comparison with crustal thickness, we have ignored the thickness of the shallow soils (with maximum depth of 50–100 m) in the first step of ground motion calculations (i.e. generation of LF and HF motions) as done by Hartzell et al. (2002). The input ground motions in DEEPSOIL are defined as ‘outcrop’ motions since they are obtained from the free surface of the spectral element model. The ‘outcrop’ motion is the motion that occurs at a free surface and is simply twice the upward propagating wave motion. Using ‘outcrop’ input ground motions is consistent with the implementation of elastic half-space as the boundary condition of soil column. The water table depth in Tehran decreases from north to the south, consistent with the topography and hydraulic gradient of the underground water. At the elevated parts of the city, which are characterized by hard and dense soils and rocks, shallow ground water is not available. However, at the central and southern parts of the city ground water depth changes from 50 to 15 m as directing southwards. Oscillations in pore-water pressure change dynamic soil properties at each time step even at small strains, for example maximum shear modulus and maximum shear stress. Consequently, implementation of the effective stress in calculations is too time consuming against total stress calculations, specifically for studies with large amounts of models. Moreover, seepage analysis parameters such as permeability coefficients are not well constrained for Tehran soils, therefore, following Roten et al. (2012) and Hartzell et al. (2002), we use the total stress approach in calculations. Annual measurements of the ground water level, attained by IWRMC (Iran Water Resource Management Company), are used for total stress calculations in this study (Keshavarz Bakhshayesh 2003). 3 RESULTS AND DISCUSSION We applied 50 s low-frequency time histories to generate broad-band synthetics at the bedrock level and then studied effects of local soil response with linear, equivalent-linear and non-linear soil column models. Linear soil response is compared with non-linear response for each station to evaluate the effect of soil non-linearities on the ground motions. 3.1 Broad-band synthetics at bedrock Low-frequency synthetic ground motions are provided for all of the 42 stations (with 1 km intervals) along the N-S profile obtained from the set of 17 low-frequency simulations in the previous study (Majidinejad et al. 2017). Each of these LF synthetics are combined with seven high-frequency scatterograms to achieve 119 broad-band 0–10 Hz ground motions. Fig. 3 depicts the LF and broad-band velocity and acceleration time histories and amplitude spectra of a typical scenario model for two stations of the footwall side along the N-S profile. One of the stations is located near the fault trace and out of the basin (Rx = + 2 km) and the other is located further above the deep parts of basin (Rx = + 14 km). As it is observable, the amplitudes of the LF model fall at frequencies above 2 Hz, while the broad-band synthetic amplitudes are considerable up to 10 Hz. As it can be seen, the general form of the broad-band time histories are consistent with real observations in the sense that ground motions first begin with HF motions and then continue with arrival of the main LF shear waves (Hartzell et al. 2004). Late surface Rayleigh waves with large amplitudes are visible at the station above the deep Tehran basin. By achieving broad-band synthetics, we are now able to calculate the HF spectral amplitudes and PGA of the expected ground motions in Tehran for an M7.2 earthquake induced by the NTF. As comprehensively discussed in Majidinejad et al. (2017), we assume that the fault-normal component is the predominant component of the ground motions due to the NTF earthquake and contains major portion of the seismic wave's energy taking into account the fault mechanism and 2-D geometry of the fault and basin. Figure 3. View largeDownload slide LF and broad-band ground motion time histories and amplitude spectra at two typical stations. Figure 3. View largeDownload slide LF and broad-band ground motion time histories and amplitude spectra at two typical stations. Fig. 4 shows the comparison of PGA and response spectra values from broad-band simulations at the Bedrock level with the Boore & Atkinson (2008) GMPEs. Results demonstrate very good agreement with the empirical relationships on the average values. Almost all the data points, except the ones corresponding to the deepest parts of the basin, fall within one standard deviation of the empirical relationships. Differences between average values of the synthetic ground motions parameters and empirical values are less than 50 per cent of the empirical standard deviations at all bins. It is noteworthy to point out that the inclusion of the high-frequency contents in broad-band synthetics improves the match between synthetics and empirical relationships, even at low-frequency SAs. Amplification of the SA values due to the deep sedimentary basin of Tehran is also present at high frequencies and HF SAs follow a similar trend to the trend of 1.5–2 Hz SAs in the LF models. This can be attributed to the scaling of the Fourier amplitudes of the HF motions to the LF amplitudes during generation of BB ground motions. However, there are more spatial variations in the HF SA values with respect to the LF SA values due to the stochastic manner of HF contents (Majidinejad et al. 2017). Figure 4. View largeDownload slide Ground motion parameters obtained from broad-band simulation at the geotechnical bedrock level compared with the empirical equation of Boore & Atkinson (2008). Black bars represent for simulation results, including the standard deviation (error bars), binned into logarithmically spaced distance ranges. The thick black line corresponds to the mean values of empirical GMPE; the dashed line show mean ± σ of the empirical GMPE correspondingly. Figure 4. View largeDownload slide Ground motion parameters obtained from broad-band simulation at the geotechnical bedrock level compared with the empirical equation of Boore & Atkinson (2008). Black bars represent for simulation results, including the standard deviation (error bars), binned into logarithmically spaced distance ranges. The thick black line corresponds to the mean values of empirical GMPE; the dashed line show mean ± σ of the empirical GMPE correspondingly. In order to take into account the soil column behaviour, for each LF model, we selected the broad-band time history that produces the nearest responses to the mean values obtained from all correspondent LF models. Therefore, we attained 17 final broad-band models—with dissimilar LF and HF contents—which are then used for investigation of the local site effects at 30 stations situated in the municipal regions of Tehran on the N-S profile (i.e. a total number of 510 simulations for each of the linear and non-linear soil columns analyses). 3.2 Linear results Soil column models are set such that minimum resolvable frequency of 30 Hz is available for all models, with layers defined at 2 m interval depths. Fig. 5 compares the average 3 and 5 Hz SAs, obtained by linear amplification of the broad-band synthetic ground motions, with the empirical GMPE of Boore & Atkinson (2008) in the N-S profile of city at the free-surface level. The required parameters for the empirical relationships, such as Vs, 30, Rx, etc., are evaluated at each site with regard to the values in the models. Although the average values of the broad-band synthetics parameters at the bedrock level appeared to be within one-half standard deviation of the empirical relationships, linearly amplified values of these parameters at the surface do not coincide with the empirical relationships. The mean values differ more than one standard deviation from the empirical observations at many parts of the N-S profile, in particular over the deep parts of the basin. This inconsistency should be attributed to the insufficient accuracy of the linear modeling in representing soil behaviour during severe ground motions. While central and southern parts of the city experience considerable amplifications because of shallow soil layers, the northern parts are not affected significantly by local site effects, due to the shallow bedrock depth on the hangingwall side (Fig. 5). Figure 5. View largeDownload slide Average spectral accelerations at the ground surface level of the stations along the N-S profile obtained by linear, equivalent-linear and fully non-linear analyses. Empirical values for each station on the profile is also calculated based on site specific values and presented for comparison. Figure 5. View largeDownload slide Average spectral accelerations at the ground surface level of the stations along the N-S profile obtained by linear, equivalent-linear and fully non-linear analyses. Empirical values for each station on the profile is also calculated based on site specific values and presented for comparison. Maximum linear SA values occur on the footwall side of the fault and reach up to ∼2.4 and ∼3.7 g on average at 3 and 5 Hz, respectively. On the hangingwall side, the maximum linear SAs at these frequencies are about ∼1.1 and ∼1.6 g which occur in the vicinity of the fault. A notable remark is that, in agreement with Dalguer & Mai (2012) and in contrast with empirical relationships, the maximum values of responses do not occur just in the vicinity of the fault, but in a distance 2–3 km from its trace. Also, it can be seen that at stations very close to the fault, results are lower than empirical predictions. 3.3 Non-linear results Many energy dissipation mechanisms occur at the shallow layers of soil during large earthquake ground motions, therefore, a non-linear model of the soil which takes into account soil's stiffness degradation and hysteresis damping may be capable of reducing the unrealistic high values of the linear amplification. In order to evaluate the soil non-linear effects on the HF ground motions, we performed both full non-linear and equivalent-linear analyses. The reference strain and damping ratios are defined at 2 m intervals of the vertical profiles. We also assumed that the soil layers with shear wave velocities above 750 m s−1 do not undergo non-linear strains, due to their stiff and dense soil compositions. Linear, non-linear and equivalent-linear results of a scenario earthquake are compared in Fig. 6 for a station at distance of 13 km from the fault trace, on the footwall side and above the basin. As it is clear in Fig. 6, linear soil response amplifies the ground motions significantly which leads to large apexes in time history, while, these peaks are lowered and also smoothed by the non-linear behaviour of the soil column. It is noteworthy that the reductions due to the soil non-linearities are more pronounced at the high amplitude peaks, which are followed by large inelastic strains. Comparison of the non-linear and equivalent-linear results shows that although equivalent-linear method predicts the deamplifications at large peaks properly, it also reduces the modest and small apexes of the time history which results in underestimation of medium peaks and a somewhat smoothed time history. This observation clearly associates with the assumption of the equivalent-linear method, which utilizes the stiffness and damping ratio calculated for the maximum peak at all time steps. In this case, the equivalent-linear method has underpredicted the spectral accelerations with respect to the fully non-linear analysis at frequencies above 3 Hz (Fig. 6, right). Figure 6. View largeDownload slide Effect of shallow soil layers on broad-band ground motions of a station located 13 km from fault's trace, on the footwall side. Left: time histories and right: response spectral accelerations. Figure 6. View largeDownload slide Effect of shallow soil layers on broad-band ground motions of a station located 13 km from fault's trace, on the footwall side. Left: time histories and right: response spectral accelerations. Comparison between linear and non-linear soil responses reveals considerable decrease in spectral acceleration values of the HF content and also increase in soil column resonant periods due to the reduced stiffness and elevated damping of the non-linear soil layers (Fig. 6, right). As stated by Hartzell et al. (2004) and also reported by Roten et al. (2012), the opposing effects of stiffness reduction and increase of damping in non-linear soil response, determine further amplification or deamplification of the responses with respect to the linear behaviour. Since damping effect is predominant at high frequencies (>2 Hz), the non-linear responses at these frequencies are usually below the linear values. On the other hand, at low frequencies (<1 Hz) damping effects are minor, therefore, non-linearities lead to higher amplifications at this range (Fig. 6, right). Fig. 5 also displays the average spectral accelerations in the case of non-linear soil response for comparison with ground motion prediction equations. Intense SA values above the basin are diminished in this case and average values of numerical results are now limited to the mean ± one standard deviation of the GMPEs at many stations. The maximum non-linear deamplification in response takes place at a station 18.5 km from the fault, on the footwall side. At this site, the non-linear deamplification ratio (i.e. $$1 - {{{\rm{S}}{{\rm{A}}_{{\rm{NL}}}}} / {{\rm{S}}{{\rm{A}}_{{\rm{Lin}}}}}}$$) reaches to 33, 60 and 63 per cent for frequencies of 2, 3 and 5 Hz, respectively. On the footwall side of the NTF, in a distance range of 0–8 km from the fault trace where the soil profiles depth to the geotechnical bedrock is not very high, due to the hysteresis behaviour in few thin soil layers, dissipation of energy is negligible. Hence, these stations experience low deamplifications or even higher amplifications as a result of stiffness reductions in the non-linear response (particularly at low frequencies). At distances within 8–22 km from the fault trace, depth of low-velocity layers to the geotechnical bedrock increases and consequently, at the same time with the reduction of dynamic stiffness during non-linear soil behaviour, the amount of hysteresis damping increases significantly. This highly elevated dissipation in non-linear response leads to considerable reduction of responses with respect to the linear behaviour. Although stations further than 22 km are located on weak soil layers, low bedrock motions at these stations are unable to develop significant plasticity in soil columns, therefore, these sites exhibit lower non-linear deamplifications. Figs 7(a) and (b) show the vertical profiles of maximum accelerations, averaged over scenario earthquakes, at stations located 9.5 and 18.5 km from the fault trace and for the case of linear and non-linear analyses. Both station are located on the footwall side and above the deep sedimentary basin (Fig. 1) but one has medium depth, low-velocity layers with average shear wave velocity of Vs, 30 = 400 m s − 1  at the top of the geotechnical bedrock (Rx = 9.5 km) and the other is located above thick very low-velocity layers, Vs, 30 = 220 m s − 1, down to the geotechnical bedrock (Rx = 18.5 km). Figs 7(c) and (d) depict the hysteresis loops of soil at depth of 30 m for these stations, which clearly show the extent of plasticity in soil columns. Limited plasticity at the station out of the basin explains similar linear and non-linear responses. Figs 7(e) and (f) are hysteresis loops at depth of 90 m and confirm the presumption of very limited plasticity at depths below 100 m. Figure 7. View largeDownload slide Effect of soil's plasticity development on deamplification caused by non-linear soil response at two representative stations. (a) and (b) vertical profiles of average PGA in linear and non-linear response. (c) and (d) hysteresis loops of the soil layer at depth of 30 m. (c) and (d) hysteresis loops at depth of 90 m. Figure 7. View largeDownload slide Effect of soil's plasticity development on deamplification caused by non-linear soil response at two representative stations. (a) and (b) vertical profiles of average PGA in linear and non-linear response. (c) and (d) hysteresis loops of the soil layer at depth of 30 m. (c) and (d) hysteresis loops at depth of 90 m. 3.4 Non-linear deamplification models Although we investigated just one N-S profile along Tehran, by considering JICA (2000) soil profiles and taking into account depth of the geotechnical bedrock, geological units and average shear wave velocity, it appears plausible to suggest three major site classes I, II and III for the Tehran plain as in Fig. 8. Sites located in each site class are roughly assumed to represent similar responses to the same bedrock excitations. Class I includes shallow bedrock profiles (depth <30 m) and coarse-grained materials with approximately high shear wave velocities (Vs, 30 above 450 m s−1). Class II consists of mixture of sand and clay and the bedrock depth at this class changes from 35 to 60 m. Average shear wave velocity (Vs, 30) at this class ranges between 300 and 450 m s−1. Finally, class III is characterized by deep bedrock (75m< depth) and fine clayey materials with low shear wave velocities (Vs, 30 < 300 m s−1). Figure 8. View largeDownload slide Map of the site classes proposed based on behaviour of soil profiles in JICA (2000). Thick black lines represent the municipal district's boundaries of the Tehran. Figure 8. View largeDownload slide Map of the site classes proposed based on behaviour of soil profiles in JICA (2000). Thick black lines represent the municipal district's boundaries of the Tehran. At this stage, we propose deamplification factors for the Tehran site classes, as done by Roten et al. (2012) for the Salt Lake City. Such results can be used for stations all over the city in order to calculate the non-linear soil responses with respect to the corresponding linear spectral accelerations. Since we aim to introduce models for any possible earthquakes in the city and regardless of the causative fault—in order to remove the source to site distance effects on input motions—we apply all the synthetic ground motions (510 time histories) to each of archetype soil profile models (10 models). Therefore, total number of 10 200 simulations are performed and results are categorized into three major site classes. Linear and non-linear simulation results are analysed to obtain statistical relations between two (linear and non-linear) sets at each class. A possible caveat here is the duration of ground motions for earthquakes with magnitudes smaller (or larger) than M7.2. Although the deamplification functions correct for distance, they do not account for shorter (or longer) duration that would be characteristic of smaller (larger) earthquakes, but which would result in different non-linear site response. Fig. 9 shows the non-linear versus linear spectral acceleration scatters for three site classes, considering all frequencies at the range of interest (top) and frequency of 5 Hz (bottom). As it can be obtained from Fig. 9, the maximum linear 5 Hz spectral accelerations reach up to 75.61, 67.56 and 50.81 m s−2 for site classes I, II and III, while the non-linear responses are limited to 31.38, 20.86 and 8.37 m s−2 correspondingly. It is notable that non-linear soil effects become evident for linear spectral accelerations larger than ∼2 m s−2. Also the dispersion of data points increases at larger SA values due to the more plasticity development and non-linearities. As anticipated site, classes I and III show the lowest and highest deamplifications to the normalized source-to-site distance ground motions, respectively, however, this conclusion may differ slightly for the earthquake taking place on the NTF as the stiffest site class I is the nearest and the loosest site class III is the furthest site class to the NTF. Figure 9. View largeDownload slide Scatterogram of non-linear SAs versus corresponding linear SAs for three site class of Tehran plain. Top: spectral accelerations at all frequencies. Bottom: spectral accelerations at frequency of 5 Hz. Black lines correspond to the correction functions calculated in the form of eq. (1) and dashed lines correspond to their standard deviations, respectively. Figure 9. View largeDownload slide Scatterogram of non-linear SAs versus corresponding linear SAs for three site class of Tehran plain. Top: spectral accelerations at all frequencies. Bottom: spectral accelerations at frequency of 5 Hz. Black lines correspond to the correction functions calculated in the form of eq. (1) and dashed lines correspond to their standard deviations, respectively. The overall trends of scattered data in Fig. 9 suggest a second-order relationship between the logarithm values of the linear and non-linear spectral accelerations, therefore, we use a relationship in the form of eq. (2) to represent the non-linear response as a function of linear spectral acceleration (Roten et al. 2012):   $$\log ({\rm{S}}{{\rm{A}}_{{\rm{NL}}}}) = b * \log ({\rm{S}}{{\rm{A}}_{{\rm{Lin}}}}) + a * {(\log ({\rm{S}}{{\rm{A}}_{{\rm{Lin}}}}))^2}.$$ (2) Table 2 lists the values and standard deviations calculated for regression parameters of eq. (2) for spectral accelerations at various frequencies in the range of interest and also PGA. Values in Table 2 can be interpreted by characteristics of the site classes. Site class I, which shows negligible deamplifications, is approximated by close to unity ‘b’ values and small values of ‘a’ parameter. Considerable deamplifications in site classes II and III, which are composed of soft materials, are reflected in high absolute values of ‘a’ parameter. Although at some frequencies, ‘a’ values of site class II are higher than corresponding values in site class III, parameter ‘b’ at these frequencies is meaningfully smaller for site class III. Small values of ‘b’ parameter imply premature initiation of non-linearities in site class III. Table 2. Parameters of regression equation (2) and their standard deviation.   Site class I  Site class II  Site class III  Frequency (Hz)  b  a  std-b  std-a  b  a  std-b  std-a  b  a  std-b  std-a  10  1.19  −0.16  0.248  0.135  1.07  −0.19  0.145  0.068  1.29  −0.52  0.313  0.216  9.5  1.12  −0.08  0.146  0.064  1.05  −0.17  0.170  0.087  0.87  −0.15  0.204  0.123  9  1.29  −0.21  0.173  0.091  1.05  −0.20  0.174  0.084  0.73  −0.11  0.126  0.057  8.5  1.18  −0.16  0.121  0.046  1.07  −0.20  0.192  0.093  0.73  −0.12  0.153  0.085  8  1.07  −0.14  0.146  0.052  1.14  −0.24  0.164  0.080  0.77  −0.13  0.182  0.107  7.5  0.96  −0.08  0.143  0.063  1.16  −0.27  0.163  0.080  0.79  −0.16  0.170  0.104  7  0.98  −0.10  0.122  0.057  1.07  −0.24  0.182  0.093  0.66  −0.09  0.184  0.106  6.5  1.01  −0.11  0.123  0.060  1.05  −0.23  0.200  0.105  0.66  −0.07  0.204  0.124  6  1.05  −0.15  0.159  0.072  1.02  −0.21  0.164  0.084  0.69  −0.10  0.130  0.065  5.5  1.02  −0.13  0.221  0.110  1.03  −0.22  0.158  0.078  0.61  −0.03  0.140  0.078  5  0.98  −0.11  0.210  0.101  1.00  −0.21  0.159  0.080  0.62  −0.06  0.138  0.077  4.5  1.04  −0.16  0.184  0.091  0.97  −0.19  0.151  0.076  0.63  −0.06  0.137  0.071  4  1.05  −0.15  0.164  0.087  0.98  −0.20  0.150  0.077  0.67  −0.08  0.149  0.085  3.5  1.02  −0.14  0.145  0.077  0.96  −0.17  0.160  0.088  0.69  −0.09  0.162  0.104  3  1.05  −0.15  0.149  0.073  0.94  −0.16  0.142  0.074  0.72  −0.12  0.157  0.097  2.5  1.00  −0.13  0.132  0.056  0.97  −0.17  0.158  0.083  0.80  −0.19  0.171  0.107  2  1.02  −0.14  0.160  0.076  0.98  −0.17  0.172  0.094  0.85  −0.22  0.155  0.091  All  0.95  −0.06  0.025  0.030  0.90  −0.10  0.024  0.035  0.80  −0.15  0.023  0.047  PGA  1.02  −0.12  0.178  0.082  0.91  −0.15  0.184  0.095  0.56  −0.03  0.170  0.119    Site class I  Site class II  Site class III  Frequency (Hz)  b  a  std-b  std-a  b  a  std-b  std-a  b  a  std-b  std-a  10  1.19  −0.16  0.248  0.135  1.07  −0.19  0.145  0.068  1.29  −0.52  0.313  0.216  9.5  1.12  −0.08  0.146  0.064  1.05  −0.17  0.170  0.087  0.87  −0.15  0.204  0.123  9  1.29  −0.21  0.173  0.091  1.05  −0.20  0.174  0.084  0.73  −0.11  0.126  0.057  8.5  1.18  −0.16  0.121  0.046  1.07  −0.20  0.192  0.093  0.73  −0.12  0.153  0.085  8  1.07  −0.14  0.146  0.052  1.14  −0.24  0.164  0.080  0.77  −0.13  0.182  0.107  7.5  0.96  −0.08  0.143  0.063  1.16  −0.27  0.163  0.080  0.79  −0.16  0.170  0.104  7  0.98  −0.10  0.122  0.057  1.07  −0.24  0.182  0.093  0.66  −0.09  0.184  0.106  6.5  1.01  −0.11  0.123  0.060  1.05  −0.23  0.200  0.105  0.66  −0.07  0.204  0.124  6  1.05  −0.15  0.159  0.072  1.02  −0.21  0.164  0.084  0.69  −0.10  0.130  0.065  5.5  1.02  −0.13  0.221  0.110  1.03  −0.22  0.158  0.078  0.61  −0.03  0.140  0.078  5  0.98  −0.11  0.210  0.101  1.00  −0.21  0.159  0.080  0.62  −0.06  0.138  0.077  4.5  1.04  −0.16  0.184  0.091  0.97  −0.19  0.151  0.076  0.63  −0.06  0.137  0.071  4  1.05  −0.15  0.164  0.087  0.98  −0.20  0.150  0.077  0.67  −0.08  0.149  0.085  3.5  1.02  −0.14  0.145  0.077  0.96  −0.17  0.160  0.088  0.69  −0.09  0.162  0.104  3  1.05  −0.15  0.149  0.073  0.94  −0.16  0.142  0.074  0.72  −0.12  0.157  0.097  2.5  1.00  −0.13  0.132  0.056  0.97  −0.17  0.158  0.083  0.80  −0.19  0.171  0.107  2  1.02  −0.14  0.160  0.076  0.98  −0.17  0.172  0.094  0.85  −0.22  0.155  0.091  All  0.95  −0.06  0.025  0.030  0.90  −0.10  0.024  0.035  0.80  −0.15  0.023  0.047  PGA  1.02  −0.12  0.178  0.082  0.91  −0.15  0.184  0.095  0.56  −0.03  0.170  0.119  View Large We also compared the amplification ratios from fully non-linear analyses with most recent empirical and semi-empirical relationships proposed by Seyhan & Stewart (2014), Sandikkaya et al. (2013) and Kamai et al. (2014), which are implemented in updated NGA2-west GMPEs. Kamai et al. (2014) and Seyhan & Stewart (2014) equations are obtained by using a data set composed of realistic data and simulations results. In order to estimate the non-linear soil response, they have performed equivalent-linear 1-D simulations. Sandikkaya et al.’s (2013) model is completely based on empirical data. The method used for normalizing the Kamai et al.’s (2014) amplification ratios is the same as the one described by Sandikkaya et al. (2013), such that the values obtained for different reference shear wave velocities become comparable. Simulations show more non-linearity in comparison with relationships, therefore, the amplification ratios for higher base ground motions drop faster than these predictions. Dependency of the non-linear soil effects on the bedrock motions is more notable at the 5 Hz frequency. As it is observable in Fig. 10, the level of amplification ratios for the Vs, 30 = 800 m s − 1 is almost constant with a very slight positive slope at large ground motions which can be due to the reduction of soil stiffness at small plastic strains, while non-linear damping is still negligible. Profiles with Vs, 30 = 480 m s − 1 show a somewhat constant amplification ratio at the frequency of 3 Hz and at small to moderate excitations of the bedrock at 5 Hz. Amplification ratios of profiles with low average velocity Vs, 30 < 300 m s − 1 exhibit a descending manner at both frequencies and their independency from the bedrock motion occurs at levels below ∼0.2 g. Results show that for large ground motions, the highest 3 Hz amplifications occur at profiles with medium shear wave velocity (i.e. Vs, 30 ∼ 500 m s − 1), while this happens at low shear wave velocities (i.e. Vs, 30 ∼ 300 m s − 1) for lower motions. The amplification ratio at 5 Hz is approximately the same in medium and loose profiles. The amplification ratios are generally higher than their corresponding values at 5 Hz which agree with the results obtained by Kamai et al. (2014). Our results are more consistent with the Sandikkaya et al.’s (2013) model, however, the general form of all studies match with each other. The agreement between the results and Kamai et al.’s (2014) model is also acceptable for almost all cases, except for the medium shear wave profiles at 5 Hz. The difference between the analyses results and relationships may originate from various reasons such as different non-linear curves implemented for prediction of soil non-linear behaviour, the difference between non-linear and equivalent-linear results and also assumptions made for reference shear wave velocities. Figure 10. View largeDownload slide Comparison of the non-linear amplification ratios with (semi-) empirical results of Kamai et al. (2014), Seyhan & Stewart (2014) and Sandikkaya et al. (2013). Figure 10. View largeDownload slide Comparison of the non-linear amplification ratios with (semi-) empirical results of Kamai et al. (2014), Seyhan & Stewart (2014) and Sandikkaya et al. (2013). 4 CONCLUSIONS In this study, we generated broad-band time histories by combining low-frequency synthetic motions, obtained by dynamic simulations, with stochastic physics based high-frequency motions calculated from wave scattering theory. Since the central and southern parts of the city of Tehran are made up of many old low-rise masonry buildings with high natural frequencies, providing such broad-band ground motions is an essential step to make deterministic low-frequency ground motions applicable to practical engineering tasks for the city. Results obtained from simulations at the bedrock level show good consistency with the empirical values in the sense that average values of ground motion parameters fall within a half standard deviation from mean values of the empirical relationships. As discussed by Hartzell et al. (2002) and Frankel & Stephenson (2000) existence of a deep basin in the vicinity of a reverse fault causes large fault-normal and vertical motions. Large amplitudes are associated with upward shear waves in the basin and formation of surface waves at the edges. In the next step, we investigated local site effects on ground motions by taking into account the soil profiles proposed for the Tehran plain by previous studies. Both linear and non-linear soil responses were calculated for stations along the N-S profile of the city. Hartzell et al. (2002) compared the effect of variations in source models and randomized soil profile parameters by considering 18 input synthetic ground motions (i.e. at the bedrock level) and 20 random soil profiles. They concluded that variations in source models are predominant on variability of the ground motion parameters. Therefore, considering 17 broad-band models for each station in our study may be appropriate for considering slight to moderate changes in soil profile characteristics (Toro 1997). Linear soil response models predicted large unrealistic amplification factors, ranging from 1.5 to more than 3, which do not seem realistic with respect to the empirical relations and observations (Boore & Joyner 1997). This large amplifications resulted in overestimation of the ground motion parameters with respect to empirical relationships at the surface. In contrast, taking into account the non-linear soil behaviour led to plausible results which follow the empirical relationships with one standard deviation. Average values of the horizontal PGA at the surface reach up to 1.09 g on the hangingwall side. On the footwall side, which contains most parts of the city, mean PGA value is about 0.93 g at the northern, 0.59 g at the central and 0.40 g at the southern parts. According to the results, the extent of plasticity development in soil column plays a major role in deamplification due to non-linear soil response. The development of plasticity is not only controlled by the soil properties but also by the ground motion intensity at the base of the column, therefore, soft soil columns with large linear ground motion response experience the highest deamplifications due to the non-linearities. This phenomenon explains why stations very close to the fault show lower amplification factors. It also explains why in some previous earthquakes such as Northridge 1994, main shocks were characterized by lower amplifications with respect to their aftershocks (Field et al. 1997). As a conclusion, since during moderate and small earthquakes soil columns remain elastic, high linear amplifications are not improbable and therefore in the purpose of performing microzonation studies, both the linear and non-linear behaviour of soil profiles should be considered in estimation of the local site effects. The spatial distribution of amplification factors in linear and non-linear cases are quite different, specifically at higher frequencies. Stations with the highest amplifications in linear case, experience the highest deamplifications due to non-linear soil response. Conclusions obtained here place a question on the values reported by Haghshenas (2005) for amplifications of the Tehran alluvium. They reported very large, up to 7, amplification ratios for stations above the Tehran basin. Unrealistically large values in that study are calculated based on microearthquake investigations, hence, do not take into account the soil non-linearities at all. However, it is interesting that our linear simulations also predict overall amplification factors (amplifications caused by deep sedimentary basin up to the bedrock multiplied by near-surface amplifications from bedrock to surface) up to 7.5–8.5 at deepest parts of the basin. Finally, it is notable that although equivalent-linear results underestimate the modest peaks of the motions, they can properly estimate the maximum and large peaks of the motion with very effective reduction in computational time cost. Results demonstrate the need for more accurate studies on properties of the soil profiles of the city. We assigned soil non-linear parameters with respect to the alluvium's structure, composition and loading history. In this regard, more boreholes with larger depths should be investigated over the city to provide a complete database for microzonation studies. Until that time, we proposed a model for three site classes in Tehran to estimate the non-linear soil column response from linear values. In order to propose a more applicable and comprehensive model, we eliminated source-to-site distance effect by considering all time histories for all soil profiles. The correction equation and corresponding parameters for frequencies in the range of 2–10 are proposed in this paper. Broad-band motions calculated here, can be used for estimation of damages and required seismic rehabilitation plans in the city due to the probable M7.2 earthquake of the NTF. The final goal of this research is to propose SA maps that include soil effects for the city of Tehran and also providing a specific synthetic-based GMPE for the Tehran region by using simulation results which will be presented in a forthcoming paper. 4.1 Data and resources The code ‘Broad-Band toolbox’ is developed by Martin Mai and Walter Imperatori at ETH Zurich and can be accessed from (https://ces.kaust.edu.sa/Pages/Software.aspx, last accessed 2017 January). Attenuation curves were plotted by OpenSHA attenuation relationship plotter (http://www.opensha.org, last accessed 2017 January). Non-linear soil analyses are performed by the academic software DEEPSOIL, developed by Youssef Hashash at University of Illinois at Urbana-Champaign and is accessible from http://deepsoil.cee.illinois.edu/. Acknowledgements The authors are grateful to Martin Mai for his kindly providing ‘Broad-Band toolbox’ software. We are also very grateful to Dr Eiichi Fukuyama and two anonymous reviewers for their insightful and constructive comments, which significantly improved the manuscript. H. Zafarani thanks the continuing support of the International Institute of Earthquake Engineering and Seismology during this research, in the framework of the Probabilities of Earthquake Ruptures in Iran (PERSIA) project. REFERENCES Amir-Faryar B., Aggour M.S., 2016. Effect of fibre inclusion on dynamic properties of clay, Geomech. Geoeng. , 11( 2), 104– 113. https://doi.org/10.1080/17486025.2015.1029013 Google Scholar CrossRef Search ADS   Anderson J.G., 1997. Benefits of scenario ground motion maps, Eng. Geol. , 48( 1–2), 43– 57. https://doi.org/10.1016/S0013-7952(97)81913-8 Google Scholar CrossRef Search ADS   Bay J.A., Sasanakul I., 2005. Measurement of nonlinear dynamic properties of bonneville silty-clays: NIW program, Department of Civil and Environmental Engineering, Utah State University, Logan, Utah, Award 04HQGR0055. Berberian M., Yeats R.S., 1999. Patterns of historical earthquake rupture in the Iranian Plateau, Bull. seism. Soc. Am. , 89, 120– 139. Boore D.M., 1983. Stochastic simulation of high-frequency ground motions based on seismological models of the radiated spectra, Bull. seism. Soc. Am. , 73, 1865– 1894. Boore D.M., 2003. Simulation of ground motion using the stochastic method, Pure appl. geophys. , 160( 3), 635– 676. https://doi.org/10.1007/PL00012553 Google Scholar CrossRef Search ADS   Boore D.M., Atkinson G.M., 2008. Ground-motion prediction equations for the average horizontal component of PGA, PGV, and 5%-samped PSA at spectral periods between 0.01 s and 10.0 s, Earthq. Spectra , 24( 1), 99– 138. https://doi.org/10.1193/1.2830434 Google Scholar CrossRef Search ADS   Boore D.M., Joyner W.B., 1997. Site amplifications for generic rock sites, Bull. seism. Soc. Am. , 87, 327– 341. Dalguer L., Mai P., 2012. Prediction of near-source ground motion exceeding 1 g at low frequencies (<2 Hz) from Mw ∼ 6.5 deterministic physics-based dynamic rupture simulations, in Proceedings of the 15th World Conference on Earthquake Engineering (15WCEE), Lisbon, Portugal. Darendeli M.B., 2001. Development of a New Family of Normalized Modulus Reduction and Material Damping Curves, PhD thesis , University of Texas. Dreger D., Tinti E., Cirella A., 2007. In Seismological Society of America 2007 Annual Meeting , 2007 April 11–13, Waikoloa, Hawaii. Emami M.H., Amini B., Jamshidi K.H., Afsharyanzadeh A.M., 1993. Geology Map of Tehran (1:100,000 scale) , Geological Survey of Iran. Farahani J.V., Zaré M., Cichowicz A., 2012. Attenuation of high-frequency P and S waves in south and southeast Tehran using blast data, Soil Dyn. Earthq. Eng. , 40, 99– 108. https://doi.org/10.1016/j.soildyn.2012.03.005 Google Scholar CrossRef Search ADS   Field E.H., Johnson P.A., Beresnev I.A., Zeng Y., 1997. Nonlinear ground-motion amplification by sediments during the 1994 Northridge earthquake, Nature , 390( 6660), 599– 602. https://doi.org/10.1038/37586 Google Scholar CrossRef Search ADS   Frankel A., Stephenson W., 2000. Three-dimensional simulations of ground motions in the Seattle region for earthquakes in the Seattle fault zone, Bull. seism. Soc. Am. , 90( 5), 1251– 1267. https://doi.org/10.1785/0119990159 Google Scholar CrossRef Search ADS   Ghorbani A., Kalantari F., Zohori M., 2014. Estimation of the drained shear strength parameters of fine grained soil of southern Tehran using nonlinear regression model, J. Eng. Geol. , 7, 1869– 1900. Haghshenas E., 2005. Conditions Géotechniques et aléa Sismique local à Téhéran 2005, PhD thesis, Joseph Fourier University (Grenoble I), Grenoble France. Hartzell S., Bonilla L., Williams R.A., 2004. Prediction of nonlinear soil effects, Bull. seism. Soc. Am. , 94( 5), 1609– 1629. https://doi.org/10.1785/012003256 Google Scholar CrossRef Search ADS   Hartzell S., Leeds A., Frankel A., Williams R.A., Odum J., Stephenson W., Silva W., 2002. Simulation of broadband ground motion including nonlinear soil effects for a magnitude 6.5 earthquake on the Seattle fault, Seattle, Washington, Bull. seism. Soc. Am. , 92( 2), 831– 853. https://doi.org/10.1785/0120010114 Google Scholar CrossRef Search ADS   Hashash Y., Groholski D., Phillips C., Park D., Musgrove M., 2012. DEEPSOIL 5.1., User Manual and Tutorial , p. 107, University of Illinois, Urbana, IL, USA. Hashash Y.M., Park D., 2001. Non-linear one-dimensional seismic ground motion propagation in the Mississippi embayment, Eng. Geol. , 62( 1–3), 185– 206. https://doi.org/10.1016/S0013-7952(01)00061-8 Google Scholar CrossRef Search ADS   Hassani B., Zafarani H., Farjoodi J., Ansari A., 2011. Estimation of site amplification, attenuation and source spectra of S-waves in the East-Central Iran, Soil Dyn. Earthq. Eng. , 31( 10), 1397– 1413. https://doi.org/10.1016/j.soildyn.2011.05.017 Google Scholar CrossRef Search ADS   Idriss I., Sun J.I., 1992. SHAKE91: A Computer Program for Conducting Equivalent Linear Seismic Response Analyses of Horizontally Layered Soil Deposits , Center for Geotechnical Modeling, Department of Civil and Environmental Engineering, University of California, Davis, CA. Irikura K., 1986. Prediction of strong acceleration motion using empirical Green's function, in Proc. 7th Japan Earthq. Eng. Symp , pp. 151– 156, IEEE. Jafari M., Razmkhah A., Sohrabi A., Keshavarz M., Pourazin K., 2001. Complementary seismic microzonation studies for South of Tehran, International Institute of Earthquake Engineering and Seismology, Tehran, Iran , p. 126. Jafari M., Saniei H., Haidari M., 2000. Structural classification of clay soils and its application in classifying Tehran City clays, Geosciences , 9, 2– 17. Jafari M.K., Shafiee A., Razmkhah A., 2002. Dynamic properties of fine grained soils in south of Tehran, J. Seismol. Earthq. Eng. , 4, 25– 35. Jica C., 2000. The study on seismic microzoning of the Greater Tehran Area in the Islamic Republic of Iran, Pacific Consultants Int. Rep., OYO Cooperation, Japan . Kamai R., Abrahamson N., Silva W.J., 2014. Nonlinear horizontal site amplification for constraining the NGA-West2 GMPEs, Earthq. Spectra , 30( 3), 1223– 1240. https://doi.org/10.1193/070113EQS187M Google Scholar CrossRef Search ADS   Keshavarz Bakhshayesh M., 2003. The effect of southern Tehran faults on alluvial deposits, in Proceeding 4th International Conference on Sismology and Earthquake Engineering (SEE 4), Tehran 12–14. Kondner R.L., Zelasko J.S., 1964. Void ratio effects on the hyperbolic stress-strain response of a sand, in Laboratory Shear Testing of Soils , ASTM International. Google Scholar CrossRef Search ADS   Liu P., Archuleta R.A., Hartzell S.H., 2006. Prediction of broadband ground-motion time histories: hybrid low/high- frequency method with correlated random source parameters, Bull. seism. Soc. Am. , 96( 6), 2118– 2130. https://doi.org/10.1785/0120060036 Google Scholar CrossRef Search ADS   Mai P.M., Beroza G., 2003. A hybrid method for calculating near-source, broadband seismograms: application to strong motion prediction, Phys. Earth planet. Inter. , 137( 1–4), 183– 199. https://doi.org/10.1016/S0031-9201(03)00014-1 Mai P.M., Imperatori W., Olsen K.B., 2010. Hybrid broadband ground-motion simulations: combining long-period deterministic synthetics with high-frequency multiple S-to-S backscattering, Bull. seism. Soc. Am. , 100( 5A), 2124– 2142. https://doi.org/10.1785/0120080194 Google Scholar CrossRef Search ADS   Majidinejad A., Zafarani H., Vahdani S., 2017. Dynamic simulation of ground motions from scenario earthquakes on the North Tehran Fault, Geophys. J. Int. , 209( 1), 434– 452. https://doi.org/10.1093/gji/ggx017 Google Scholar CrossRef Search ADS   Masing G., 1926. Eigenspannungen und verfestigung beim messing, in Proceedings of the 2nd International Congress of Applied Mechanics , pp. 332– 335, Zürich. Matasović N., Vucetic M., 1995. Generalized cyclic-degradation-pore-pressure generation model for clays, J. Geotech. Eng. , 121( 1), 33– 42. https://doi.org/10.1061/(ASCE)0733-9410(1995)121:1(33) Google Scholar CrossRef Search ADS   Mena B., Mai P.M., Olsen K.B., Purvance M.D., Brune J.N., 2010. Hybrid broadband ground-motion simulation using scattering green's functions: application to large-magnitude events, Bull. seism. Soc. Am. , 100( 5A), 2143– 2162. https://doi.org/10.1785/0120080318 Google Scholar CrossRef Search ADS   Motazedian D., Atkinson G.M. 2005. Stochastic finite-fault modeling based on a dynamic corner frequency, Bull. seism. Soc. Am. , 95( 3), 995– 1010. https://doi.org/10.1785/0120030207 Google Scholar CrossRef Search ADS   Pedrami M., 1981. Pasadenian orogeny and geology of last 700,000 years of Iran, Int. Rep., Geological Survey of Iran (in Persian) , 273. Phillips C., Hashash Y.M., 2009. Damping formulation for nonlinear 1D site response analyses, Soil Dyn. Earthq. Eng. , 29( 7), 1143– 1158. https://doi.org/10.1016/j.soildyn.2009.01.004 Google Scholar CrossRef Search ADS   Rieben H., 1955. The geology of the Teheran plain, Am. J. Sci. , 253( 11), 617– 639. https://doi.org/10.2475/ajs.253.11.617 Google Scholar CrossRef Search ADS   Ritz J.F. et al.  , 2012. Paleoearthquakes of the past 30,000 years along the North Tehran Fault (Iran), J. geophys. Res. , 117, doi:10.1029/2012JB009147. Roten D., Olsen K., Pechmann J., 2012. 3D simulations of M 7 Earthquakes on the wasatch fault, utah, Part II: broadband (0–10 Hz) ground motions and nonlinear soil behavior, Bull. seism. Soc. Am. , 102( 5), 2008– 2030. https://doi.org/10.1785/0120110286 Google Scholar CrossRef Search ADS   Sandikkaya M.A., Akkar S., Bard P., 2013. A nonlinear site-amplification model for the next pan-european ground-motion prediction equations, Bull. seism. Soc. Am. , 103( 1), 19– 32. https://doi.org/10.1785/0120120008 Google Scholar CrossRef Search ADS   Seed H.B., Idriss I.M., 1969. Influence of soil conditions on ground motions during earthquakes, J. Soil Mech. Found. Div. , 95, 99– 138. Seyhan E., Stewart J.P., 2014. Semi-empirical nonlinear site amplification from NGA-West2 data and simulations, Earthq. Spectra , 30( 3), 1241– 1256. https://doi.org/10.1193/063013EQS181M Google Scholar CrossRef Search ADS   Shafiee A., Azadi A., 2007. Shear-wave velocity characteristics of geological units throughout Tehran City, Iran, J. Asian Earth Sci. , 29( 1), 105– 115. https://doi.org/10.1016/j.jseaes.2006.02.005 Google Scholar CrossRef Search ADS   Standard A., 2003. D3999-91. ‘‘Standard test methods for the determination of the modulus and damping properties of soils using the cyclic triaxial apparatus’’, in Annual Book of ASTM Standard . Stewart J.P., Afshari K., Hashash Y.M., 2014. Guidelines for performing hazard-consistent one-dimensional ground response analysis for ground motion prediction, in PEER Rep. , p. 16. Thráinsson H., Kiremidjian A.S., Winterstein S.R., 2000. Modeling of Earthquake Ground Motion in the Frequency Domain , John A. Blume Earthquake Engineering Center. Tiniti E., Fukuyama E., Piatanesi A., Cocco M., 2005. A kinematic source-time function compatible with earthquake dynamics, Bull. seism. Soc. Am. , 95( 4), 1211– 1223. https://doi.org/10.1785/0120040177 Google Scholar CrossRef Search ADS   Toro G., 1997. Probabilistic Models of Shear-Wave Velocity Profiles at the Savannah River Site , Aiken, South Carolina, Risk Engineering, Inc. Tech. Rep., Boulder, CO, p. 72. Towhata I., 2008. Geotechnical earthquake engineering, in Springer Series in Geomechanics and Geoengineering , Springer, Berlin. doi:10.1007/978-3-540-35783-4. Google Scholar CrossRef Search ADS   Vucetic M., Dobry R., 1991. Effect of soil plasticity on cyclic response, J. Geotech. Eng. , 117, 89– 107. https://doi.org/10.1061/(ASCE)0733-9410(1991)117:1(89) Google Scholar CrossRef Search ADS   Zafarani H., Vahidifard H., Ansari A., 2013. Prediction of broadband ground-motion time histories: the case of Tehran, Iran, Earthq. Spectra , 29, 633– 660. https://doi.org/10.1193/1.4000150 Google Scholar CrossRef Search ADS   Zafarani H., Vahidifard H., Ansari A., 2012. Sensitivity of ground-motion scenarios to earthquake source parameters in the Tehran metropolitan area, Iran, Soil Dyn. Earthq. Eng. , 43, 342– 354 https://doi.org/10.1016/j.soildyn.2012.07.007 Google Scholar CrossRef Search ADS   Zeng Y., 1993. Theory of scattered P-and S-wave energy in a random isotropic scattering medium, Bull. seism. Soc. Am. , 83, 1264– 1276. Zeng Y., Anderson J.G., Su F., 1995. Subevent rake and random scattering effects in realistic strong ground motion simulation, Geophys. Res. Lett. , 22, 17– 20. https://doi.org/10.1029/94GL02798 Google Scholar CrossRef Search ADS   Zeng Y., Anderson J.G., Yu G., 1994. A composite source model for computing realistic synthetic strong ground motions, Geophys. Res. Lett. , 21, 725– 728. https://doi.org/10.1029/94GL00367 Google Scholar CrossRef Search ADS   Zeng Y., Su F., Aki K., 1991. Scattering wave energy propagation in a random isotropic scattering medium: 1. Theory, J. geophys. Res. , 96, 607– 619. https://doi.org/10.1029/90JB02012 Google Scholar CrossRef Search ADS   © The Author(s) 2018. Published by Oxford University Press on behalf of The Royal Astronomical Society.

### Journal

Geophysical Journal InternationalOxford University Press

Published: May 1, 2018

## You’re reading a free preview. Subscribe to read the entire article.

### DeepDyve is your personal research library

It’s your single place to instantly
that matters to you.

over 18 million articles from more than
15,000 peer-reviewed journals.

All for just $49/month ### Explore the DeepDyve Library ### Search Query the DeepDyve database, plus search all of PubMed and Google Scholar seamlessly ### Organize Save any article or search result from DeepDyve, PubMed, and Google Scholar... all in one place. ### Access Get unlimited, online access to over 18 million full-text articles from more than 15,000 scientific journals. ### Your journals are on DeepDyve Read from thousands of the leading scholarly journals from SpringerNature, Elsevier, Wiley-Blackwell, Oxford University Press and more. All the latest content is available, no embargo periods. DeepDyve ### Freelancer DeepDyve ### Pro Price FREE$49/month
\$360/year

Save searches from
PubMed

Create lists to

Export lists, citations