S31E-01 INVITED
Shear wave splitting intensity tomography beneath southwestern Japan and coupling with numerical flow models
The inversion of shear wave splitting measurements for anisotropic structure is not often attempted due to the limitations imposed by sparse data and the difficulty of inverting for laterally varying general anisotropy. However, carefully constructed splitting data sets from dense broadband arrays can potentially illuminate regions of the upper mantle and facilitate a tomographic approach. In addition, the application of simplifying assumptions designed to study anisotropy in specific tectonic settings, such as a subduction zone, can reduce the complexity of the problem. We have developed a method for 2.5-D splitting intensity tomography that incorporates constraints from numerical models of geodynamical processes and applied this technique to study upper mantle anisotropy beneath southwestern Japan. We calculate wave-equation splitting intensity sensitivity kernels for the orientation and strength of anisotropy using the Born approximation. We focus on computing sensitivity kernels in heterogeneous, anisotropic starting models taken from a numerical modeling study of anisotropy development in a subduction zone mantle wedge. The calculation of such kernels in realistic heterogeneity is essential for true multi-scale, wave-equation splitting tomography. We use the results of the inversion to refine our flow models with control parameters that best fit the splitting observations, and further improve the models by using perturbed flow models as new starting models for the inversion. In this way, we identify tomographic models of upper mantle anisotropy beneath southwestern Japan that are consistent both with constraints from shear wave splitting observations and with constraints from geodynamics.
S31E-02
Non-linear 3D Born Shear Wave Tomography in Southeastern Asia
We have developed a 3D radially anisotropic shear velocity model of the upper mantle in southeastern Asia from the inversion of long period seismic multimode waveforms. Our approach is based on normal mode perturbation theory, specifically, on a recent modification of the Born approximation, which we call "N-Born", and which includes a non-linear term that allows the accurate inclusion of accumulated phase shifts which arise when the wavepath traverses a spatially extended region with a smooth velocity anomaly of constant sign. We apply the N-Born approximation in the forward modeling part and calculate linear 3D Born kernels in the inverse part. Our starting model is a 3D radially anisotropic model which we derived from a large dataset of teleseismic multimode long period waveforms in the period range 60 to 400 s, using a finite-frequency 2D approximation (NACT, Li and Romanowicz, 1995). This model covered a larger region of East Asia (longitude 30 to 150 degrees and latitude -10 to 60 degrees), while our N-Born model is restricted to a smaller subregion (longitude 75 to 150 degrees and latitude 0 to 45 degrees) for computational efficiency. In this subregion, our N-Born isotropic and anisotropic models are both parameterized at relatively short wavelengths corresponding to a spherical spline level 6 (~200km). Our N-Born model can fit waveforms as well as the NACT model, with up to ~ 83% variance reduction. While the models agree in general, the N-Born isotropic model shows a stronger fast velocity anomaly beneath the Tibetan plateau in the depth range of 150 km to 250 km, which disappears at greater depth, consistent with other studies. More importantly, the N-Born anisotropic model can recover well the downwelling structure associated with subducted slabs. Beneath the Tibet plateau, radial anisotropy shows VSH>VSV, which is indicative of horizontal rather than vertical flow and may help distinguish between end member models of the tectonics of Tibet.
S31E-03
Arbitrary-resolution global sensitivity kernels
Extracting observables out of any part of a seismogram (e.g. including diffracted phases such as Pdiff) necessitates the knowledge of 3-D time-space wavefields for the Green functions that form the backbone of Fréchet sensitivity kernels. While known for a while, this idea is still computationally intractable in 3-D, facing major simulation and storage issues when high-frequency wavefields are considered at the global scale. We recently developed a new "collapsed-dimension" spectral-element method that solves the 3-D system of elastodynamic equations in a 2-D space, based on exploring symmetry considerations of the seismic-wave radiation patterns. We will present the technical background on the computation of waveform kernels, various examples of time- and frequency-dependent sensitivity kernels and subsequently extracted time-window kernels (e.g. banana- doughnuts). Given the computationally light-weighted 2-D nature, we will explore some crucial parameters such as excitation type, source time functions, frequency, azimuth, discontinuity locations, and phase type, i.e. an a priori view into how, when, and where seismograms carry 3-D Earth signature. A once-and-for-all database of 2-D waveforms for various source depths shall then serve as a complete set of global time-space sensitivity for a given spherically symmetric background model, thereby allowing for tomographic inversions with arbitrary frequencies, observables, and phases.
S31E-04
Global Attenuation Tomography and Implications for Upper-Mantle Thermal Structure
Observation of seismic-wave attenuation provides a direct measure of the Earth's anelasticity. The sensitivity of attenuation to temperature, composition, partial melt, and water content is different from that of seismic velocity, and joint interpretation of elastic and anelastic models may be used to improve constraints on these properties throughout the Earth. Historically, the development of attenuation models has lagged behind velocity models. However, the availability of large seismic datasets and improved techniques to treat these data have recently led to better and higher-resolution attenuation models. We have developed a new 3-D global model of shear attenuation in the upper mantle. This new model, QRFSI12, is derived from > 30,000 fundamental-mode Rayleigh wave amplitude measurements at each period (period range 50-250 s). The amplitudes are inverted simultaneously for the coefficients of the 3-D model as well as frequency-dependent amplitude correction factors for each source and receiver. We have found that focusing by elastic heterogeneity can significantly influence surface-wave amplitudes and that this effect can be modeled at long periods using ray-theoretical approximations. We therefore subtract focusing effects from the data prior to inversion by using phase-velocity maps determined from jointly inverting amplitude and phase-delay datasets. In the shallow mantle, QRFSI12 exhibits a strong correlation with tectonic features, and different tectonic provinces are characterized by distinct attenuative properties. At depths > 250 km, the model is dominated by high attenuation beneath the southeastern Pacific and eastern Africa and low attenuation associated with subduction zones in the western Pacific. Comparison of QRFSI12 with global shear-velocity models shows a strong anti-correlation throughout the upper mantle. At 100-km depth, a clear trend of increasing velocity and decreasing attenuation with increasing age of the seafloor is apparent, and tectonically active continental areas are associated with slower velocities and higher attenuation than stable continental interiors. At depths of 150 and 200 km, oceanic regions exhibit a larger decrease in attenuation per fractional increase in velocity than stable continental regions do, suggesting differences in the mechanisms that influence the seismic properties within these two regions. Comparison with recent laboratory measurements (Faul and Jackson, 2005) of attenuation and velocity for olivine helps to quantify the extent to which temperature alone can explain the observed variability. We find that the mineral-physics predictions agree well with the global seismic models for the oceanic regions between 150- and 250-km depth, but that the cratonic areas cannot be fit.
S31E-05
Can we go From Tomographically Determined Seismic Velocities to Composition? Amplitude Resolution Issues in Local Earthquake Tomography
There have been a number of recent papers (i.e. Lee (2003), James et al. (2004), Hacker and Abers (2004), Schutt and Lesher (2006)) which calculate predicted velocities for xenolith compositions at mantle pressures and temperatures. It is tempting, therefore, to attempt to go the other way ... to use tomographically determined absolute velocities to constrain mantle composition. However, in order to do this, it is vital that one is able to accurately constrain not only the polarity of the determined velocity deviations (i.e. fast vs slow) but also how much faster, how much slower relative to the starting model, if absolute velocities are to be so closely analyzed. While much attention has been given to issues concerning spatial resolution in seismic tomography (i.e. what areas are fast, what areas are slow), little attention has been directed at the issue of amplitude resolution (how fast, how slow). Velocity deviation amplitudes in seismic tomography are heavily influenced by the amount of regularization used and the number of iterations performed. Determining these two parameters is a difficult and little discussed problem. I explore the effect of these two parameters on the amplitudes obtained from the tomographic inversion of the Chile Argentina Geophysical Experiment (CHARGE) dataset, and attempt to determine a reasonable solution space for the low Vp, high Vs, low Vp/Vs anomaly found above the flat slab in central Chile. I then compare this solution space to the range in experimentally determined velocities for peridotite end-members to evaluate our ability to constrain composition using tomographically determined seismic velocities. I find that in general, it will be difficult to constrain the compositions of normal mantle peridotites using tomographically determined velocities, but that in the unusual case of the anomaly above the flat slab, the observed velocity structure still has an anomalously high S wave velocity and low Vp/Vs ratio that is most consistent with enstatite, but inconsistent with the predicted velocities of known mantle xenoliths.
S31E-06
Global P-wave Tomography Using Finite-Frequency Modeling and Finite-Frequency Data
We use our matched-filtering measurements of finite-frequency, teleseismic body-wave traveltimes and amplitudes as input data for global-scale "banana-doughnut" tomography. We present a new mantle model of Vp anomalies obtained from joint inversion of P-wave traveltimes and amplitudes. A Qs model is also obtained since amplitudes are inverted simultaneously for focusing and attenuation. In accordance with their definition in the Born approximation, traveltimes and amplitudes were measured by cross-correlating observed and predicted seismograms. Observables are consistent across frequency bands since for any source-receiver path, finite-frequency measurements are derived by filtering a single broadband seismogram into distinct sub-bands. Our bandpass filters have central frequencies between 0.03 and 1~Hz. In this range, predicted waveforms depend heavily on accurate estimates of source time functions, a challenge explicitly addressed by our method. Frechet kernels are computed using the exact frequency responses of the bandpass filters used. Our data set comprises teleseismic P-wave data from 1999-2007, available from IRIS DMC. We have processed all events for which mb ≥ 6.0, and most events for which 5.7 ≤ mb < 6.0, for a total of >1,500 earthquakes. We have thus far concentrated on the inversion of teleseismic P-wave data observed in North America. The obtained models predict traveltime dispersion due to diffraction effects on the order of 0.1-0.5 sec, consistent with what we observe in the data. Amplitude dispersion is on the order of 20%. Our models explain roughly half of the variance associated with the amplitude data. For North America, we observe a close correspondence of upper mantle Vp patterns with known surface tectonics. Thanks to the excellent USArray data, features on the order of 200~km can be resolved in large parts of the upper mantle under the Western U.S. Deeper down and further east, the presumed fragments of the subducting Farallon plate are clearly visible. We expect to include the full data set with global source-receiver coverage and to present preliminary results for a global inversion in addition to the North American tomography.
S31E-07
Multi-scale issues and two scale homogenization solutions for the direct and inverse problems in seismology for layered earth model and beyond.
In many cases, in the seismic wave propagation modeling context, scales much smaller than the minimum wavelength are present in the earth model we wish to propagate in. For many numerical methods these small scales are a challenge leading to high numerical cost. The purpose of the work presented here is to understand and to build the effective medium and equations allowing to average the small scales of the original medium without losing the accuracy of the wavefield computation. We show that high order two scale homogenization provides a promising solution to this king of problem. The presentation will be mostly limited to the layered model case. In that case, it appears that the order 0 homogenization gives the result that was obtained by Backus in 1962 which implies that order 0 homogenized model is transversely isotropic even though the original model is isotropic. It appears the order 0 is not enough to obtain surface wave with correct group and phase velocities and that higher order homogenization terms up to 2 are often required, which implies to modify the wave equation and the boundary conditions. We show how to extend the theory from the periodic case to the non-periodic case. Examples in periodic and non-periodic medium are given. Applications to numerical method like the spectral element method as well as preliminary results of homogenization in media varying rapidly in all direction will be presented. The implication of this work on tomography parameterization will be discussed.
S31E-08 INVITED
Full waveform seismic inversion: Prospects for scaling to petaflops computers
The U.S., Japanese, and European governments are all pursuing petaflops computing programs, and the first systems capable of a peak petaflops performance are expected to appear in 2008. Problems in the geosciences have been among the important drivers for the development of such systems. One such problem is the seismic inverse problem of determining the distribution of earth properties from surface observations of earthquake-induced high-frequency ground motion in large heterogeneous elastic regions. Scalability of inverse methods on supercomputers requires both algorithmic scalability (scaling to large problems sizes) and parallel scalability (scaling to large numbers of processors). Here we focus on the deterministic inverse problem, and we address both parallel efficiency and algorithmic efficiency of a class of optimization methods for solution of full elastic waveform seismic inverse problems as resolution limits are pushed. We discuss algorithmic choices designed to assure scalability for various components of the inverse method, including inexact Newton for nonlinear iterations, conjugate gradients for linear iterations, primal-dual active sets for treatment of parameter bounds, and primal-dual treatment of total variation regularization. We conclude with an assessment of the prospects of scalability of high-resolution seismic inversion to upcoming petascale computing systems.