Full-waveform Modeling Research Articles

Full-Waveform Modeling of Complex Media Seismic Waves for Irregular Topography and Its Application in Metal Ore Exploration

Full-Waveform Modeling of Complex Media Seismic Waves for Irregular Topography and Its Application in Metal Ore Exploration

Forward modeling of quake’s infrasound recorded in the stratosphere on board balloon platforms

Acoustic waves generated by seismic waves contain information on the internal structure of planets, and can be sensed by pressure sensors onboard high-altitude balloons. To identify the various contributions (infrasound signal, noise, balloon response, etc.) in such pressure records, a full waveform modeling is implemented and completed by infrasound ray tracing and additional data analysis. Here, we analyze the Stratéole-2 pressure data associated with two earthquakes (Garcia et al. Geophys. Res. Lett. 49(15):e98844, 2022) and compared these to full waveform simulations by SPECFEM2D-DG-LNS software. Even if our simulations do not precisely reproduce the waveform observed in the frequency range [0.05, 0.3] Hz, we show that the waveform presents more sensitivity to quake and internal structure parameters than to atmospheric structure, and that seismic surface wave dispersion is observed in balloon pressure records. The long-duration pressure oscillations observed after the main infrasonic signal cannot be fully reproduced by our one-dimensional input model even when source time function complexity and aftershocks are considered. These features are ascribed mainly to the complex vertical ground movements below the balloon and partly to late secondary infrasound arrivals excited by the interactions of seismic waves with the topography. These results enhance the advantages and limitations of quake-related infrasound observations on board terrestrial and planetary balloon platforms.Graphical

Generative adversarial network-enhanced directional seismic wavefield decomposition and its application in reverse time migration

Reverse time migration (RTM) is an accurate method for imaging complex geologic structures without imposing any dip limitations. However, a large amount of high-amplitude low-frequency noise, which is mainly generated by the crosscorrelation of source and receiver wavefields propagating in the same directions, seriously contaminates the image quality. The causal imaging condition with separated up- and downgoing wavefields is an effective approach to reduce these low-frequency artifacts. Explicit up- and downgoing wavefield decomposition based on the Hilbert transform is computationally expensive due to additional wavefield extrapolation and storage for the imaginary parts. Directionally propagating wavefields have distinctive kinematic patterns such as traveltime and wavefront curvature, which provide us an opportunity to implement the wavefield decomposition using the statistical neural network method. Using extrapolated wavefields as the input and the decomposed up-, down-, left-, and rightgoing wavefields as the labeled data, we train a pair of generative adversarial networks to predict the directional wavefields. The training data sets are generated using seismic full-waveform modeling and explicit wavefield decomposition based on the Hilbert transform. Then, the decomposed directional wavefields are incorporated into a novel imaging condition that depends on the subsurface dip angles to compute the reflectivity perpendicular to the reflectors. Numerical experiments demonstrate that our method can produce accurate directional wavefield decomposition results and high-quality reflectivity images without low-wavenumber artifacts.

Robust quantitative estimation of the seismic attenuation from shallow geotechnical borehole VSP data

Estimating seismic attenuation from shallow geotechnical borehole surveys can be a delicate task. Measured data are often collected within the source near-field domain and the classic inverse-distance correction of the geometric spreading (GS) is inappropriate at this scale. We develop a novel approach based on a 3D full-waveform modeling to substitute the inverse-distance correction. It consists of scaling the picked amplitudes using their counterparts obtained from an elastic simulation carried out under conditions mimicking the data acquisition. The seismic attenuation may be inferred from the corrected amplitudes using either a piecewise regression or a ray-based inversion. Numerical experiments involving P- and S-wave synthetic data indicate that our correction better compensates for the GS effect than the inverse-distance correction. For a synthetic example with 5% noise, the Q-factor values derived from amplitude corrected via the proposed approach have a relative error of approximately 10% compared with 40% for the traditional correction. We investigate the effect of the velocity and density uncertainty upon the calculated correction terms and show that our approach is unbiased and stable. Finally, the robustness of our workflow is assessed on a real case study involving a P-wave data set acquired in a geotechnical borehole.

Full waveform model for axisymmetric black hole mergers

Full waveform model for axisymmetric black hole mergers

Seismoelectric wave conversions at an interface: a quantitative comparison between laboratory data and full-waveform modelling

SUMMARYSeismo-electromagnetic phenomena, electrokinetic in nature, take place whenever a seismic wave propagates in fluid-bearing media, its energy depending mainly on the electrical properties of the fluid and the hydraulic properties of the porous medium. They result from a conversion of mechanical into electromagnetic (EM) energy due to the transient ionic interactions occurring at the pore scale. Two of these phenomena are usually studied: the electric field accompanying seismic waves, and the EM field that travels independently, generated at discontinuities of physicochemical properties in the porous medium. Although the first event is sensitive to physical parameters of the surrounding medium, the second catches information about interfaces in the subsurface, with the resolution of seismic methods, making it very attractive to near surface exploration. In this context, we propose a new experimental setup where both phenomena can be simultaneously studied. At first, we use a porous medium composed of homogeneous water-saturated sand and study the characteristics of the coseismic electric field. Afterwards, a thin layer of Vosges sandstone is inserted into the sand, which allows the study of the EM waves generated at the two closely spaced sand-sandstone interfaces. We record the seismic displacement field at the upper surface of the sand volume using a laser vibrometer, and use stainless steel electrodes buried in the sand to acquire individual electric potentials rather than electric fields, seeking to favour the measurement of the EM interface-generated signals. With the help of direct numerical simulations, we compare experimental measurements and theoretical predictions, based on a well established set of seismoelectric governing equations. In both types of experiments, this comparison shows very good agreements between experimental and numerical waveforms, thus confirming the relevant theory. The electric potential data also show that the EM signals generated at interfaces are clearly recorded at distances of about 10 seismic P wavelengths away from the interface. By contrast, the same events are barely noticeable near the inserted layer when measured using classical electric dipolar arrays.

3D Sensitivity Kernels With Full Attenuation Computed by a Combination of the Strong Stability Preserving Runge‐Kutta Method and the Scattering‐Integral Method

AbstractTo date developments of seismic attenuation models of the Earth have lagged those of velocity models. This is partly due to difficulties in isolating waveform perturbations caused by attenuation and velocity heterogeneities and partly due to different theories (e.g., ray theory, finite frequency theory) used in most previous and current studies to approximate sensitivity kernels and invert for the attenuation structure. We present in this paper a new method for computing the 3D waveform Fréchet kernels that account for full physical‐dispersion and dissipation attenuation. For solving the 3D isotropic anelastic wave equation described by the generalized Maxwell body model, we extend our previously proposed full waveform modeling method in Cartesian coordinates to that in spherical coordinates, which can provide stable numerical solutions even in the presence of strong attenuation. Then, we apply the scattering‐integral method for calculating 3D travel time and amplitude sensitivity kernels with respect to velocity and attenuation structures. We demonstrate the accuracy of our forward method and the effectiveness of the implementation of absorbing and free surface boundary conditions through numerical tests. Moreover, by choosing the Northwestern United States region as a realistic example, we verify the accuracy of the computed 3D sensitivity kernels through comparing the waveform measurements with predictions from the kernels. Finally, we discuss the importance of calculating full anelastic sensitivity kernels including both effects of physical dispersion and dissipation, where we specially explore the effect of scattering due to random velocity and attenuation heterogeneities on waveform measurements.

Hypocenter, Fault Plane, and Rupture Characterization of Australian Earthquakes: Application to the September 2021 Mw 5.9 Woods Point Earthquake

AbstractThe Australian Seismometers in Schools (AuSIS) network operates 50 broadband seismic stations across Australia that are hosted at schools. The instruments augment the Australian National Seismograph Network providing valuable data from urban and regional Australia. The network coverage is quite sparse, but these vital records of rare, moderate Australian earthquakes can improve our understanding of the deformation within the stable continental region of Australia, especially for events with no surface rupture. In this study, we present the feasibility of identifying the fault plane of moderate earthquakes on the Australian continent, using data from the AuSIS network. We examine the fault plane of the September 2021 Mw 5.9 Woods Point earthquake that occurred about 130 km northeast of the Melbourne metropolitan area. We estimate the hypocenter and the centroid moment tensor (CMT) to identify the fault plane from the auxiliary plane in the focal mechanism. We explore a range of 1D models and a 3D Earth model to simulate seismic arrivals and full waveform data. The hypocenter is resolved using P- and S-wave arrivals in a probabilistic framework and the CMT is derived from full waveform modeling through grid search over a set of trial points around the hypocenter. Our solution suggests the mainshock ruptured the depth of 15 ± 4 km, with a strike-slip mechanism striking 348° north on a nearly vertical plane. The high double-couple percentage of this event indicates a simple rupture that propagated from the south (hypocenter) toward the north (centroid) and remained subsurface. This indicates that the causative fault had a deeper structure than the previously known shallow, northwest–southeast-striking faults of the region. The P and T axes deduced from our fault model are notably aligned with the maximum horizontal crustal stress in the region.

Rapid Seismic Waveform Modeling and Inversion With Neural Operators

Seismic waveform modeling is a powerful tool for determining earth structure models and unraveling earthquake rupture processes, but it is usually computationally expensive. We introduce a scheme to vastly accelerate these calculations with a recently developed machine learning paradigm called the neural operator. Once trained, these models can simulate a full wavefield at negligible cost. We use a U-shaped neural operator to learn a general solution operator to the 2D elastic wave equation from an ensemble of numerical simulations performed with random velocity models and source locations. We show that full waveform modeling with neural operators is nearly two orders of magnitude faster than conventional numerical methods, and more importantly, the trained model enables accurate simulation for velocity models, source locations, and mesh discretization distinctly different from the training dataset. The method also enables convenient full-waveform inversion with automatic differentiation.

3-D Full-Waveform Modeling and Analysis of Induced Polarization and Magnetic Viscosity Effect in Time-Domain Electromagnetic Method

3-D Full-Waveform Modeling and Analysis of Induced Polarization and Magnetic Viscosity Effect in Time-Domain Electromagnetic Method

An analysis of core–mantle boundary related seismic waves using full-waveform modelling and adjoint methods

SUMMARYThe core–mantle boundary (CMB) is the most abrupt internal discontinuity in the Earth, marking the solid–fluid boundary between mantle and outer core that strongly affects the dynamics of the Earth’s interior. However, good agreement between models of CMB topographic variations is still lacking. This is probably due to difficulties relating to observations on seismograms and to the lack of good models of lowermost mantle velocity structure. Using spectral-element synthetic seismograms and adjoint methods, we perform traveltime analyses of seismic waves interacting with the CMB. We focus on reflected and refracted P and S waves. We select some of the most important and routinely used seismic phases: ScS, SKS, SKKS, PcP, PKP, PKKP and PcS, given their path through mantle and core and their interaction with the CMB. These seismic waves have been widely deployed by seismologists trying to image CMB topography and lowermost mantle structure. To analyse the reliability of measuring their traveltimes to infer CMB topography, we perform experiments in two ways. First, we compute synthetic seismograms with a dominant period of T ≈ 11s, for computational efficiency, using existing models of CMB topography. We compare traveltime perturbations measured by cross-correlation on the synthetics to those predicted using ray theory. We find deviations from a perfect agreement between ray theoretical predictions of time shifts and those measured on synthetics with and without CMB topography. Second, we calculate Fréchet sensitivity kernels of traveltimes with respect to shear and compressional wave speeds. We also explicitly compute boundary sensitivities with respect to the CMB interface. We observe that the overall sensitivity of the traveltimes is mostly due to volumetric velocity structure and that imprints of CMB on traveltimes are less pronounced. Our study explains the observed difficulties relating to inferring CMB topography using traveltimes and provides a suite of finite frequency sensitivity kernels computed with the adjoint method. The kernels allow us to qualitatively explain the behaviour of measured traveltimes and understand the trade-off between velocity and CMB topography. They can also serve as reference of finite frequency effects on traveltimes of observed seismic phases. From our analyses we conclude that: i) traveltime anomalies measured on Swaves are more in accord with ray theoretical predictions, ii) PcP, PKP, ScS and SKS phases have more pronounced sensitivity to the boundary and iii) separating the greater effects of velocity from those due to the boundary structure is difficult, as they intricately affect the traveltime. We propose that jointly inverting for CMB topography and lowermost mantle velocity structure using full-waveform synthetics and adjoint sensitivity kernels can progress our understanding of deep Earth structure and finite frequency effects on observed waveforms.

Adaptive parameter-related implementation of the free surface in elastic anisotropic full-waveform modeling

The near-surface anisotropy is one of the factors that lead to errors in static corrections, seismic imaging, and inaccurate surface velocity building. Therefore, it is important to study the propagation features of seismic waves in the near-surface with anisotropy. The finite-difference method is commonly used in seismic wave forward modeling, however, the numerical implementation of the free surface is nontrivial in this method, in particular to its implementation in the anisotropic media which is less studied. Here, we propose an adaptive parameter-related free surface implementation. This method is developed from the physical viewpoint with the average medium theory and limit theorem, in which the stress-free condition is implemented by modifications of the original constitutive relation and density at the free-surface points without additional computational requirements. Numerical experiments with the increasing complexity of anisotropy in both 2D and 3D demonstrate that this method can produce precise simulation results with a spatial sampling of 10 grid points per minimum wavelength. In addition, a time-space discontinues grid strategy can be introduced in this implementation to allow the usage of the coarse grid in the deeper part where the oversampling for the surface waves is not required, which increases the flexibility and the computational efficiency of the proposed method.

Predicting infrasound transmission loss using deep learning

SUMMARY Modelling the spatial distribution of infrasound attenuation (or transmission loss, TL) is key to understanding and interpreting microbarometer data and observations. Such predictions enable the reliable assessment of infrasound source characteristics such as ground pressure levels associated with earthquakes, man-made or volcanic explosion properties, and ocean-generated microbarom wavefields. However, the computational cost inherent in full-waveform modelling tools, such as parabolic equation (PE) codes, often prevents the exploration of a large parameter space, that is variations in wind models, source frequency and source location, when deriving reliable estimates of source or atmospheric properties—in particular for real-time and near-real-time applications. Therefore, many studies rely on analytical regression-based heuristic TL equations that neglect complex vertical wind variations and the range-dependent variation in the atmospheric properties. This introduces significant uncertainties in the predicted TL. In the current contribution, we propose a deep learning approach trained on a large set of simulated wavefields generated using PE simulations and realistic atmospheric winds to predict infrasound ground-level amplitudes up to 1000 km from a ground-based source. Realistic range dependent atmospheric winds are constructed by combining ERA5, NRLMSISE-00 and HWM-14 atmospheric models, and small-scale gravity-wave perturbations computed using the Gardner model. Given a set of wind profiles as input, our new modelling framework provides a fast (0.05 s runtime) and reliable (∼5 dB error on average, compared to PE simulations) estimate of the infrasound TL.

Fully scalable solver for frequency-domain visco-elastic wave equations in 3D heterogeneous media: A controllability approach

We develop a controllability strategy for the computation of frequency-domain solutions of the 3D visco-elastic wave equation, in the perspective of seismic imaging applications. We generalize the controllability results for such equations beyond the sound-soft scattering (obstacle) problem. We detail the conjugate gradient implementation and show how an inner elliptic problem needs to be solved to compute the Riesz representative of the gradient at each iteration. We select a spectral-element spatial discretization and a fourth-order Runge-Kutta time discretization. We implement the controllability method in the framework of the SEM46 full waveform modeling and inversion software, to inherit for its excellent scalability which relies on an efficient domain decomposition algorithm. We perform a series of numerical experiments to validate the approach and illustrate its scalability up to more than fifteen hundred cores. In this case, with an elapsed time of less than 50 minutes, we solve a problem on a cubic domain containing up to 160 wavelengths in each direction, involving more than 1.7 billion unknowns.

Suppressing position approximation imaging artifacts in full waveform modeling and inversion of marine seismic data

Full-waveform inversion (FWI) is an advanced seismic imaging method for determining high resolution velocity models of different scales, such as in subseafloor geo-hazard monitoring and risk analysis. Among all two-way wave-equation based modeling methods, a mixed-grid finite-difference modeling (FDM) in frequency domain has the advantage of low numerical dispersion errors for mesh size up to 1/4 wave-length. However, irregular source and receiver positions situated away from the FDM nodes introduce position approximation error (PAE) as large as half mesh size. When a source or receiver is situated in strong heterogeneity or near a free surface, the PAE induces imaging artifacts in high frequencies. We propose an imaging artifact suppression method for multi-scale FDM and FWI in marine seismic survey. The delta function of the source or receiver position is approximated by a discrete Kaiser-windowed sinc function in unbounded medium. The free-surface effects to sources and receivers are compensated with the mirror source concept. The final coefficients for delta function approximation and wavefield interpolation are spatially confined tensor products. The results using ocean bottom cable and vertical cable survey layouts in synthetic models validate the accuracy and feasibility of this new FDM and FWI strategy. The FWI tests shows that this strategy has reduced small scale artifacts at high frequencies for survey layouts of good and poor illuminations. The extension of our method to 3-D visco-acoustic model or reverse time migration is straightforward.

3D Microseismic Monitoring Using Machine Learning

AbstractMicroseismic source localization is important for inferring the dynamic status of the subsurface stress field during hydraulic fracturing. Traditional deterministic methods for 3D microseismic source localization require either ray tracing or full waveform modeling, thus are computationally expensive. We propose a very efficient (e.g., within 1 s) microseismic source localization method based on machine learning. First, three‐dimensional (3D) ray tracing is performed with hypothetical event locations and realistic acquisition geometry to calculate the theoretical travel‐times. The theoretical travel‐time differences and the spatial locations of the stations are treated as the input features of the training data set, and the corresponding source locations are used as the labels. The manually or automatically picked arrival time differences between different stations and a reference station after a microseismic event and the actual station locations are fed into the well‐trained model for a fast and accurate location prediction. The proposed method is efficient enough to be widely applied for the real‐time monitoring of hydraulic fracturing. The machine learning model is analogous to 3D grid search but performs 3D ray tracings before the actual localization, and needs to be retrained when applied to a new study area. We use several synthetic tests and a real data example from the Weiyuan hydraulic fracturing experiment in Sichuan, China, to demonstrate the effectiveness of the proposed method. The application of the proposed method to earthquake localization is also demonstrated to be straightforward.

3-D multiparameter full-waveform inversion for ocean-bottom seismic data using an efficient fluid–solid coupled spectral-element solver

SUMMARYOcean-bottom seismic acquisition systems deployed on the seabed give access to three-component geophone data and hydrophone data. Compared with conventional streamer acquisitions, the separation of sources and receivers makes it possible to increase the maximum offset and azimuth coverage for improving the illumination at depth. Furthermore, the three-component geophones naturally capture elastic wave propagation effects. While this information is mostly overlooked up to now, reconstructing jointly P- and S-wave velocities would significantly improve the subsurface characterization. To achieve a 3-D high-resolution multiparameter reconstruction, we design an efficient 3-D fluid–solid coupled full waveform modelling and inversion engine. In this engine, fluid and solid domains are divided explicitly and handled with the acoustic and elastic wave equations, respectively. The numerical implementation is based on a time-domain spectral-element method (SEM) with a flexible 3-D Cartesian-based hexahedral mesh, which contributes to an accurate coupling of the acoustic and elastic wave equations and high computational efficiency through domain-decomposition based parallelization. We select the best acoustic–elastic coupled formulations among 4 possibilities with criteria based on numerical accuracy and implementation efficiency. Moreover, we propose a specific hybrid approach for the misfit gradient building so as to use a similar modelling solver for both forward and adjoint simulations. Synthetic case studies on a 3-D extended Marmousi-II model and a 3-D deep-water crustal-scale model illustrate how our modelling and inversion engine can efficiently extract information from ocean-bottom seismic data to simultaneously reconstruct both P- and S-wave velocities within a full waveform inversion framework.

Empirical Acoustic Source Model for Chemical Explosions in Air

ABSTRACT Chemical explosions generate pressure disturbances in air that radiate as nonlinear shock waves near the source and transition into acoustic waves with distance. Because low-frequency acoustic waves generally travel large distances without significant loss of energy, they are often used for explosion monitoring and yield estimation. However, quantitative relationships between acoustic energy and explosion yields are required for accurate yield estimation. Here, we develop an empirical acoustic source model for chemical explosions from experimental data. The empirical model returns the acoustic pressure waveform for the detonation of 1 kg of trinitrotoluene, which is conventionally used to represent the explosive release of 4.184 MJ of explosion energy. The full-waveform model can be used to predict acoustic signals for an arbitrary yield of a high-explosive detonation based on the standard scaling law and to estimate acoustic energies in a specific frequency range. We evaluate the accuracy of the acoustic source model independently by estimating the yield of other explosive events that are not included in the model development. Statistical characteristics of the model and their implications for the uncertainty quantification of estimated yields are discussed.

Source Characterization for Two Small Earthquakes in Dartmouth, Nova Scotia, Canada: Pushing the Limit of Single Station

AbstractA pair of small earthquakes (MN 2.4 and 2.6, Earthquakes Canada) hit the city of Dartmouth, Nova Scotia, Canada, in early March 2020. The events were recorded by three seismic stations within 200 km, but only one station (HAL, <10 km) is close enough to offer high-quality broadband signals. In this study, we explore their source parameters using the nearest station through waveform modeling. A nearby quarry blast (MN 2.0) with known Global Positioning System coordinates is adopted as a reference for regional velocity model building and location calibration. We first build a half-space velocity model by estimating the P-S travel-time difference of the blast and determine the near-surface velocity through full-waveform modeling (i.e., comparing a set of synthetic waveforms with the observed blast). The velocity model is then used to evaluate the pair of earthquakes, in which waveform fitting and Rg/S amplitude ratios suggest source depths of ∼0.7 km. The epicenters of these two earthquakes are situated in a recently constructed commercial development. Lastly, single-station template matching finds no similar earthquakes near the hypocenters of the two events in the past decade and only three aftershocks in the following four months. Taking advantage of a ground-truth blast and waveform modeling, our study demonstrates the potential to construct a detailed regional velocity model and determine accurate earthquake source parameters in regions where only a single station is available.

Anomalous azimuthal variations with 360° periodicity of Rayleigh phase velocities observed in Scandinavia

SUMMARY We use the recently deployed ScanArray network of broad-band stations covering most of Norway and Sweden as well as parts of Finland to analyse the propagation of Rayleigh waves in Scandinavia. Applying an array beamforming technique to teleseismic records from ScanArray and permanent stations in the study region, in total 159 stations with a typical station distance of about 70 km, we obtain phase velocities for three subregions, which collectively cover most of Scandinavia (excluding southern Norway). The average phase dispersion curves are similar for all three subregions. They resemble the dispersion previously observed for the South Baltic craton and are about 1 per cent slower than the North Baltic shield phase velocities for periods between 40 and 80 s. However, a remarkable sin(1θ) phase velocity variation with azimuth is observed for periods >35 s with a 5 per cent deviation between the maximum and minimum velocities, more than the overall lateral variation in average velocity. Such a variation, which is incompatible with seismic anisotropy, occurs in northern Scandinavia and southern Norway/Sweden but not in the central study area. The maximum and minimum velocities were measured for backazimuths of 120° and 300°, respectively. These directions are perpendicular to a step in the lithosphere–asthenosphere boundary (LAB) inferred by previous studies in southern Norway/Sweden, suggesting a relation to large lithospheric heterogeneity. In order to test this hypothesis, we carried out 2-D full-waveform modeling of Rayleigh wave propagation in synthetic models which incorporate a steep gradient in the LAB in combination with a pronounced reduction in the shear velocity below the LAB. This setup reproduces the observations qualitatively, and results in higher phase velocities for propagation in the direction of shallowing LAB, and lower ones for propagation in the direction of deepening LAB, probably due to the interference of forward scattered and reflected surface wave energy with the fundamental mode. Therefore, the reduction in lithospheric thickness towards southern Norway in the south, and towards the Atlantic ocean in the north provide a plausible explanation for the observed azimuthal variations.