High-order Finite-difference Time-domain Research Articles

SMIwiz: An integrated toolbox for multidimensional seismic modelling and imaging

This paper contributes an open source software - SMIwiz, which integrates seismic modelling, reverse time migration (RTM), and full waveform inversion (FWI) into a unified computer implementation. SMIwiz has the machinery to do both 2D and 3D simulation in a consistent manner. The package features a number of computational recipes for efficient calculation of imaging condition and inversion gradient: a dynamic evolving computing box to limit the simulation cube and a well-designed wavefield reconstruction strategy to reduce the memory consumption when dealing with 3D problems. The modelling in SMIwiz runs independently: each shot corresponds to one processor in a bijective manner to maximize the scalability. A batchwise job scheduling strategy is designed to handle large 3D imaging tasks on computer with limited number of cores. The viability of SMIwiz is demonstrated by a number of applications on benchmark models. Program summaryProgram Title: SMIwizCPC Library link to program files:https://doi.org/10.17632/tygszns27k.1Developer's repository link:https://github.com/yangpl/SMIwizLicensing provisions: GNU General Public License v3.0Programming language: C, Shell, FortranExternal dependencies: MPI [1], FFTW [2]Nature of problem: Seismic modelling and imaging (FWI and RTM)Solution method: High-order finite-difference time-domain (FDTD) for modelling on staggered grid; Quasi-Newton LBFGS algorithm for nonlinear optimization; line search to estimate step length based on Wolfe condition References[1]https://www.mpich.org/[2]http://fftw.org/

LibEMM: A fictious wave domain 3D CSEM modelling library bridging sequential and parallel GPU implementation

This paper delivers a software - libEMM - for 3D controlled-source electromagnetics (CSEM) modelling in fictitious wave domain, based on the newly developed high-order finite-difference time-domain (FDTD) method on non-uniform grid. The numerical simulation can be carried out over a number of parallel processors using MPI-based high performance computing architecture. The FDTD kernel coded in C has been parallelized with OpenMP for speedup using local shared memory. In addition, the software features a GPU implementation of the same algorithm based on CUDA programming language, which can be cross-validated and compared in terms of efficiency. A perspective of libEMM on the horizon is its application to 3D CSEM inversion in land and marine environment. Program summaryProgram Title: libEMMCPC Library link to program files:https://doi.org/10.17632/p769t7c5bk.1Developer's repository link:https://github.com/yangpl/libEMMLicensing provisions: GNU General Public License v3.0Programming language: C, CUDA, Fortran, ShellExternal dependencies: MPI [1], FFTW [2], CUDA [3]Nature of problem: Controlled-source electromagnetics (CSEM)Solution method: High-order finite-difference time-domain (FDTD) on non-uniform grid by fictious wave domain transformation References[1]https://www.mpich.org/[2]http://fftw.org/[3]https://developer.nvidia.com/cuda-toolkit

Controlled-source electromagnetic modeling using a high-order finite-difference time-domain method on a nonuniform grid

Simulation of 3D low-frequency electromagnetic (EM) fields propagating in the earth is computationally expensive. We have developed a fictitious wave-domain high-order finite-difference time-domain (FDTD) modeling method on nonuniform grids to compute frequency-domain 3D controlled-source EM data. The method overcomes the inconsistency issue widely present in the conventional second-order staggered-grid finite-difference scheme over a nonuniform grid, achieving high accuracy with an arbitrarily high-order scheme. The finite-difference coefficients adaptive to the node spacings can be accurately computed by inverting a Vandermonde matrix system using an efficient algorithm. A generic stability condition applicable to nonuniform grids is established, revealing the dependence of the time step and these finite-difference coefficients. A recursion scheme using fixed-point iterations is designed to determine the stretching factor to generate the optimal nonuniform grid. The grid stretching in our method reduces the number of grid points required in the discretization, making it more efficient than the standard high-order FDTD with a densely sampled uniform grid. Instead of stretching in vertical and horizontal directions, better accuracy of our method is observed when the grid is stretched along the depth without horizontal stretching. The efficiency and accuracy of our method are demonstrated by numerical examples.

Numerical simulations of sonic boom propagation over urban areas

Acceptability of supersonic transportation by population requires an accurate prediction of ground noise levels generated by sonic boom. This study aims at predicting sonic boom propagation over urban areas. For this purpose, numerical simulations are performed; the full 2D Euler equations are solved using high-order finite-difference time-domain techniques. First, the case of an isolated building is considered. From a geometrical analysis, two characteristic zones are highlighted: an illumated region in front of the building and a shadow zone at its rear. The sonic boom waveforms at the ground are composed of several arrivals, related to reflection at the building facades and diffraction at the building corners. The evolution of the noise levels is then shown to follow closely the geometrical analysis, with an amplification in the illuminated region and a large reduction in the shadow zone. Second, the case of two identical buildings is investigated. The acoustic field inside the street canyon is examined. In particular, the boom waveforms exhibit low-frequency oscillations, in addition to the geometrical arrivals. They are related to resonant modes of the canyon. Finally, an urban geometry representative of European city centres is considered. The variability of the boom waveforms and the noise levels is shown.

Crosstalk Analysis of Printed Circuit Board Traces With Right-Angled Bent Corners via Time Domain Hybrid Method

Discontinuous impedance exists when bent corners are introduced into the traces of printed circuit board (PCB). By combining transmission line (TL) equations, higher order finite-difference time-domain (FDTD) method, state variable method, and Norton's theorem, an efficient time domain hybrid method is proposed to realize the fast crosstalk simulation of the traces with right-angled bent corners (RABCs) in this work. First, the whole structure of the traces is decomposed into a cascade network, which consists of multiple equal length parallel microstrip lines (EL-PMSLs) and the RABC structures. Then, the crosstalk responses on these EL-PMSLs can be quickly solved by applying TL equations combined with the higher order FDTD method. With the application of Norton's theorem, the equivalent circuit model of whole PCB traces is constructed. Finally, the equivalent circuits of the RABC structures are extracted, and the port voltages of these RABCs are calculated via state variable method and fed back to these EL-PMSLs as boundaries, thus the round-trip transmission of interference signals on the PCB traces is realized. Numerous crosstalk simulations of PCB traces with typical RABC structures, as well as extracted from a practical PCB, are presented and compared with the results of commercial software CST to verify the correctness and efficiency of the proposed method.

A Novel High-Order Symplectic Compact FDTD Schemes for Optical Waveguide Simulation

As a 2-D full-wave numerical algorithm in the time domain, the compact Finite Difference Time Domain (FDTD) is an efficient algorithm for eigenvalue analysis of optical waveguide system. However, the numerical dispersion accuracy and stability of fast algorithm need to be improved while simulating at high frequency. A novel high-order symplectic compact FDTD scheme is developed and validated for optical waveguide modal analysis. The stability condition and the numerical dispersion of schemes with fourth-order accuracy in temporal and spatial using the symplectic integrator and compact scheme are analyzed. By comparisons with other time-domain schemes, their stable and accurate performance is qualitatively verified. The proposed high-order SC-FDTD method can be used for efficiently simulating electrically large and longitudinally invariant optical devices since the reduction of simulation dimensionality and the novel high-order symplectic algorithm can greatly reduce the memory cost and the numerical dispersive errors.

A high-order finite-difference scheme for time-domain modeling of time-varying seismoelectric waves

Simulation of the seismoelectric effect serves as a useful tool to capture the observed seismoelectric conversion phenomenon in porous media, thus offering promising potential in underground exploration activities to detect pore fluids such as water, oil, and gas. The static electromagnetic (EM) approximation is among the most widely used methods for numerical simulation of the seismoelectric responses. However, the static approximation ignores the accompanying electric field generated by the shear wave, resulting in considerable errors when compared to analytical results, particularly under high-salinity conditions. To mitigate this problem, we have adopted a spatial high-order finite-difference time-domain method based on Maxwell’s full equations of time-varying EM fields to simulate the seismoelectric response in 2D mode. To improve the computational efficiency influenced by the velocity differences between seismic and EM waves, different time steps are set according to the stability conditions, and the seismic feedback values of EM time nodes are obtained by linear approximation within the seismic unit time step. To improve the simulation accuracy of the seismoelectric response with the time-varying EM calculation method, finite-difference coefficients are obtained by solving the spatial high-order difference approximation based on the Taylor expansion. Our method yields consistent simulation results compared to those obtained from the analytical method under different salinity conditions, thus indicating its validity for simulating seismoelectric responses in porous media. We further apply our method to layered and anomalous body models and extend our algorithm to three dimensions. Results indicate that the time-varying EM calculation method can effectively capture the reflection and transmission phenomena of the seismic and EM wavefields at the interfaces of contrasting media. This may allow for the identification of abnormal locations, thus highlighting the capability of seismoelectric response simulation to detect subsurface properties.

Characterization of topographic effects on sonic boom reflection by resolution of the Euler equations.

The influence of topography on sonic boom propagation is investigated. The full two-dimensional Euler equations in curvilinear coordinates are solved using high-order finite-difference time-domain techniques. Simple ground profiles, corresponding to a terrain depression, a hill, and a sinusoidal terrain, are examined for two sonic boom waves: a classical N-wave and a low-boom. Ground reflection of the sonic boom is affected by elevation variations: a concave ground profile induces compression, which tends to increase the peak pressure in particular, while the opposite is true for convex elevation variations, which lead to expansion and a reduction in peak pressure. The reflected boom is then strongly altered. Furthermore, a sufficiently concave topography can cause focal zones, which generate extra contributions at ground level in the form of U-waves in addition to the reflected wave. This mechanism has the largest effect on waveforms at ground level. The variations of standard metrics are of a few dBs compared to a flat ground for both sonic boom waves, and they are notably greater for the terrain depression than for the hill. Finally, in the case of a sinusoidal terrain, the pressure waveforms are composed of multiple arrivals due to successive focal zones.

Sonic boom reflection over irregular terrain

An overview of recent studies on sonic boom reflection over irregular terrain is presented. Numerical simulations based on the two-dimensional Euler equations are performed using high-order finite-difference time-domain techniques, allowing for an accurate prediction of diffraction. Studies are done for two sonic boom waves: a conventional N-wave and a low-boom wave. First, ground elevation effects on sonic boom are examined. An academic ground profile, corresponding to a terrain depression, is considered to highlight mechanisms at play. In particular, caustics due to ground elevation variations induce another contribution at the ground level, in the form of a U-wave. Results for real ground profiles are then presented and cumulative effects, due to the succession of terrain irregularities, are shown. Second, propagation of sonic boom over an urban profile is investigated. Possible resonance effects inside street canyons are discussed. For all cases, statistics of noise metrics are examined. [This work was performed within the framework of the LabexCeLyA of the University of Lyon, within the program ‘‘Investissements d’Avenir” (ANR-10-LABX-0060/ANR-16-IDEX-0005) operated by the French National Research Agency and has received funding from the European Union's Horizon 2020 research and innovation programme (Grant Agreement No. 769896 (RUMBLE).]

Quantifying Sub-gridding Errors in Standard and Hybrid Higher Order 2D FDTD Simulations

Sub-gridding errors for a 2D Finite-Difference Time-Domain (FDTD) simulation are compared for both the standard FDTD and Hybrid higher order FDTD cases. Subgridding contrast ratios of 1:3, 1:9, 1:15, and 1:27 are considered and analyzed. A correlation is seen between the increase of contrast ratio with the increase of sub-gridding errors for both standard and hybrid cases. However, a trend of errors reduction when using hybrid formulations over standard formulations is apparent for each contrast ratio.

Nonuniform and Higher-order FDTD Methods for the Schrödinger Equation

Two Finite-Difference Time-Domain (FDTD) methods are developed for solving the Schrödinger equation on nonuniform tensor-product grids. The first is an extension of the standard second-order accurate spatial differencing scheme on uniform grids to nonuniform grids, whereas the second utilizes a higher-order accurate spatial scheme using an extended stencil. Based on discrete-time stability theory, an upper bound is derived for the time step of both proposed schemes. It is shown that the time step derived in this way can be larger compared to the known stability criterion. Furthermore, the numerical dispersion error is investigated as a function of the time step, the spatial step and the propagation direction. Numerical experiments are compared with analytical solutions and demonstrate the increased accuracy of the higher-order scheme as well as the advantageous properties of nonuniform gridding.

Fictitious wave domain modelling and analysis of marine CSEM data

SUMMARY Modelling marine controlled-source electromagnetic (CSEM) responses in the fictitious time domain is a novel approach, which facilitates the full exploration of EM diffusive properties in the fictitious wave domain (FWD). Concepts, such as reflections, refractions, diffractions and transmissions, which are used for the analysis of elastic wave propagation can thus be adopted in FWD for interpreting CSEM data. In this paper, we use a high-order finite difference time domain (FDTD) algorithm for modelling marine CSEM responses in both the fictitious time domain and the diffusive frequency domain. A complex frequency shifted perfectly matched layer (CFS–PML) boundary condition is adopted to the FDTD modelling. We demonstrate the performance of the CFS–PML boundary condition and validate the high-order FDTD code in the FWD with the half-space sea water model and in the frequency domain with the 1-D canonical reservoir model. We investigate and analyse the propagation characteristics of electromagnetic fields in the FWD, where we apply wave propagation concepts to interpret marine CSEM data. Similarities between wave and field propagations relevant for marine CSEM data are demonstrated through several 1-D to 3-D numerical examples.

A 3-D High-Order Reverse-Time Migration Method for High-Resolution Subsurface Imaging With a Multistation Ultra-Wideband Radar System

A three-dimensional high-order reverse-time migration (3-D HO-RTM) method is proposed to perform subsurface electromagnetic imaging with an ultra-wideband radar (UWBR) system consisting of a multi-input and multi-output antenna array. By using a UWBR system to collect temporal scattering signals, subsurface targets can be detected, and the image of targets can be obtained by imaging methods such as the back-propagation method, frequency-wavenumber migration technique, time-reversal mirror, and reverse-time migration method. The proposed HO-RTM method is based on the high-order finite-difference time-domain (HO-FDTD) method to significantly reduce the computational cost in the conventional RTM method. The measured data from an experimental lunar exploration system Chang'E-5 have been collected on a 7 m × 2.5 m × 2.5 m laboratory model with volcanic ash and validated by the 3-D HO-RTM method. Results show that all buried objects can be effectively identified by the HO-RTM, and its computer memory and CPU time are only 3.87% and 0.7128% of the conventional RTM method, respectively.

Numerical modeling for karst cavity sonar detection beneath bored cast in situ pile using 3D staggered grid finite difference method

Karst cavities beneath bored cast in situ piles constitute substantial dangers to the stability and quality of underground construction projects. Therefore, it is important to detect and image these karst cavities before the construction of piles for such projects. New sonar detection methods, which have advantages such as low cost and high efficiency, can ensure the quality of the bedrock for underground construction projects. In this study, a three-dimensional (3D) hole-bedrock-cave model for the sonar detection of karst cavities was developed using the finite difference time domain (FDTD) method. The modeling parameters were calibrated with measured data, and the numerical results were then compared with measured field detection data to validate the reliability of the proposed method. The numerical results revealed the following findings. (1) The waves transmitted to the center of the pile hole bottom can be divided into three types: the waves traveling upward into the slurry that are absorbed by a perfectly matched layer (PML) and can be ignored, the waves traveling horizontally along the hole bottom that constitute the primary interference during the detection process, and the waves traveling downward into the bedrock that represent the effective detection signal. To minimize the negative effects of multiple reflections from the pile hole wall, we deployed a receiver at the midpoint between the transmission point and the hole wall, at which location the signals from the multiple reflections are symmetrical and easy to identify. (2) The reflected signals generated from different depths, sizes and directions of cavities have different travel time features that make it possible to estimate the depth, size and direction of a cavity beneath a pile. (3) The velocities of multiple surface wave reflections can be predicted in the frequency spectrum of a test signal and then used to predict the arrival times of multiple surface wave reflections and identify reflected P-waves. (4) Cavities have a focusing effect on waves when they are completely filled with a low-velocity medium. The energy of a wave reflected from the floor of a karst cavity after focusing is theoretically stronger than that reflected from the roof of a karst cavity. (5) PMLs can effectively absorb waves reflected from artificial boundaries to avoid false interference. High-order FDTD methods can eliminate numerical dispersion and reduce the computational costs. In general, the sonar detection of karst cavities beneath bored piles can be simulated using a 3D high-order staggered grid finite difference method and PMLs. This numerical modeling scheme is reliable and can improve the accuracy and feasibility of practical detection experiments.

Spectral ADI‐HO‐FDTD method for periodic structures

In this Letter, the spectral technique that has been applied to the finite-difference-time-domain (FDTD) method is applied to the alternating-direction-implicit higher-order FDTD (ADI-HO-FDTD) method for solving periodic structures, resulting in the periodic spectral (PS) ADI-HO-FDTD. The algorithm for inversion of the block periodic tri-diagonal matrices is presented. The effectiveness and computational efficiency of the proposed method are also shown. Numerical results indicate that the PS-ADI-HO-FDTD algorithm can save computation time by more than 27% in comparison with the conventional HO-FDTD method.

High-order accurate FDTD schemes for dispersive Maxwell’s equations in second-order form using recursive convolutions

Abstract We propose a novel finite-difference time-domain (FDTD) scheme for the solution of Maxwell’s equations in which linear dispersive effects are present. The method uses high-order accurate approximations in space and time for the dispersive Maxwell’s equations written as a second-order vector wave equation with a time-history convolution term. The modified equation approach is combined with the recursive convolution (RC) method to develop high-order approximations accurate to any desired order in space and time. High-order-accurate centered approximations of the physical Maxwell interface conditions are derived for the dispersive setting in order to fully restore accuracy at discontinuous material interfaces. Second- and fourth-order accurate versions of the scheme are presented and implemented in two spatial dimensions for the case of the Drude linear dispersion model. The stability of these schemes is analyzed. Finally, our approach is also amenable to curvilinear numerical grids if used with an appropriate generalized Laplace operator.

Reverse Time Migration Using the Pseudospectral Time-Domain Algorithm

We propose a pre-stack reverse time migration (RTM) seismic imaging method using the pseudospectral time-domain (PSTD) algorithm. Traditional pseudospectral method uses the fast Fourier transform (FFT) algorithm to calculate the spatial derivatives, but is limited by the wraparound effect due to the periodicity assumed in the FFT. The PSTD algorithm combines the pseudospectral method with a perfectly matched layer (PML) for acoustic waves. PML is a highly effective absorbing boundary condition that can eliminate the wraparound effect. It enables a wide application of the pseudospectral method to complex models. RTM based on the PSTD algorithm has advantages in the computational efficiency compared to traditional methods such as the second-order and high order finite difference time-domain (FDTD) methods. In this work, we implement the PSTD algorithm for acoustic wave equation based RTM. By applying the PSTD-RTM method to various seismic models and comparing it with RTM based on the eighth-order FDTD method, we find that PSTD-RTM method has better performance and saves more than 50% memory. The method is suitable for parallel computation, and has been accelerated by general purpose graphics processing unit.

A Hybrid Method of Higher-Order FDTD and Subgridding Technique

A hybrid method of higher-order finite-difference time-domain (FDTD) method and subgridding technique is proposed in this letter. The whole computational domain is decomposed into two kinds of regions. One is coarse-grid higher-order region, and another is fine-grid low-order region. An interlayer is introduced to connect the two kinds of regions. The update schemes for fields in each region and in the interlayer are given. Different time-steps are used in the higher-order region and in the subgridding region to further improve computational efficiency. Numerical properties in terms of nonphysical interface reflection and late-time stability are numerically studied. Numerical examples are also given to validate the efficiency and accuracy of this method.

Implementation and application of ADE‐PML for higher‐order FDTD method

A new implementation of a perfect matched layer with an auxiliary differential equation (ADE-PML) is presented for the higher-order finite-difference time-domain (FDTD) method. The implementation of the ADE-PML is more straightforward and more memory saving than that of the traditional PML. A two-dimensional (2D) case and two 3D cases with different microwave circuit components are designed to analyse the absorbing performance of the ADE-PML and the simulation results show that the ADE-PML can provide a quite satisfactory absorbing boundary condition.

Improved Hybrid FDTD Method for Studying Tunable Graphene Frequency-Selective Surfaces (GFSS) for THz-Wave Applications

One hybrid finite-difference time-domain (FDTD) method is proposed for studying transmission characteristics of a THz wave through some tunable dispersive graphene frequency-selective surfaces (GFSSs) biased by an electrostatic or a magnetostatic field, respectively. It integrates auxiliary differential equation (ADE)-FDTD with high-order conformal FDTD (2, 4) so as to handle atomically thin periodic structures on an electrically anisotropic substrate with high flexibility and accuracy. The dispersion error and computational efficiency of the developed algorithm is validated in comparison with the piecewise linear recursive convolution (PLRC)-FDTD and commercial software HFSS, respectively. Numerical studies are further performed to demonstrate its capability for characterizing the effects of chemical potential of graphene, biasing electrostatic and magnetostatic fields, and geometrical parameters of different GFSSs on the central frequencies and 3-dB bandwidths of their passbands and stopbands. In particular, both frequency selectivity and tunability of GFSSs made of single- and double-layer periodic split-ring resonators are predicted, which can be exploited for developing some novel tunable THz structures.