We provide an introduction to the use of the spectral-element method (SEM) in seismology. Following a brief review of the basic equations that govern seismic wave propagation, we discuss in some detail how these equations may be solved numerically based upon the SEM to address the forward problem in seismology. Examples of synthetic seismograms calculated based upon the SEM are compared to data recorded by the Global Seismographic Network. Finally, we discuss the challenge of using the remaining differences between the data and the synthetic seismograms to constrain better Earth models and source descriptions. This leads naturally to adjoint methods, which provide a practical approach to this formidable computational challenge and enables seismologists to tackle the inverse problem. 

We review the methods of simulating elastic wave propagation in a borehole. We considered two different approaches: a quasi-analytic approach using the Discrete Wavenumber Summation Method, and the purely numerical Finite Difference Method. We consider the special geometry of the borehole and discuss the problem in cylindrical coordinates. We point out some numerical difficulties that are particularly unique to this problem in cylindrical coordinates.

We present a semi-analytic method based on the propagation matrix formulation of indirect boundary element method to compute response of elastic (and acoustic) waves in multi-layered media with irregular interfaces. The method works recursively starting from the top-most free surface at which a stress-free boundary condition is applied, and the displacement-stress boundary conditions are then subsequently applied at each interface. The basic idea behind this method is the matrix formulation of the propagation matrix (PM) or more recently the reflectivity method as wide used in the geophysics community for the computation of synthetic seismograms in stratified media. The reflected and transmitted wave-fields between arbitrary shapes of layers can be computed using the indirect boundary element method (BEM, sometimes called IBEM). Like any standard BEM, the primary task of the BEM-based propagation matrix method (thereafter called PM-BEM) is the evaluation of element boundary integral of the Green's function, for which there are standard method that can be adapted. In addition, effective absorbing boundary conditions as used in the finite difference numerical method is adapted in our implementation to suppress the spurious arrivals from the artificial boundaries due to limited model space. To our knowledge, such implementation has not appeared in the literature. We present several examples in this paper to demonstrate the effectiveness of this proposed PM-BEM for modelling elastic waves in media with complex structure.

This study presents stochastic methods to simulate wave envelopes in layered random media. High-frequency Seismograms of small earthquakes are so complex due to lithospheric inhomogeneity that seismologists often analyze wave envelopes rather than wave traces to quantify the subsurface inhomogeneity. Since the statistical properties of the inhomogeneity vary regionally, it is important to develop and examine direct envelope simulation methods for non-uniform random media. As a simple example,  this study supposes plane wave propagation through two-layer random media in 2-D composed of weak and strong inhomogeneity zones. The characteristic spatial-scale of the inhomogeneity is supposed to be larger than the wavelength, where small-angle scattering around the forward direction dominates large-angle scattering. Two envelope simulation methods based on the small-angle scattering approximation are examined. One method is to solve a differential equation for the two-frequency mutual coherence function with the Markov approximation. The other is to solve the stochastic ray bending process by using the Monte Carlo method based on the Markov approximation for the mutual coherence function. The resultant wave envelopes of the two methods showed excellent coincidence both for uniform and for two-layer random media. Furthermore, we confirmed the validity of the two methods comparing with the envelopes made from the finite difference simulations of waves. The two direct envelope simulation methods presented in this study can be a mathematical base for the study of high-frequency wave propagation through randomly inhomogeneous lithosphere in seismology.

Constraining numerical geodynamo models with surface geomagnetic observations is very important in many respects: it directly helps to improve numerical geodynamo models, and expands their geophysical applications beyond geomagnetism. A successful approach to integrate observations with numerical models is data assimilation, in which Bayesian algorithms are used to combine observational data with model outputs, so that the modified solutions can then be used as initial conditions for forecasts of future physical states. In this paper, we present the first geomagnetic data assimilation framework, which comprises the MoSST core dynamics model, a newly developed data assimilation component (based on ensemble covariance estimation and optimal interpolation), and geomagnetic field models based on paleo, archeo, historical and modern geomagnetic data. The overall architecture, mathematical formulation, numerical algorithms and computational techniques of the framework are discussed. Initial results with 100-year geomagnetic data assimilation and with synthetic data assimilation are presented to demonstrate the operation of the system.

We present theoretical and experimental studies on the effects of formation properties on seismoelectric conversions in fluid-filled boreholes. First, we derive the theoretical formulations for seismoelectric responses for an acoustic source in a borehole. Then, we compute the electric fields in boreholes penetrating formations with different permeability and porosity, and then we analyze the sensitivity of the converted electric fields to formation permeability and porosity. We also describe the laboratory results of the seismoelectric and seismomagnetic fields induced by an acoustic source in borehole models to confirm our theoretical and numerical developments qualitatively. We use a piezoelectric transducer to generate acoustic waves and a point electrode to receive the localized seismoelectric fields in layered boreholes and the electric component of electromagnetic waves in a fractured borehole model. Numerical results show that the magnitude ratio of the converted electric wave to the acoustic pressure increases with the porosity and permeability increases in both fast and slow formations. Furthermore, the converted electric signal is sensitive to the formation permeability for the same source frequency and formation porosity. Our experiments validate our theoretical results qualitatively. An acoustic wave at a fracture intersecting a borehole induces a radiating electromagnetic wave.

We apply the CIP (Cubic Interpolated Profile) scheme to the numerical simulation of the acoustic wave propagation based on characteristic equations.The CIP scheme is based on a concept that both the wavefield and its spatial derivative propagate along the same characteristic curves derived from a hyperbolic differential equation. We describe the derivation of the characteristic equations for the acoustic waves from the basic equations by means of the directional splitting and the diagonalization of the coefficient matrix, and establish geophysical boundary conditions. Since the CIP scheme calculates both the wavefield and its spatial derivatives, it is easy to realize the boundary conditions theoretically. We also show some numerical simulation examples and the CIP can simulate acoustic wave propagation with high stability and less numerical dispersion. The method of characteristics with the CIP scheme is a very powerful technique to deal with the wave propagation in complex geophysical problems. 

We study the 2-parameter acoustic Born series for an actual medium with constant velocity and a density distribution. Using a homogeneous background we define a perturbation, the difference between actual and reference medium (we use background and reference as synonyms), which exhibits an anisotropic behavior due to the density distribution. For an actual medium with a constant velocity, the reference velocity can be selected so that the waves in the actual medium travel with the same speed as the waves in the background medium. Scattering theory decomposes the actual wave field into an infinite series where each term contains the perturbation and the propagators in the background medium. Hence, in this formalism, all propagations occur in the background medium and the actual medium is included only through the perturbations which scatter the propagating waves. The density-only perturbation has an isotropic and an anisotropic component. The anisotropic component is dependent on the incident direction of the propagating waves and behaves as a purposeful perturbation in the sense that it annihilates the part of the Born series that acts to correct the time to build the actual wave field, an unnecessary activity when the reference velocity is equal to the one in the actual medium. This means that the forward series is not attempting to correct for an issue that does not exist. We define the purposeful perturbation concept as the intrinsic knowledge of precisely what a given term is designed to accomplish. This is a remarkable behavior for a formalism that predicts the scattered wave field with an infinite series. At each order of approximation the output of the series is consistent with the fact that the time is correct because the velocity is always constant. In the density-only perturbation, the forward series only seeks to predict the correct amplitudes. Finally, we extend the analysis to a wave propagating in a medium where both density and velocity change. By selecting a convenient set of parameters, we find a conceptual framework for the multiparameter Born series. This framework provides an insightful analysis that can be mapped and applied to the concepts and algorithms of the inverse scattering series.

We develop a new computational method for modeling and inverting frequency domain airborne electromagnetic (EM) data. Our method is based on the contraction integral equation method for forward EM modeling and on inversion using the localized quasi-linear (LQL) approximation followed by the rigorous inversion, if necessary. The LQL inversion serves to provide a fast image of the target. These results are checked by a rigorous update of the domain electric field, allowing a more accurate calculation of the predicted data. If the accuracy is poorer than desired, rigorous inversion follows, using the resulting conductivity distribution and electric field from LQL as a starting model. The rigorous inversion iteratively solves the field and domain equations, converting the non-linear inversion into a series of linear inversions. We test this method on synthetic and field data. The results of the inversion are very encouraging with respect to both the speed and the accuracy of the algorithm, showing this is a useful tool for airborne EM interpretation.

We present the application of the theory of Padé approximants to extending the perturbative solutions of acoustic wave equation for a three dimensional vertically varying medium with one interface. These type of solutions have limited convergence properties depending on either the degree of contrast between the actual and the reference medium or the angle of incidence of a plane wave component. We show that the sequence of Padé approximants to the partial sums in the forward scattering series for the 3D wave equation is convergent for any contrast and any incidence angle. This allows the construction of any reflected waves including phase-shifted post-critical plane waves and, for a point-source problem, refracted events or headwaves, and it also provides interesting interpretations of these solutions in the scattering theory formalism. 

The WKBJ solution for the one-way wave equations in media with smoothly varying velocity variation with depth, c(z), is reformulated from the principle of energy flux conservation for acoustic media. The formulation is then extended to general heterogeneous media with local angle domain methods by introducing the concepts of Transparent Boundary Condition (TBC) and Transparent Propagator (TP). The influence of the WKBJ correction on image amplitudes in seismic imaging, such as depth migration in exploration seismology, is investigated in both smoothly varying c(z) and general heterogeneous media. We also compare the effect of the propagator amplitude compensation with the effect of the acquisition aperture correction on the image amplitude. Numerical results in a smoothly varying c(z) medium demonstrate that the WKBJ correction significantly improves the one-way wave propagator amplitudes, which, after compensation, agree very well with those from the full wave equation method. Images for a point scatterer in a smoothly varying c(z) medium show that the WKBJ correction has some improvement on the image amplitude, though it is not very significant. The results in a general heterogeneous medium (2D SEG/EAGE salt model) show similar phenomena. When the acquisition aperture correction is applied, the image improves significantly in both the smoothly varying c(z) medium and the 2D SEG/EAGE salt model. The comparisons indicate that although the WKBJ compensation for propagator amplitude may be important for forward modeling (especially for wide-angle waves), its effect on the image amplitude in seismic imaging is much less noticeable compared with the acquisition aperture correction for migration with limited acquisition aperture in general heterogeneous media.

This article develops a new, granular approach to tsunami runup and inundation. The small grains employed here are not fluid, but bits, or balls, of tsunami energy. By careful formulation of the ball accelerations, both wave-like and flood-like behaviors are accommodated so tsunami waves can be run seamlessly from deep water, through wave breaking, to the final surge onto shore and back again. In deep water, tsunami balls track according long wave ray theory. On land, tsunami balls behave like a water landslide. In shallow water, the balls embody both deep water and on land elements. In modeling several 2-D and 3-D cases, we find that wave breaking generally causes relative runup to increase with beach slope and wave period and decrease with input wave amplitude. Because of their highly non-linear nature, runup and inundation are best considered to be random processes rather than deterministic ones. Models and observations hint that for uniform input waves, normalized runup statistics everywhere follow a single skewed distribution with a spread between 1/2 and 2 times its mean.

Multiple scattering of elastic waves in realistic media makes that average field intensities or energy densities follow diffusive processes. In such regime the successive P to S energy conversions by distributed random inhomogeneities give rise to equipartition which means that in the phase space the available elastic energy is distributed in average with equal amounts among the possible states of P and S waves. In such diffusive regime the P to S energy ratio equilibrates in an universal way independent of the particular details of the scattering. It has been demonstrated that averaging the cross correlations at any two points of an elastic medium subjected to diffuse elastic wavefields leads to the emergence of the Green function, which is the wave field that would be observed at one position if an impulsive load is applied at the other. In this work we study the problem of the retrieval of the 2D tensor elastodynamic Green function in an infinite elastic space containing a circular cylinder inclusion. We illuminate isotropically the elastic space with plane waves. We assume the spectra for both P and S waves uniform but such that the energy ratio E_S/E_P =(α/β)², which is the one predicted by equipartition theory in two-dimensions. We then show that the Fourier transform of azimuthal average of the cross-correlation of motion between two points within an elastic medium is proportional to the imaginary part of the exact Green tensor function between these points. The numerical results presented here point out the possibility of detection and imaging of diffractors and resonant diffractors by cross correlation even in presence of attenuation exists.