Real-time 3-D space numerical shake prediction for earthquake early warning

In earthquake early warning systems, real-time shake prediction through wave propagation simulation is a promising approach. Compared with traditional methods, it does not suffer from the inaccurate estimation of source parameters. For computation efficiency, wave direction is assumed to propagate on the 2-D surface of the earth in these methods. In fact, since the seismic wave propagates in the 3-D sphere of the earth, the 2-D space modeling of wave direction results in inaccurate wave estimation. In this paper, we propose a 3-D space numerical shake prediction method, which simulates the wave propagation in 3-D space using radiative transfer theory, and incorporate data assimilation technique to estimate the distribution of wave energy. 2011 Tohoku earthquake is studied as an example to show the validity of the proposed model. 2-D space model and 3-D space model are compared in this article, and the prediction results show that numerical shake prediction based on 3-D space model can estimate the real-time ground motion precisely, and overprediction is alleviated when using 3-D space model.     

Scaling behavior and the effects of heterogeneity on shallow seismic imaging of mineral deposits: A case study from Brunswick No. 6 mining area,Canada

We have studied the scaling behavior of compressional-wave velocity and density logs from an exploration borehole that extends down to about 700 m depth in the Brunswick No. 6 mining area, Bathurst Mining Camp, Canada. Using statistical methods, vertical and horizontal scale lengths of heterogeneity were estimated. Vertical scale length estimates from the velocity, density and calculated acoustic impedance are 14 m, 33 m, and about 20 m, respectively. Although the estimated scale length for the acoustic impedance implies a weak scattering environment, elastic finite difference modeling of seismic wave propagation in 2D heterogeneous media demonstrates that even this weak scattering medium can mask seismic signals from small, but yet economically feasible, massive sulfide deposits. Further analysis of the synthetic seismic data suggests that in the presence of heterogeneity, lenticular-shaped targets may only exhibit incomplete diffraction signals whereby the down-dip tails of these diffractions are mainly visible on the stacked sections. Therefore, identification of orebody generated diffractions is much easier on the unmigrated stacked sections than on migrated stacked sections. The numerical seismic modeling in 2D heterogeneous media indicates that in the presence of large horizontal, but small vertical scale lengths (structural anisotropy), identification of massive sulfide deposits is possible, but their delineation at depth requires detailed velocity modeling and processing algorithms which can handle the anisotropy.     

基于广义正交多项式褶积微分算子的地震波场数值模拟方法

地震波场模拟方法研究对于与波动现象有关的地震学问题的重要性是不言而喻的.就目前现有的各种正演算法来说,精度较高的算法（如有限元法、谱元法、高阶有限差分法等）,其计算速度较慢;计算速度较快的算法（如低阶有限差分法、付氏伪谱法等）计算精度却比较低.为了兼顾地震波场模拟的精度与速度,本文推出了一种快速的、高精度地震波场模拟方法（基于Forsyte广义正交多项式的褶积微分算子法）,该方法是以计算数学中的Forsyte广义正交多项式插值函数为基础,构建一个新的褶积微分算子,并将该算子引入到地震波动方程的一阶速度-应力方程的空间微分运算中去,采用时间交错网格有限差分算子替代普通的差分算子以匹配高精度的褶积微分算子,从而构造一种全新的地震波场数值模拟方法.该方法同时具有广义正交多项式方法的高精度和短算子低阶有限差分算法的高速度.通过对算子长度的调节及算子系数的优化,可同时兼顾波场解的全局信息与局部信息.复杂非均匀介质模型中的波场数值模拟实验证实了该方法的可行性及优越性.     

Numerical Simulation of Fault Zone Guided Waves: Accuracy and 3-D Effects

-- Fault zones are thought to consist of regions with reduced seismic velocity. When sources are located in or close to these low-velocity zones, guided seismic head and trapped waves are generated which may be indicative of the structure of fault zones at depth. Observations above several fault zones suggest that they are common features of near fault radiation, yet their interpretation may be highly ambiguous. Analytical methods have been developed to calculate synthetic seismograms for sources in fault zones as well as at the material discontinuities. These solutions can be used for accurate modeling of wave propagation in plane-parallel layered fault zone structures. However, at present it is not clear how modest deviations from such simplified geometries affect the generation efficiency and observations of trapped wave motion. As more complicated models cannot be solved by analytical means, numerical methods must be employed. In this paper we discuss 3-D finite-difference calculations of waves in modestly irregular fault zone structures. We investigate the accuracy of the numerical solutions for sources at material interfaces and discuss some dominant effects of 3-D structures. We also show that simple mathematical operations on 2-D solutions generated with line sources allow accurate modeling of 3-D wave propagation produced by point sources. The discussed simulations indicate that structural discontinuities of the fault zone (e.g., fault offsets) larger than the fault zone width affect significantly the trapping efficiency, while vertical properly gradients, fault zone narrowing with depth, small-scale structures, and moderate geometrical variations do not. The results also show that sources located with appropriate orientations outside and below a shallow fault zone layer can produce considerable guided wave energy in the overlying fault zone layer.     

Fourier–Hilbert versus Hartley–Hilbert transforms with some geophysical applications

An intimate mathematical relation between Hartley and Hilbert transforms is given here in contrast with the well known Fourier and Hilbert transform relations. It is interesting to note that the Fourier–Hilbert and Hartley–Hilbert transforms while possessing the same magnitude differ in phase by 270°. The inverse Hartley–Hilbert transform returns the original function unlike the Fourier–Hilbert transform which results the negative of the original function. Further, it may be realized that the envelope defined here of the analytic signal in both Fourier–Hilbert and Hartley–Hilbert domains numerically remain the same while differing in polarity. The feasibility of Hartley–Hilbert transform for a straight forward interpretation, total magnetic anomaly due to a thin plate from Tejpur, India and self potential data of the Sulleymonkey anomaly in the Ergani Copper district, Turkey are illustrated in contrast with the Fourier–Hilbert transform. This pair of transforms have innumerable geophysical applications.     

The Numerical Modeling of 3-D Elastic Wave Equation Using a High-Order,Staggered-Grid, Finite Difference Scheme

This article provides the application of the high-order, staggered-grid, finite-difference scheme to model elastic wave propagation in 3-D isotropic media. Here, we use second-order, temporal-and high-order spatial finite-difference formulations with a staggered grid for discretization of the 3-D elastic wave equations of motion. The set of absorbing boundary conditions based on paraxial approximations of 3-D elastic wave equations are applied to the numerical boundaries. The trial resuits for the salt model show that the numerical dispersion is decreased to a minimum extent, the accuracy high and diffracted waves abundant. It also shows that this method can be used for modeling wave propagation in complex media with the lateral variation of velocity.     

Acoustic energy transfer to the upper atmosphere from sinusoidal sources and a role of nonlinear processes

In the framework of solving the problem of acoustic energy transfer from near-surface sources through the upper atmosphere, the propagation of sinusoidal signals of different origin is studied. All calculations are made by means of a model that takes into account the inhomogeneity of the atmosphere, nonlinear effects, absorption, divergence of wave front due to long-range acoustic wave propagation, etc., but does not include the effect of gravity. Infrasonic waves of various periods and their absorption at various heights of the atmosphere are investigated. The calculations show that a sinusoidal signal is destroyed by nonlinear processes during its upward propagation; it transforms into two, initial and final, impulses. The location of the “transformation zone” depends on frequency; its height increases with decreasing frequency. The acoustic waves can heat the upper atmosphere, for example, waves with a period of 3 min generated by thunderstorms can heat the atmosphere by up to ΔT_a=13.08 K/day in the region of 323–431 km. The efficiency of a point artificial emitter is too weak to heat the atmosphere significantly.     

WAVE EQUATION MIGRATION WITH THE ACCURATE SPACE DERIVATIVE METHOD*

A stacked seismic section represents a wave-field recorded at regularly spaced points on the surface. The seismic migration process transforms this recorded data into a reflectivity display. In recent years, Jon F. Claerbout and his co-workers developed migration techniques based on the numerical approximation of the wave equation by finite difference methods. This paper describes an alternative method, termed ASD (for Accurate Space Derivative), and its application to the wave equation migration problem. In this approach to the numerical solution of partial differential equations, partial derivatives are computed by finite Fourier transform methods. This migration method can accommodate media with vertical as well as horizontal velocity variations.     

Parallel PSM/FDM Hybrid Simulation of Ground Motions from the 1999 Chi-Chi,Taiwan, Earthquake

-- A new technique for the parallel computing of 3-D seismic wave propagation simulation is developed by hybridizing the Fourier pseudospectral method (PSM) and the finite-difference method (FDM). This PSM/FDM hybrid offers a good speed-up rate using a large number of processors. To show the feasibility of the hybrid, a numerical 3-D simulation of strong ground motion was conducted for the 1999 Chi-Chi, Taiwan earthquake (M_w 7.6). Comparisons between the simulation results and observed waveforms from a dense strong ground motion network in Taiwan clearly demonstrate that the variation of the subsurface structure and the complex fault slip distribution greatly affect the damage during the Chi-Chi earthquake. The directivity effect of the fault rupture produced large S-wave pulses along the direction of the rupture propagation. Slips in the shallow part of the fault generate significant surface waves in Coastal Plain along the western coast. A large velocity gradient in the upper crust can propagate seismic waves to longer distances with minimum attenuation. The S waves and surface waves were finally amplified further by the site effect of low-velocity sediments in basins, and caused the significant disasters.     

Seismic wavefield modeling based on time-domain symplectic and Fourier finite-difference method

Seismic wavefield modeling is important for improving seismic data processing and interpretation. Calculations of wavefield propagation are sometimes not stable when forward modeling of seismic wave uses large time steps for long times. Based on the Hamiltonian expression of the acoustic wave equation, we propose a structure-preserving method for seismic wavefield modeling by applying the symplectic finite-difference method on time grids and the Fourier finite-difference method on space grids to solve the acoustic wave equation. The proposed method is called the symplectic Fourier finite-difference (symplectic FFD) method, and offers high computational accuracy and improves the computational stability. Using acoustic approximation, we extend the method to anisotropic media. We discuss the calculations in the symplectic FFD method for seismic wavefield modeling of isotropic and anisotropic media, and use the BP salt model and BP TTI model to test the proposed method. The numerical examples suggest that the proposed method can be used in seismic modeling of strongly variable velocities, offering high computational accuracy and low numerical dispersion. The symplectic FFD method overcomes the residual qSV wave of seismic modeling in anisotropic media and maintains the stability of the wavefield propagation for large time steps.     

Infrastructure and evolution of ocean-ridge discontinuities in Iceland

There are two main ocean-ridge discontinuities in Iceland: the Tjörnes Fracture Zone (TFZ) and the South Iceland Seismic Zone (SISZ). The TFZ is a 120-km-long and as much as 70-km-wide WNW-trending zone of high seismicity. It has three main seismic lineaments: the Husavik-Flatey Fault (HFF), the Dalvik lineament, and the Grimsey lineament. The HFF, a dextral strike-slip fault and active as a transform fault for about 9 Ma, has a cumulative transform-parallel displacement of some 60 km. Offshore, the HFF is marked by a transform (fracture-zone) valley, 5–10 km wide and 3–4 km deep. Onshore the Flateyjarskagi Peninsula the HFF is marked by a 3–5-km-wide zone of intense crustal deformation with numerous strike-slip and normal faults, transform-parallel dykes, dense sets of mineral veins, and subzones of completely crushed rocks, that is, fault cores. Where the HFF comes on land on Tjörnes there is a similar, but much thinner, zone of crushed rocks. The seismic lineaments are located a few tens of kilometres south (Dalvik) and north (Grimsey) of, and run subparallel with, the HFF. Both lineaments are composed of sets of NNW-trending sinistral faults arranged en echelon.The SISZ is a 70-km-long and 10–20-km wide zone of almost continuous seismicity located between the overlapping West and East Volcanic Zones. It produces the largest earthquakes in Iceland, some of which exceed M7, during which the N–S width of the zone may be as great as 50–60 km. The SISZ is partly covered with Holocene lava flows where the seismogenic faults occur as dextral NNE-trending and sinistral ENE-trending conjugate arrays with push-ups between their nearby ends. The same fault-segment trends occur in the Pleistocene pile north of the Holocene lava flows.The HFF is neither perpendicular to the nearby ridge segments nor parallel with the spreading vector. As a consequence, the North Volcanic Zone has propagated to the north and the Kolbeinsey Ridge to the south during the past 1 Ma, resulting in the development of the Grimsey and Dalvik lineaments. Similarly, the tip of the East Volcanic Zone has been propagating rapidly to the southwest during the past 3 Ma. The tip has been at its present location for no more than several hundred thousand years, thus making the SISZ less stable than the HFF. If the propagation of the tip of the East Volcanic Zone continues, it will eventually reach the Reykjanes Ridge, whereby either the West or the East Volcanic Zone becomes extinct. Then the SISZ dies out as a major seismic zone.     

Elasticity of antigorite,seismic detection of serpentinites,and anisotropy in subduction zones

Serpentinization of the mantle wedge is an important process that influences the seismic and mechanical properties in subduction zones. Seismic detection of serpentines relies on the knowledge of elastic properties of serpentinites, which thus far has not been possible in the absence of single-crystal elastic properties of antigorite. The elastic constants of antigorite, the dominant serpentine at high-pressure in subduction zones, were measured using Brillouin spectroscopy under ambient conditions. In addition, antigorite lattice preferred orientations (LPO) were determined using an electron back-scattering diffraction (EBSD) technique. Isotropic aggregate velocities are significantly lower than those of peridotites to allow seismic detection of serpentinites from tomography. The isotropic V_P/V_S ratio is 1.76 in the Voigt–Reuss–Hill average, not very different from that of 1.73 in peridotite, but may vary between 1.70 and 1.86 between the Voigt and Reuss bonds. Antigorite and deformed serpentinites have a very high seismic anisotropy and remarkably low velocities along particular directions. V_P varies between 8.9 km s^? 1 and 5.6 km s^? 1 (46% anisotropy), and 8.3 km s^? 1 and 5.8 km s^? 1 (37%), and V_S between 5.1 km s^? 1 and 2.5 km s^? 1 (66%), and 4.7 km s^? 1 and 2.9 km s^? 1 (50%) for the single-crystal and aggregate, respectively. The V_P/V_S ratio and shear wave splitting also vary with orientation between 1.2 and 3.4, and 1.3 and 2.8 for the single-crystal and aggregate, respectively. Thus deformed serpentinites can present seismic velocities similar to peridotites for wave propagation parallel to the foliation or lower than crustal rocks for wave propagation perpendicular to the foliation. These properties can be used to detect serpentinite, quantify the amount of serpentinization, and to discuss relationships between seismic anisotropy and deformation in the mantle wedge. Regions of high V_P/V_S ratios and extremely low velocities in the mantle wedge of subduction zones (down to about 6 and 3 km.s^?1 for V_P and V_S, respectively) are difficult to explain without strong preferred orientation of serpentine. Local variations of anisotropy may result from kilometer-scale folding of serpentinites. Shear wave splittings up to 1–1.5 s can be explained with moderately thick (10–20 km) serpentinite bodies.     

Combined three-dimensional electric and seismic tomography study on the Åknes rockslide in western Norway

We present a combined 3-D geoelectric and seismic tomography study conducted on the large Åknes rockslide in western Norway. Movements on the slope are strongly influenced by water infiltration, such that the hydrogeological regime is considered as a critical factor affecting the slope stability. The aim of our combined geophysical study was to identify and visualize the main shallow tension fractures and to determine their effect on hydraulic processes by comparing the geophysical results with information from borehole logging and tracer tests. To resolve the complex subsurface conditions of the highly fractured rock mass, a three-dimensional set-up was chosen for our seismic survey. To map the water distribution within the rock mass, a pattern of nine intersecting 2-D geoelectric profiles covered the complete unstable slope. Six of them that crossed the seismic survey area were considered as a single data set in a 3-D inversion. For both methods, smoothing-constraint inversion algorithms were used, and the forward calculations and parameterizations were based on unstructured triangular meshes. A pair of parallel shallow low-velocity anomalies (< 1400 m/s) observed in the final seismic tomogram was immediately underlain by two anomalies with resistivities <13 kΩm in the resistivity tomogram. In combination with borehole logging results, the low-velocity and resistivity anomalies could be associated with the drained and water-filled part of the tension fractures, respectively. There were indications from impeller flowmeter measurements and tracer tests that such tension fractures intersected several other water-filled fractures and were responsible for distinct changes of the main groundwater flow paths.     

退化的Fourier偏移算子及其在复杂断块成像中的应用

波动方程宽角抛物逼近得到的通常是非常系数的单程波传播算子，其系数是速度横向变化的函数，因此需要利用有限差分(FD)进行数值实施. 通过对Lippmann Schwinger单程波动积分方程的退化核逼近,本文研究了一类宽角退化算子的偏移成像. 这种退化偏移算子只用快速Fourier变换进行波场延拓，将常规的Fourier分裂步地震偏移方法(SSF)推广适应强速度横向变化介质和大角度传播波场. 退化的Fourier偏移算子通过在两个分裂步项之间作波数域线性插值来实现波场延拓，每延拓一层需要比常规的SSF地震偏移方法多一次快速Fourier变换(FFT). 通过SEG/EAGE盐丘模型和实际地震资料的应用表明，退化Fourier偏移算子能很好地对盐下的陡倾角断层和实际地震剖面上的复杂小断块和大断裂地质构造成像.     

A qualitative study of seismic-ionospheric precursor phenomena monitored by a very close to the epicenter VLF and LF receiver

This work investigates the occurrence of disturbances across a wide range of VLF and LF frequencies received prior to a seismic event (Mw = 4), that took place on May 12th 2012, the epicenter of which was very close (14 km) to the VLF/LF station. The signals analyzed were emitted from five VLF and five LF European transmitters. This seismic event produced precursory ionospheric disturbances, identified as spectral distortion, three days before its occurrence, providing a distinct pattern open to further investigation. Although the basis of the ionosphere interaction with seismic phenomena has been well documented in previous studies, the close proximity of the receiver to the seismic event provides a new perspective to this study. The monitored signals have undergone normalization and then they have been processed by means of the Hilbert-Huang Transform. Diagrams of the signals relevant to the phenomena are presented and the disturbances that are present in the raw data are accentuated through further processing.     

Wave Propagation, Scattering and Imaging Using Dual-domain One-way and One-return Propagators

— Dual-domain one-way propagators implement wave propagation in heterogeneous media in mixed domains (space-wavenumber domains). One-way propagators neglect wave reverberations between heterogeneities but correctly handle the forward multiple-scattering including focusing/defocusing, diffraction, refraction and interference of waves. The algorithm shuttles between space-domain and wavenumber-domain using FFT, and the operations in the two domains are self-adaptive to the complexity of the media. The method makes the best use of the operations in each domain, resulting in efficient and accurate propagators. Due to recent progress, new versions of dual-domain methods overcame some limitations of the classical dual-domain methods (phase-screen or split-step Fourier methods) and can propagate large-angle waves quite accurately in media with strong velocity contrasts. These methods can deliver superior image quality (high resolution/high fidelity) for complex subsurface structures. One-way and one-return (De Wolf approximation) propagators can be also applied to wave-field modeling and simulations for some geophysical problems. In the article, a historical review and theoretical analysis of the Born, Rytov, and De Wolf approximations are given. A review on classical phase-screen or split-step Fourier methods is also given, followed by a summary and analysis of the new dual-domain propagators. The applications of the new propagators to seismic imaging and modeling are reviewed with several examples. For seismic imaging, the advantages and limitations of the traditional Kirchhoff migration and time-space domain finite-difference migration, when applied to 3-D complicated structures, are first analyzed. Then the special features, and applications of the new dual-domain methods are presented. Three versions of GSP (generalized screen propagators), the hybrid pseudo-screen, the wide-angle Padé-screen, and the higher-order generalized screen propagators are discussed. Recent progress also makes it possible to use the dual-domain propagators for modeling elastic reflections for complex structures and long-range propagations of crustal guided waves. Examples of 2-D and 3-D imaging and modeling using GSP methods are given.     

球坐标系下谱元法三维地震波场模拟

谱元法已成为区域性乃至大陆性尺度地震波场模拟的重要工具.对于区域或大陆尺度层析成像而言,地球曲率不可忽略,此时模拟地震波传播采用球坐标系更为合适.本文从球坐标系下弹性波动方程弱形式出发,基于球坐标系变分原理给出了球坐标系下求解三维地震波方程的谱元法.另一方面,计算Fréchet敏感核是进行全波形反演的关键,本文借助伴随原理,推导了全波走时层析成像三维Fréchet敏感核表达式.为了验证球坐标系下谱元法的精度,我们将数值模拟结果与normal mode方法得到的解析解在1-D PREM模型下进行了对比.同时,我们将此方法应用到华北克拉通区域,以期获得地球内部结构精确成像.基于3-D全球径向各向异性地幔模型S362ANI和3-D地壳模型Crust1.0,我们建立了华北克拉通初始3-D背景模型,并将数值模拟结果与实际观测台站记录波形资料进行对比分析,利用互相关方法提取走时残差,最后给出了Fréchet敏感核在3-D空间中的分布,这些工作为下一步开展球坐标系下三维大尺度全波形反演奠定了基础.     

高精度有限差分地震波正演方法

三维地震模拟不仅可以更精确地研究地震波的传播规律,而且是三维地震资料处理和解释的工具。本文采用精细积分法用于求解波动方程.对波动方程在空间采用差分格式,时间域采用积分法求解析的方法。文中详细论述了精细积分法的数值方法,并给出了计算公式。理论分析和数值算例的结果表明了用这种混合方法得到的解与精确解十分吻合,比有限差分法具有更高的精度。文中给出的地震波正演模拟算例说明了该方法适用于复杂地表和复杂构造地质体。     

Earthquake site effect modeling in the Granada basin using a 3-D indirect boundary element method

In this work, we develop the indirect boundary element method (IBEM) to simulate the seismic site response in a realistic, large-scale 3-D sedimentary basin. Most previous applications of boundary element method have used full-space Green’s functions for wave propagation between element points. We use half-space Green’s functions, which include the seismic wavefield interactions at the free surface and require only the boundary elements of the basin interface. In this way, the size of the matrix equation for solution in the IBEM can be reduced to approximately a quarter of that using full-space Green’s functions. The site response modeling of the Granada basin in southern Spain using the IBEM shows that the basin-induced scattering waves were identified as propagating back and forth inside the basin. The scattered waves also generate surface waves that are weakly propagated outside of the basin. The wave propagation inside and outside of the basin shows different patterns. We observe that the scattered wave is locally amplified, and its propagation direction deviates from that of the incident waves propagation direction. Therefore, the computed seismic response in the basin could provide us with good estimates of the seismic motion.     

有限元算法在声波方程数值模拟中的频散分析

针对有限元算法在地震波数值模拟中的数值频散问题,利用集中质量矩阵双线性插值有限元算法,推导了二维声波方程的频散函数.在此基础上采用定量分析方法,对比分析了网格纵横长度比变化时的入射方向、空间采样间隔、地震波频率以及地层速度对数值频散的影响.数值算例和模型正演结果表明:当采用集中质量矩阵双线性插值有限元算法时,为了有效地压制数值频散,在所使用震源子波的峰值频率对应的波长内,采样点数目应不少于20个;减小网格长度的纵横比可以有效地抑制入射角(波传播方向与z轴的夹角)较小的地震波的数值频散;地震波频率越高,传播速度越慢,频散越严重,尤其是当相速度与其所对应的频率比值小于2倍空间采样间隔时,不仅会出现严重的数值频散,还会出现假频现象.