University of California, Berkeley



COMPARISON OF VALIDATIONS OF GROUND MOTION SIMULATION PROCEDURES

Report to the PEER-Lifelines Program, Project 1C02d

Nancy Collins, Robert Graves, and Paul Somerville

URS Corporation, Pasadena Office

April 13, 2006

introduction

The idea of simulating broadband strong ground motion time histories is not new, and dates back at least to the pioneering work of Hartzell (1978) and Irikura (1978). These early studies proposed a method of summing recordings of small earthquakes (empirical Green’s functions) to estimate the response of a larger earthquake. Since then, the simulation techniques have been extended to include stochastic representation of source and path effects (e.g., Boore, 1983), theoretical full waveform Green’s functions (e.g., Zeng, 1994), or various combinations of these approaches (e.g., Hartzell, 1989). Over the years, a large number of investigators have made significant contributions and refinements to these methodologies. Hartzell (1999) provides a detailed and comprehensive review of many of these existing simulation methodologies.

VALIDATION EXERCISES

This report compares the results of validations of three ground motion simulation procedures against recorded data. The simulation procedures are those of PEA (Silva et al., 2002), UNR (Zeng and Anderson, 1994; Zeng, 2002), and URS (Graves and Pitarka, 2004). The three simulation procedures analyzed in this report use significantly different approaches to the representation of the earthquake source, seismic wave propagation, and site response. These procedures are outlined in the three appendixes at the end of this report.

The earthquakes and recording stations used for the validation of the ground motion simulation procedures are described in Chiou et al. (2003). The earthquakes include the 1979 Imperial Valley, 1989 Loma Prieta, 1992 Landers, 1994 Northridge, 1995 Kobe, and 1999 Kocaeli earthquakes.

We quantify the goodness of fit between the response spectra of the recorded and simulated ground motions by calculating the residuals between the recorded and simulated values. We analyze the goodness of fit by examining the dependence of these residuals on distance, magnitude, and directivity parameter. The goodness of fit is indicated by the absence of trends in the residuals.

DISTANCE RESIDUALS

Figures 1, 2 and 3 show the residuals of the individual data points for each earthquake as a function of distance, together with a quadratic fit to the residuals for each earthquake and an average of these quadratic fits over the earthquakes, for a set of response spectral periods. Figure 4 shows a summary of the quadratic curves for the three simulation methods. In general, all three methods are unbiased for distances less than 70 km, but underpredict the data at distances larger than about 70 km at periods of 3 seconds and less.

MAGNITUDE RESIDUALS

The magnitude residuals for close distances (less than 20 km) are shown by event for each of the three methods in Figures 5, 6 and 7. The residuals for individual recording stations are shown by open circles, and the median value of the residuals for each even is shown by a solid square. The solid lines are quadratics drawn through the median residual for each event. Figure 8 shows a summary of the quadratic curves for the three simulation methods. The PEA and UNR methods tend to underpredict the magnitude scaling at short periods.

Equivalent results are shown for the whole distance range (0 – 500 km) in Figures 9 through 12. As expected, the magnitude scaling is generally better behaved when averaged over this distance range. Finally, results are shown for the distance range of 70 – 500 km in Figures 13 through 16. The tendency for large distance residuals beyond 70 km, noted in the discussion of distance residuals, is clearly evident in these magnitude residuals.

DIRECTIVITY RESIDUALS

Figures 17, 18 and 19 show the residuals of the individual data points for each earthquake as a function of directivity parameter, together with a quadratic fit to the residuals for each earthquake and an average of these quadratic fits over the earthquakes, in the distance range of 0 – 20 km for a set of response spectral periods. The directivity parameter quantifies the degree to which the ground motion contains forward rupture directivity effects, as formulated by the model of Somerville et al. (1997), ranging from zero (for full backward directivity) to 1.0 (for full forward directivity). Figure 20 shows a summary of the quadratic curves for the three simulation methods. The Somerville et al. (1997) model indicates that directivity effects are systematically present at periods longer than 0.5 seconds, and become larger with increasing period. At a period of 2 seconds, the PEA method underpredicts directivity effects and the URS model overpredicts them. For periods of 4 and 5 seconds, the UNR method overpredicts directivity, and the URS method predicts the expected directivity effect.

CONCLUSIONS

The three simulation procedures analyzed in this report use significantly different procedures for the representation of the earthquake source, seismic wave propagation, and site response. Nevertheless, all three simulation procedures are reasonably successful in modeling recorded ground motions. This gives rise to the expectation that they should produce ground motions having similar scaling characteristics when used in the simulation of ground motions from future earthquakes. However, for this expectation to be realized, it will be necessary to remove differences in the representation of the earthquake source used in the different simulation procedures (Collins et al., 2006).

REFERENCES

Abrahamson, N.A., Somerville, P.G., Cornell, C.A. (1990). "Uncertainty in numerical strong motion predictions" Proc. Fourth U.S. Nat. Conf. Earth. Engin., Palm Springs, CA., 1, 407-416.

Boore D. "Stochastic simulation of high frequency ground motions based on seismological models of the radiated spectra." Bull. Seism. Soc. Am. 1983; 73: 1865-1894.

Borcherdt R. "Estimates of site-dependent response spectra for design (methodology and justification)." Earthquake Spectra 1994; 10(4): 617-653.

Chiou, B. et al. (2003). Validation Guidelines for Numerical Simulation of Ground Motion on Rock Conditions. July 10, 2003.

Collins, N., R. Graves, and P. Somerville (2006). Revised Analysis of 1D Rock Simulations for The NGA-E Project. Report to the Peer-Lifelines Program, Project 1C02d, April 13, 2006.

Graves, R.W. and A. Pitarka (2004). Broadband time history simulation using a hybrid approach. Proceedings of the 13th World Conference on Earthquake Engineering, Vancouver, Canada, August 1-6, 2004, Paper No. 1098.

Irikura K. "Semi-empirical estimation of strong ground motions during large earthquakes." Bull. Disast. Prev. Res. Inst., Kyoto Univ. 1978; 33: 63-104.

Hartzell S. "Earthquake aftershocks as Green's functions." Geophys. Res. Lett. 1978; 5: 1-4.

Hartzell S. "Comparison of seismic waveform inversion results for the rupture history of a finite fault: application to the 1986 North Palm Springs, California, earthquake." J. Geophys. Res. 1989; 94: 7515-7534.

Hartzell S, Harmsen S, Frankel A, Larsen S. "Calculation of broadband time histories of ground motion: comparison of methods and validation using strong ground motion from the 1994 Northridge earthquake." Bull. Seism. Soc. Am. 1999; 89: 1484-1504.

Silva, W., N. Gregor and R. Darragh (2002). Validation of 1-d numerical simulation procedures. Final Report to PEER-Lifelines Project.

Somerville, P.G., N.F. Smith, R.W. Graves, and N.A. Abrahamson (1997). Modification of empirical strong ground motion attenuation relations to include the amplitude and duration effects of rupture directivity, Seismological Research Letters, 68, 180-203.

Somerville, P.G., K. Irikura, R. Graves, S. Sawada, D. Wald, N. Abrahamson, Y. Iwasaki, T. Kagawa, N. Smith and A. Kowada (1999). Characterizing earthquake slip models for the prediction of strong ground motion. Seismological Research Letters, 70, 59-80.

Zeng, Y. (2002). Final Technical Report on Validation of 1-D Numerical Simulation Procedures, Final Technical Report, PEER Project 1C02, Task 1: Earthquake ground motion, Seismological Lab, University of Nevada – Reno.

Zeng Y, Anderson JG, Yu G. "A composite source model for computing synthetic strong ground motions." Geophys. Res. Lett. 1994; 21: 725-728.

Zeng, Y., J. G. Anderson and G. Yu (1994). A composite source model for computing realistic synthetic strong ground motions, J. Res. Lett., 21, 725-728.

[pic]

Figure 1. Residuals between response spectra of data and simulations and quadratic fits to the residuals for each of six earthquakes as a function of distance for the PEA simulation procedure.

[pic]

Figure 2. Residuals between response spectra of data and simulations and quadratic fits to the residuals for each of six earthquakes as a function of distance for the UNR simulation procedure.

[pic]

Figure 3. Residuals between response spectra of data and simulations and quadratic fits to the residuals for each of six earthquakes as a function of distance for the URS simulation procedure.

[pic]

Figure 4. Quadratic fit to residuals between response spectra of data and simulations averaged over six earthquakes as a function of distance for each of three simulation procedures.

[pic]

Figure 5. Residuals between response spectra of data and simulations for each of six earthquakes as a function of magnitude for the PEA simulation procedure for the distance range 0 – 20 km.

[pic]

Figure 6. Residuals between response spectra of data and simulations for each of six earthquakes as a function of magnitude for the UNR simulation procedure for the distance range 0 – 20 km.

[pic]

Figure 7. Residuals between response spectra of data and simulations for each of six earthquakes as a function of magnitude for the URS simulation procedure for the distance range 0 – 20 km.

[pic]

Figure 8. Quadratic fit to residuals between response spectra of data and simulations averaged over six earthquakes as a function of magnitude for each of three simulation procedures for the distance range 0 – 20 km.

[pic]

Figure 9. Residuals between response spectra of data and simulations for each of six earthquakes as a function of magnitude for the PEA simulation procedure for the distance range 0 – 500 km.

[pic]

Figure 10. Residuals between response spectra of data and simulations for each of six earthquakes as a function of magnitude for the UNR simulation procedure for the distance range 0 – 500 km.

[pic]

Figure 11. Residuals between response spectra of data and simulations for each of six earthquakes as a function of magnitude for the URS simulation procedure for the distance range 0 – 500 km.

[pic]

Figure 12. Quadratic fit to residuals between response spectra of data and simulations averaged over six earthquakes as a function of magnitude for each of three simulation procedures for the distance range 0 – 500 km.

[pic]

Figure 13. Residuals between response spectra of data and simulations for each of six earthquakes as a function of magnitude for the PEA simulation procedure for the distance range 70 – 500 km.

[pic]

Figure 14. Residuals between response spectra of data and simulations for each of six earthquakes as a function of magnitude for the UNR simulation procedure for the distance range 70 – 500 km.

[pic]

Figure 15. Residuals between response spectra of data and simulations for each of six earthquakes as a function of magnitude for the URS simulation procedure for the distance range 70 – 500 km.

[pic]

Figure 16. Quadratic fit to residuals between response spectra of data and simulations averaged over six earthquakes as a function of magnitude for each of three simulation procedures for the distance range 70 -500 km.

[pic]

Figure 17. Residuals between response spectra of data and simulations and quadratic fits to the residuals for each of six earthquakes as a function of the directivity parameter for the PEA simulation procedure for distances of 0 – 20 km.

[pic]

Figure 18. Residuals between response spectra of data and simulations and quadratic fits to the residuals for each of six earthquakes as a function of the directivity parameter for the UNR simulation procedure for distances of 0 – 20 km.

[pic]

Figure 19. Residuals between response spectra of data and simulations and quadratic fits to the residuals for each of six earthquakes as a function of the directivity parameter for the URS simulation procedure for distances of 0 – 20 km.

[pic]

Figure 20. Quadratic fit to residuals between response spectra of data and simulations averaged over six earthquakes as a function of directivity parameter for each of three simulation procedures for the distance range 0 – 20 km.

APPENDIX I. STOCHASTIC FINITE-FAULT MODEL (PEA)

The simulation procedure implemented here is termed the stochastic finite fault model (Silva et al., 1990; 1997) in which near surface (top 100 to 1,000 ft) materials are considered in detail through an equivalent-linear formulation. The stochastic finite-fault model implemented here is quite simple in concept, using a single-corner-frequency omega-square source spectrum (M = 5.0) for each subfault. Large earthquakes are simulated by simply delaying and summing contributions from the M 5 subfaults. The model is simple, includes a frequency domain random vibration theory equivalent-linear site response implemented for both rock and soil sites (Silva et al., 1997). The model, including site effects, has recently been validated at about 500 sites for 15 earthquakes (M 5.2 to 7.4) over fault distances ranging from 1 km to 470 km (Silva et al., 1997) and for subduction zone earthquakes for M up to 8.1. In general, the model is unbiased over the frequency range of recorded motions (spectral acceleration averaging from about 0.3 to 100 Hz).

To model the effects of the shallow geotechnical layer on the predicted motions, generic rock and deep firm soil profiles were placed on top of the regional models (Table 1). Small strain kappa values were set to 0.04 sec, base case values for both deep soil and soft rock implied by the Abrahamson and Silva (1997) and Sadigh et al. (1997) empirical attenuation relations (Silva et al., 1997; EPRI, 1993). For the Kocaeli and Duzce, Turkey earthquakes, many of the recording sites had shallow measured shear-wave velocity profiles, greatly facilitating site classification and the development of realistic generic profiles to basement depths. The measured profiles generally fell into soft rock and deep firm soils, very similar to Northern California surficial geology based profiles for Franciscan (Fr) rock and Quaternary alluvium (Qal) (Silva et al., 1998). As a result, the Northern California Franciscan and Quaternary alluvium generic profiles were adopted to model the geotechnical layer for the Turkey earthquakes recording sites. For both the Taiwan and Turkey sites, Northern California nonlinear dynamic material properties (Silva et al., 1998; EPRI, 1993) were assumed. Final assessment of the most appropriate G/Gmax and hysteretic damping curves requires refinement of source parameters (slip model and rise time) and better characterization of the Chi-Chi, Taiwan earthquake site profiles.

Slip models and nucleation points used for all three earthquakes were taken from Somerville (personal communication). The slip models are shown in Figure 3 and were developed using inversions which accommodate both spatially varying rise times as well as rupture velocities. Since the stochastic finite fault model implemented here uses a fixed rise time and rupture velocity as well as a simpler approach to wave propagation, some inconsistencies likely exist between these slip models and ones optimized for the stochastic model.

For the Kocaeli, Turkey earthquake, the bias is near zero when averaged over all distances, but shows an underprediction at long periods, exceeding about 3 seconds. For the closer sites (within a 50 km rupture distance) there is a general overprediction at periods shorter than about 3 seconds. The uncertainty computed over all the sites is higher than typical (Silva et al., 1997), ranging from about 0.75 at peak acceleration to about 1 at 10 seconds. Inversions for more appropriate slip models, rupture velocities, and rise times are required to resolve the over prediction for the sites within a 50 km rupture distance.

The general procedure followed when the model was validated using a large suite of earthquakes (Silva et al., 1997) was to find the best fitting rupture velocities, rise times, and subevent stress drop as well as Q(f) models, small strain kappa values, and nonlinear dynamic material properties (Silva et al., 1998). Q(f) and kappa models are first estimated using point-source inversions (for stress drop, Q(f), and kappa values). The Q(f) and small stress initial kappa models are then used in finite fault simulations where rise time, rupture velocity, and subevent stress drop are varied for each earthquake to minimize the bias. Average values over all earthquakes are then estimated (a rise time verses moment relation developed) and the validations redone with the “global” average values for rupture velocity, rise time, and subevent stress drop. This approach results in a larger uncertainty over all earthquakes but eliminates the need for developing and defending parametric distributions for these parameters for the next earthquake. This process should be repeated (updated) as new earthquakes are added to the validation set. Once reliable site conditions and generic profiles are available for the Chi-Chi, Taiwan earthquake and it is considered to be representative of California earthquakes and the near source sites for the Duzce, Turkey earthquake are judged to be useable, this update process will be implemented for the stochastic finite source model.

References

Abrahamson, N.A. and Silva, W.J. (1997). "Empirical response spectral attenuation relations for shallow crustal earthquakes." Seism. Res. Lett., 68(1), 94-127.

Abrahamson, N.A., Somerville, P.G., Cornell, C.A. (1990). "Uncertainty in numerical strong motion predictions" Proc. Fourth U.S. Nat. Conf. Earth. Engin., Palm Springs, CA., 1, 407-416.

Bouchon, M., M-P. Boui, H. Karabulut, M. N. Toksoz, M. Dietrich, and A. J. Rosakis (2001). “How fast is rupture during an earthquake? New insights from the 1999 Turkey earthquakes.” Geophys. Res. Let., 28(14) 2723-2726.

Boore, D. M. (2001). “Comparisons of ground motions from the 1999 Chi-Chi earthquake with empirical predictions largely based on data from California.” Bull. Seism. Soc. Am., 91(5), 1212-1217.

Electric Power Research Institute (1993). "Guidelines for determining design basis ground motions." Palo Alto, Calif: Electric Power Research Institute, vol. 1-5, EPRI TR-102293. vol. 1: Methodology and guidelines for estimating earthquake ground motion in eastern North America.

vol. 2: Appendices for ground motion estimation.

vol. 3: Appendices for field investigations.

vol. 4: Appendices for laboratory investigations.

vol. 5: Quantification of seismic source effects.

Lee, C.-T., C.-T., Cheng, C.-W., Liao, and Y.-B., Tsai (2001). “Site classification of Taiwan free-field strong-motion stations”. Bull. Seism. Soc. Am., 91(5), 1283-1297.

Neugebauer, J., M. Loeffler, H. Berckhemer, and A. Yatman (1997). “Seismic observations at an overstep of the western North Anatolian Fault (Abant-Sapanca region, Turkey). Geol. Rundsch., 86, 93-102.

Roecker, S.W., Y. H. Yeh, and Y. B. Tsai (1987). “Three-dimensional P and S wave velocity structures beneath Taiwan: Deep structure beneath an ARC-Continent collision”. J. of Geoph. Research, 92(B10), 10,547-10,570.

Sadigh, C.-Y. Chang, J.A. Egan, F. Makdisi, and R.R. Youngs (1997). "Attenuation relationships for shallow crustal earthquakes based on California strong motion data." Seism. Soc. Am., 68(1), 180-189.

Silva, W.J. Costantino, C. Li, Sylvia (1998). AQuantification of nonlinear soil response for the Loma Prieta, Northridge, and Imperial Valley California earthquakes.@ The effects of Surface Geology on Seismic Motion, Irikura, Kudo, Okada & Sasatani (eds.).

Silva, W.J., N. Abrahamson, G. Toro and C. Costantino. (1997). "Description and validation of the stochastic ground motion model." Report Submitted to Brookhaven National Laboratory, Associated Universities, Inc. Upton, New York 11973, Contract No. 770573.

Sokolov, V. Y., C.-H. Loh, K.-L. Wen (2001). “Empirical models for site- and region-dependent ground-motion parameters in the Taipei Area: A unified approach” Earthquake Spectra, 17(2), 313-333.

APPENDIX II. COMPOSITE SOURCE MODEL (UNR)

Over the past years, we have focused efforts to develop and improve a numerical simulation procedure to compute synthetic strong motion seismogram using a composite source model (Zeng et al., 1994). The method has been successful in generating realistic strong motion seismograms. The realism is demonstrated by comparing synthetic strong motions with observations from the recent California earthquakes at Landers, Loma Prieta (Su et al., 1994a,b) and Northridge (Zeng and Anderson, 1996; Anderson and Yu, 1996; Su et al., 1998), earthquakes in the eastern US (Ni et al., 1999) and earthquakes in Guerrero, Mexico (Zeng et al., 1994; Johnson, 1999), Turkey (Anderson et al., 1997) and India (Khattri et al, 1994; Zeng et al, 1995). We have also successfully applied the method for earthquake engineering applications to compute the ground motion of scenario earthquakes. During the process of continuing development, we have included scattering waves from small scale heterogeneity structure of the earth, site specific ground motion prediction using weak motion site amplification, and nonlinear soil response using the geotechnical engineering model. We have evaluated the numerical procedure for simulating near-fault long-period ground motions and rupture directivity, revisiting some of the above earthquake events, including Loma Prieta, Landers and Northridge. We also tested its ability to predict the near-fault ground motion observation from the 1979 Imperial Valley, California earthquake and the 1995 Kobe event (Zeng and Anderson, 2000).

We used the Green's function representation theorem to compute synthetic strong motion seismograms. For a given slip distribution, the synthetics can be obtained by taking a spatial and temporal convolution of the Green's function with the fault slip function (Aki and Richards, 1980). We modeled the source slip processes using the composite source model (Zeng et al., 1994). This model assumes a large earthquake is a superposition of smaller subevents that all break during the earthquake rupture processes. The number and radius of the subevents follow the Guttenberg and Richter frequency-magnitude relation given in form of a power law distribution of radii, [pic] where p is the fractal dimension. Given the source parameters of a large earthquake, we can numerically generate the subevents following the power law relation mentioned earlier. We then place these events within the fault plane and allow them to overlap. The random nature of the heterogeneities on a complex fault is achieved by distributing the subevents randomly on the fault plane. Rupture propagates from the hypocenter, and each subevent radiates a displacement pulse of a rupture crack (alternative models are by Brune, 1970, Sato and Hirasawa, 1973, or Zeng, 2001) when the rupture front reaches the subevent. The asymmetrical subevent rupture model of Zeng (2001) is the best of these three source functions to simulate the effect of rupture directivity as characterized by Somerville et al. (1997).

Once the source has been specified, we can propagate the motion generated at the source to the site using layered crustal model (Luco and Apsel, 1983) or 3-D inhomogeneity structure using finite difference model. At high frequencies, contribution from scattering becomes important. It extends the ground motion duration, redistributes energy among different component of seismic waves, and contributes to the spatial incoherence of ground shaking. We modeled the scattering waves as signal generated noise with its envelope functions consistent with the prediction based on the multiple scattering theory of Zeng et al. (1991) and then add these scattering waves to the Green's function computed deterministically using generalized reflection/transmission method (Zeng et al., 1995, Zeng, 1995). The solution is then convolved with a plane wave propagation function through a near surface 1-D velocity layering as complex as that suggested by sonic well logs.

The Kocaeli, Turkey and Chi-Chi, Taiwan earthquakes were the two largest events in 1999. Both events were well recorded by dense strong motion arrays, especially the Chi-Chi earthquake. The wealth of ground motion data from these two earthquakes has provided us the opportunity to understand the source rupture processes of large earthquakes and their high-frequency radiation. The composite source model was used to study the earthquake sources and to predict the high frequency near-field strong ground motion. The genetic algorithm procedure developed by Zeng and Anderson (1996) was used to compute a joint GPS static deformation and strong motion waveform inversion for the Kocaeli and the Chi-Chi, earthquake to determine their composite source rupture processes. For both earthquakes, a fault model with a curved fault plane based on the surface rupture data from field mapping was used. The total moment is 2.1x1027 dyne-cm for the Kocaeli earthquake and 2.9x1027 dyne-cm for the Chi-Chi earthquake.

We found evidence of supershear for the Kocaeli earthquake fault rupture. For both events, the subevents of the composite sources show significant slow rupture or long rise time compared with other earthquakes we have previously studied. By extrapolating the inversion solution to higher frequency waveform simulation, we found while slow subevent rupture velocities fit the low frequency observations, higher subevent rupture velocities are required to fit the high frequency ground accelerations. This indicates that after the rupture front passes along the rupture plane, the fault continues slipping slowly and then stops abruptly. Station site effects were considered by propagating the synthetics through a 30 meters low velocity soil layer on top of the original velocity models using DESRA2 (Lee and Finn, 1982). Our final ground motion simulations show unbiased results through the entire frequency band in comparison with the observed acceleration response spectra for both earthquakes. The average errors of the logarithmic spectral acceleration between the synthetics and the observations are 0.52 and 0.55 for the horizontal and vertical components, respectively for the Kocaeli earthquake. For the Chi-Chi earthquake, the errors are 0.60 and 0.65 for the horizontal and vertical components, respectively. In general, the synthetic accelerations agree with the observation in their peak values and duration.

References

Aki, K., and P. G. Richards (1980). Quantitative Seismology Theory and Methods. W. H. Freeman and Co., San Francisco.

Anderson, J. G. and G. Yu (1996). Predictability of strong motions from the Northridge, California, earthquake, Bull. Seism. Soc. Am. No. 86, 1B, S100-S114.

Anderson, J. G., Y. Zeng, H. Sucuoglu (1997). Analysis accelerations from the Dinar, Turkey earthquake, presented at the Eighth International Conference on Soil Dynamics and Earthquake Engineering at Istanbul, Turkey.

Brune, J. N. (1970). Tectonic stress and the spectra of seismic shear waves from earthquake. J. Geophys. Res. 75, 4997-5009.

Johnson, M. (1999). Composite source model parameters for large earthquakes (M>5.0) in the Mexican subduction zone, M. S. thesis, University of Nevada, Reno.

Khattri, K. N., Y. Guang, J. G. Anderson, J. N. Brune and Y. Zeng (1994). Seismic hazard estimation using modelling of earthquake strong ground motions: A brief analysis of 1991 Uttarkashi earthquake, Himalaya and prognostication for a great earthquake in the region, Current Science, 67, 343-353

Lee, M. K. W. and W. D. L. Finn (1982). Dynamic effective stress response analysis of soil deposites with energy transmitting boundary including assessment of liquefaction potential, Rev., Dept. of Civil Eng., Soil Mechanics Serirs No. 38, the Univ. of British Columbia, Vancouver, Canada.

Luco, J. E. and R. J. Apsel (1983). On the Green's function for a layered half-space, part I, Bull. Seism. Soc. Am. 73, 909-929.

Ni., S.-D., J. G. Anderson, Y. Zeng (1999). Comparison of strong ground motions from the 1988 Seguenay earthquake with the synthetic simulations using the composite source model, in preparation.

Sato, T. and T. Hirasawa (1973). Body wave spectra from propagating shear cracks, J. of Phys. of the Earth, 21, 415-431.

Somerville, P. G., N. F. Smith, R. W. Graves, and N. A. Abrahamson (1997). Modification of empirical strong ground motion attenuation relations to include the amplitude and duration effects of rupture Directivity, Seism. Res. Lett., 68, 199-222.

Su, Feng, Y. Zeng and J. G. Anderson (1994a). Simulation of the Loma Prieta earthquake strong ground motion using a composite source model, EOS, Trans. A.G.U., 75, 44, p448.

Su, F., Y. Zeng and J. G. Anderson (1994b). Simulation of Landers earthquake strong ground motion using a composite source model, Seism. Res. Lett., 65, p52.

Su, F., J. G. Anderson and Y. Zeng (1998). Study of weak and strong motion including nonlinearity in the Northridge, California, earthquake sequence, submitted to Bull. Seis. Soc. Am. 88, 1411-1425.

Yu, G. (1994). Ph.D. Thesis, University of Nevada, Reno.

Zeng, Y., F. Su and K. Aki (1991). Scattered Wave Energy Propagation in Random Isotropic Scattering Medium - Part I: Theory, J. Geophy. Res. 96, 607-619.

Zeng, Y., J. G. Anderson and G. Yu (1994). A composite source model for computing realistic synthetic strong ground motions, J. Res. Lett., 21, 725-728.

Zeng, Y., J. G. Anderson and Feng Su (1995). Subevent rake and random scattering effects in realistic strong ground motion simulation, Geophy. Res. Lett, 22, 17-20.

Zeng, Y. (1995). A realistic synthetic Green's function calculation using a combined deterministic and stochastic modeling approach, EOS, Trans. A.G.U., 76, F357.

Zeng Y. and J. G. Anderson (1996). A composite source modeling of the 1994 Northridge earthquake using Genetic Algorithm, Bull. Seism. Soc. Am. 86, 71-83.

Zeng, Y. and J. G. Anderson (2000). Earthquake source and near-field directivity modeling of several large earthquakes, EERI Proceedings for the Sixth International Conference on Seismic Zonation.

Zeng, Y. (2001). Validation and modeling of earthquake strong ground motion using a composite source model, EOS, Trans. A.G.U., 82, F869.

APPENDIX III. HYBRID BROADBAND GROUND MOTION SIMULATION PROCEDURE (URS)

The broadband ground motion simulation procedure is a hybrid technique that computes the low frequency and high frequency ranges separately and then combines the two to produce a single time history (Hartzell, 1999). At frequencies below 1 Hz, the methodology is deterministic and contains a theoretically rigorous representation of fault rupture and wave propagation effects, and attempts to reproduce recorded ground motion waveforms and amplitudes. At frequencies above 1 Hz, it uses a stochastic representation of source radiation, which is combined with a simplified theoretical representation of wave propagation and scattering effects. The use of different simulation approaches for the different frequency bands results from the seismological observation that source radiation and wave propagation effects tend to become stochastic at frequencies of about 1 Hz and higher.

Our methodology offers a significant enhancement over previous broadband simulation techniques through the use of frequency-dependent non-linear site amplification factors. These factors are incorporated by first restricting the computational velocity model in both the deterministic and stochastic bandwidths to have an average near-surface shear wave velocity between 600 and 1000 m/s. We then apply site-specific amplification factors, which are derived using the empirical relations of Borcherdt (1994). This approach significantly reduces the numerical computational burden, particularly for the deterministic domain, and also provides an efficient mechanism for including detailed site-specific geologic information in the ground motion estimates.

In the sections that follow, we first provide detailed descriptions of the deterministic and stochastic simulation methodologies. Next, we discuss the derivation and implementation of the non-linear site amplification factors.

Determinstic Simulation Methodology (f < 1 Hz)

The low frequency simulation methodology uses a deterministic representation of source and wave propagation effects and is based on the approach described by Hartzell (1983). The basic calculation is carried out using a 3D viscoelastic finite-difference algorithm, which incorporates both complex source rupture as well as wave propagation effects within arbitrarily heterogeneous 3D geologic structure. The details of the finite-difference methodology are described by Graves (1996) and Pitarka (1998). Anealsticity is incorporated using the coarse-grain approach of Day (2001).

The earthquake source is specified by a kinematic description of fault rupture, incorporating spatial heterogeneity in slip, rupture velocity and rise time. Following Hartzell (1983), the fault is divided into a number of subfaults. The slip and rise time are constant across each individual subfault, although these parameters are allowed to vary from subfault to subfault. We use a slip velocity function that is constructed using two triangles as shown in Figure 1. This functional form is based on results of dynamic rupture simulations (e.g., Guaterri, 2003). We constrain the parameters of this function as follows:

[pic] (1)

where M0 is the seismic moment, Tr is the rise time and A is normalized to give the desired final slip.

The expression for Tr comes from the empirical analysis of Somerville (1999). In general, Tr may vary across the fault; however, in practice we only allow a depth dependent scaling such that Tr increases by a factor of 2 if the rupture is between 0 and 5 km depth. This is consistent with observations of low slip velocity on shallow fault ruptures (Kagawa, 2001).

The rupture initiation time (Ti) is determined using the expression

[pic] (2)

where R is the rupture path length from the hypocenter to a given point on the fault surface, Vr is the rupture velocity and is set at 80% of the local shear wave velocity (Vs), and (t is a timing perturbation that scales linearly with slip amplitude such that [pic] where the slip is at its maximum and [pic] where the slip is at the average slip value. We typically set [pic] This scaling results in faster rupture across portions of the fault having large slip as suggested by source inversions of past earthquakes (Hisada, 2001).

For scenario earthquakes, the slip distribution can be specified using randomized spatial fields, constrained to fit certain wave number properties (e.g., Somerville, 1999; Mai, 2002). In the simulation of past earthquakes, we use smooth representations of the static slip distribution determined from finite-fault source inversions. Typically, these inversions will also include detailed information on the spatial variation of rupture initiation time and slip velocity function, either by solving for these parameters directly or by using multiple time windows. However, we do not include these in our simulations, but rather rely on equations (1) and (2) to provide them. Our philosophy is that the level of detailed resolution of these parameters provided by the source inversions will generally not be available a priori for future earthquakes. Furthermore, since the inversions determine these parameters by optimally fitting the selected observations, there are no guarantees that they will produce an optimal waveform fit at sites not used in the inversion. Hopefully, an improved understanding of dynamic rupture processes will help to provide better constraints on these parameters in the future.

Stochastic Simulation Methodology (f > 1 Hz)

The high frequency simulation methodology is a stochastic approach that sums the response for each subfault assuming a random phase, an omega-squared source spectrum and simplified Green’s functions. The methodology follows from the procedure that was first presented by Boore (1983) with the extension to finite-faults given by Beresnev (1997). We have incorporated several modifications of the original finite-fault methodology of Beresnev (1997), which are described below.

In our approach, each subfault is allowed to rupture with a subfault moment weighting that is proportional to the final static slip amount given by the prescribed rupture model. The final summed moment release is then scaled to the prescribed target mainshock moment. This alleviates the problem of requiring that each of the subfaults scale to an integer multiple of [pic] (where [pic] is the stress parameter and [pic] is the subfault dimension), which tends to make many of the subfaults have zero moment release. The subfault dimensions are determined using the scaling relation of Beresnev (2001).

Beresnev (1998) define a radiation-strength factor (s), which is used as a free parameter in the specification of the subfault corner frequency (fc)

[pic] (3)

where z is a scaling factor relating fc to the rise time of the subfault source. In our approach, instead of allowing this to be a free parameter, we set [pic] and let

[pic] [pic] (4)

where Df is a depth scaling factor, [pic]km, [pic]km and h is the depth of the subfault center in km, and Af is a dip scaling factor, [pic], [pic] and δ is the subfault dip. The constants c0 and c1 are set at 0.4 and 0.35, respectively, based on calibration experiments. This parameterization follows from the observation in crustal earthquakes that slip velocity is relatively low for shallow near-vertical ruptures and increases with increasing rupture depth and decreasing fault dip (Kagawa, 2001). Since corner frequency is proportional to slip velocity, this formulation replicates the trend of the observations. We note that although this formulation reduces the number of free parameters, it certainly is not unique and probably has trade-offs with other parameters in the stochastic model. In particular, allowing the subfault stress parameter (σp) to be variable across the fault would accommodate a similar type of slip velocity scaling. Instead, we fix [pic]in our simulations.

Our formulation also allows the specification of a plane layered velocity model from which we calculate simplified Green’s functions (GFs) and impedance effects. The GFs are comprised of the direct and Moho-reflected rays, which are traced through the specified velocity structure. Following Ou (1990), each ray is attenuated by 1/Rp where Rp is the path length traveled by the particular ray. For each ray and each subfault, we calculate a radiation pattern coefficient by averaging over a range of slip mechanisms and take-off angles, varying [pic] about their theoretical values. Anelasticity is incorporated using a travel-time weighted average of the Q values for each of the velocity layers and using a kappa operator set at [pic]. Finally, gross impedance effects are included using quarter wavelength theory (Boore, 1997) to derive amplification functions that are consistent with the specified velocity structure.

Site Specific Amplification Factors

Borcherdt (1994) derived empirically based amplification functions for use in converting response spectra from one site condition to a different site condition. The general form of these functions is given by

[pic] (5)

where Vsite is the 30 m travel-time averaged shear wave velocity (Vs30) at the site of interest, Vref is the corresponding velocity measure at a reference site where the ground response is known, and mx is an empirically determined factor. Borcherdt (1994) specified one set of factors at short periods (centered around 0.3 s) and one set at mid-periods (centered around 1.0 s). Furthermore, non-linear effects are included since the mx decrease as the reference ground response PGA increases. The decrease in the mx is sharper for the short period factors than the mid-period ones, reflecting the observed increase of non-linear effects at shorter periods.

In our simulation methodology, we restrict the computational velocity models in both the deterministic and stochastic calculations to have Vs30 values between 600 and 1000 m/s. This is our Vref. To obtain an amplification function for a given site velocity, we first determine the short- and mid-period factors from equation (5) using the tabulated mx from Borcherdt (1994) given the reference PGA from the stochastic response. Next, we construct a smoothly varying function in the frequency domain by applying a simple taper to interpolate the factors between the short- and mid-period bands. The function tapers back to unity at very short and very long periods. An example set of these functions is shown in Figure 2

In practice, we apply these amplification functions to the amplitude spectra of the Fourier transformed simulated time histories. This process is done to the deterministic and stochastic results separately since these may have different computational reference velocities. Although the amplification factors are strictly defined for response spectra, the application in the Fourier domain appears to be justified since the functions vary slowly with frequency. Finally, the individual responses are combined into broadband response using a set of matched butterworth filters. The filters are 4th-order and zero-phase with a lowpass corner at 1 Hz for the deterministic response and a highpass corner at 1 Hz for the stochastic response. The important properties of the matched filters are 1) they do not alter the phase of the response and 2) they sum to unity for all frequencies. After applying the filters to the individual responses, they are summed together to produce a single broadband time history.

Application to Validation

In order to test the adequacy of our simulation methodology, we compare our computed synthetic strong motion time histories with those recorded during past earthquakes. The only earthquake specific parameters we use are seismic moment, overall fault dimensions and geometry, hypocenter location, and a generalized model of the final slip distribution. For future earthquakes, these are the parameters that we feel can either be reliably estimated (e.g., seismic moment, fault dimensions) or parametrically assessed using multiple realizations (e.g., hypocenter location, slip distribution). All other source parameters are determined using the scaling relations described in the previous section. Since we have not optimized the rupture models for these exercises, we cannot hope to match all the details of the recorded waveforms. However, our goal is to reproduce the overall characteristics of the observed motions over a broad frequency range throughout the region surrounding the fault. This includes matching the trends and levels of common ground motion parameters such as PGA, PGV, SA and duration of shaking, adequately capturing near-fault phenomena such as rupture directivity and footwall / hanging wall effects, and reproducing region or site specific effects such as basin response and site amplification.

The broadband simulation methodology presented here provides a general framework for synthesizing ground motion time histories for future scenario earthquakes. One of the main enhancements of our approach over previous techniques is the use of frequency-dependent non-linear site amplification factors. Our methodology produces quite favorable results when compared against the strong ground motions recorded during the 1989 Loma Prieta and 1994 Northridge earthquakes.

In developing this methodology, we have tried to incorporate as much detail as possible in describing the source, path and site effects in order to adequately capture the main characteristics of the expected ground motions. For the path and site effects, this stresses the importance of developing detailed 3D seismic velocity models for earthquake prone regions. However, we recognize that extremely detailed descriptions of the earthquake rupture process will generally not be available a priori for future events. Thus, our methodology uses simple, yet flexible, rules to parameterize the slip, slip velocity function and rupture velocity. More robust constraints on these parameters from detailed source inversion studies and dynamic rupture analyses will be used in future enhancements of the methodology.

References

Beresnev I, Atkinson G. "Modeling finite fault radiation from the ωn spectrum." Bull. Seism. Soc. Am. 1997; 87: 67-84.

Beresnev I, Atkinson G. "Subevent structure of large earthquakes - A ground motion perspective." Geophys. Res. Lett. 2001; 28(1): 53-56.

Beresnev I, Atkinson G. "Stochastic finite fault modeling of ground motions from the 1994 Northridge, California earthquake. I. validation on rock sites." Bull. Seism. Soc. Am. 1998; 88: 1392-1401.Borcherdt R. "Estimates of site-dependent response spectra for design (methodology and justification)." Earthquake Spectra 1994; 10(4): 617-653.

Boore D. "Stochastic simulation of high frequency ground motions based on seismological models of the radiated spectra." Bull. Seism. Soc. Am. 1983; 73: 1865-1894.

Boore D, Joyner W, "Site amplification for generic rock sites." Bull. Seism. Soc. Am. 1997; 87: 327-341.

Brocher T, Brabb R, Catchings R, Fuis G, Fumal T, Jachens R, Jayko A, Kayen R, McLaughlin R, Parsons T, Rymer M, Stanley R, Wentworth C. "A crustal-scale 3D seismic velocity model for the San Francisco Bay area, California." Eos 1997; 78; F435-F436.

Day S, Bradley C. "Memory efficient simulation of anelastic wave propagation." Bull. Seism. Soc. Am. 2001; 91: 520-531.

Graves R. "Simulating seismic wave propagation in 3D elastic media using staggered grid finite differences." Bull. Seism. Soc. Am. 1996; 86: 1091-1106.

Guatteri M, Mai PM, Beroza G, Boatwright J. "Strong ground motion prediction from stochastic-dynamic source models." Bull. Seism. Soc. Am. 2003; 93: 301-313.

Hartzell S. "Earthquake aftershocks as Green's functions." Geophys. Res. Lett. 1978; 5: 1-4.

Hartzell S. "Comparison of seismic waveform inversion results for the rupture history of a finite fault: application to the 1986 North Palm Springs, California, earthquake." J. Geophys. Res. 1989; 94: 7515-7534.

Hartzell S, Harmsen S, Frankel A, Larsen S. "Calculation of broadband time histories of ground motion: comparison of methods and validation using strong ground motion from the 1994 Northridge earthquake." Bull. Seism. Soc. Am. 1999; 89: 1484-1504.

Hartzell S, Heaton T. "Inversion of strong ground motion and teleseismic waveform data for the fault rupture history of the 1979 Imperial Valley, California earthquake." Bull. Seism. Soc. Am. 1983; 73: 1553-1583.

Hartzell S, Liu P, Mendoza C. "The 1994 Northridge, California earthquake: Investigation of rupture velocity, risetime, and high frequency radiation." J. Geophys. Res. 1996; 101; 20,091-20,108. Irikura K. "Semi-empirical estimation of strong ground motions during large earthquakes." Bull. Disast. Prev. Res. Inst., Kyoto Univ. 1978; 33: 63-104.

Hisada Y. "A theoretical omega-square model considering spatial variation in slip and rupture velocity. part 2. case for a two-dimensional source model." Bull. Seism. Soc. Am. 2001; 91: 651-666.

Kagawa T, Irikura K, Somerville P. "A study on ground motion and fault rupture due to subsurface faults." Eos. Trans. AGU 2001; 82(47): Fall Meet. Suppl., Abstract S31C-06.

Magistrale H, Day S, Clayton R, Graves R. "The SCEC Southern California reference three-dimensional seismic velocity model version 2." Bull. Seism. Soc. Am. 2000; 90(6B): S77-S94.

Mai P, Beroza G. "A spatial random field model to characterize complexity in earthquake slip." J. Geophys. Res. 2002; 107(B11): doi:10.1029/2001JB000588.

Ou GB, Herrmann R. "A statistical model for ground motion produced by earthquakes at local and regional distances." Bull. Seism. Soc. Am. 1990; 80: 1397-1417.

Pitarka A. "3D elastic finite difference modeling of seismic wave propagation using staggered grid with non-uniform spacing." Bull. Seism. Soc. Am. 1998; 88: 54-68.

Schneider J, Silva W, Stark C. "Ground motion model for the 1989 M6.9 Loma Prieta earthquake including effects of source, path and site." Earthquake Spectra 1993; 9: 251-287.

Somerville P, Irikura K, Graves R, Sawada S, Wald D, Abrahamson N, Iwasaki Y, Kagawa T, Smith N, Kowada A. "Characterizing crustal earthquake slip models for the prediction of strong ground motion." Seism. Res. Lett. 1999; 70: 59-80.

Stidham C, Antolik M, Dreger D, Larsen S, Romanowicz B. "Three-dimensional structure influences on the strong motion wavefield of the 1989 Loma Prieta earthquake" Bull. Seism. Soc. Am. 1999; 99: 1184-1202.

Wald D, Helmberger D, Heaton T. "Rupture history of the 1989 Loma Prieta, California earthquake." Bull. Seism. Soc. Am. 1991; 81: 1540-1572.

Wills C, Petersen M, Bryant W, Reichle M, Saucedo G, Tan S, Taylor G, Treiman J. "A site conditions map for California based on geology and shear wave velocity." Bull. Seism. Soc. Am. 2000; 90(6B): S187-S208.

Zeng Y, Anderson JG, Yu G. "A composite source model for computing syntheticstrong ground motions." Geophys. Res. Lett. 1994; 21: 725-728.

-----------------------

[pic]

Figure 1. Slip velocity function used in the deterministic simulations [see equation (1)].

[pic]

Figure 2. Frequency dependent amplification functions with an input PGA of 20% g, Vref = 620 m/s and various site velocities.

................
................

In order to avoid copyright disputes, this page is only a partial summary.

Google Online Preview   Download